Using a dynamic formula with parrallel kringing in R

Question

I am trying to achieve parallel Kriging in R on several variables using a loop. Here is a reproducible example using data meuse and a code for parallel kriging that I found here. In the code below, each kriging is overwritten, but that's not the issue, since I can't even run the first kriging.

library(parallel)
library(sp)
library(gstat)
data(meuse)
    names <- colnames(meuse)[3:6]
    coordinates(meuse) = ~x+y
    data(meuse.grid)
    gridded(meuse.grid) = ~x+y
    m <- vgm(.59, "Sph", 874, .04)
    no_cores <- 7
    # ordinary kriging:
    for (i in 1:length(names)) {
      parts <- split(x = 1:length(meuse.grid), f = 1:no_cores)
cl <- makeCluster(no_cores)          
clusterExport(cl = cl, varlist = c("meuse", "meuse.grid", "m", "parts"), envir = .GlobalEnv)
          clusterEvalQ(cl = cl, expr = c(library('sp'), library('gstat')))
          parallelX <- parLapply(cl = cl, X = 1:no_cores, fun = function(x) krige(formula = log(get(names[i]))~1, locations = meuse, newdata = meuse.grid[parts[[x]],], model = m))
          stopCluster(cl)
          # Merge all the predictions    
          mergeParallelX <- maptools::spRbind(parallelX[[1]], parallelX[[2]])
          mergeParallelX <- maptools::spRbind(mergeParallelX, parallelX[[3]])
          mergeParallelX <- maptools::spRbind(mergeParallelX, parallelX[[4]])
          mergeParallelX <- maptools::spRbind(mergeParallelX, parallelX[[5]])
          mergeParallelX <- maptools::spRbind(mergeParallelX, parallelX[[6]])
          mergeParallelX <- maptools::spRbind(mergeParallelX, parallelX[[7]])
          # Create SpatialPixelsDataFrame from mergeParallelX
          mergeParallelX <- SpatialPixelsDataFrame(points = mergeParallelX, data = mergeParallelX@data)
        }

I keep having issues which are, I think, related to the dynamic formula. I tried a lot of other functions, such as paste, paste0, as.formula, objects, etc. Nothing works, I can't paste the variable names into looped parallel kriging, no matter how I code the dynamic formula. Maybe that's because of how the parallel function is coded?

Any idea about a dynamic formula that works with the above code?

Answer 1

Guessing that your problem is constructing a formula dynamically, you can build one as character and convert to a formula with as.formula:

So this:

> k1 = krige(log(copper)~1, meuse, meuse.grid, model = m)
[using ordinary kriging]

is the same as:

> k2 = krige(as.formula("log(copper)~1"), meuse, meuse.grid, model = m)
[using ordinary kriging]

which you can construct

> names
[1] "cadmium" "copper"  "lead"    "zinc"   
> i
[1] 2
> k3 = krige(as.formula(paste0("log(",names[i],")~1")), meuse, meuse.grid, model = m)
[using ordinary kriging]

And then:

> identical(k1, k2)
[1] TRUE
> identical(k1, k3)
[1] TRUE
> identical(k2, k3)
[1] TRUE