Передать функцию dataframe в аргумент формулы gstat - PullRequest
Новая
       98

Передать функцию dataframe в аргумент формулы gstat

0 голосов
/ 05 октября 2018

Я пытаюсь создать растры с обратным расстоянием, используя gstat() и raster::interpolate().У меня проблемы с передачей имени столбца в аргумент формулы функции gstat.

Жесткое кодирование имени столбца работает просто отлично: 

gs <- gstat(formula=v1~1, locations = data)
r <- raster(shape, res=1000, crs = crs(data))
idw <- raster::interpolate(r, gs)
idwr <- mask(r, gs)
plot(idwr)

enter image description here

Оборачивая это в функцию, чтобы я мог просмотреть несколько имен столбцов, выдается ошибка: 

apply_gstat <- function(col_name, data = data, shape = shape) {
    gs <- gstat(formula=col_name~1, locations = data)
    r <- raster(shape, res=1000, crs = crs(data))
    ## interpolate() throws an error because of issue with gstat
    idw <- raster::interpolate(r, gs)
    idwr <- mask(r, gs)
    plot(idwr)
}

col_names <- c("v1", "v2", "v3")
lapply(col_names, function(x) {
    gstat_apply(col_name = x, data = data, shape = shape)
}


Error in predict.gstat(model, blockvals, debug.level = debug.level, ...) : 
  too many spatial dimensions: 18
In addition: Warning message:
In predict.gstat(model, blockvals, debug.level = debug.level, ...) :
  NAs introduced by coercion

1 Ответ

0 голосов
/ 06 октября 2018

Вероятно, это будет работать, если вы используете as.formula(paste(col_name, '~1')) Вот отдельный пример: 

library(raster)
library(gstat)
data(meuse)
r <- raster(system.file("external/test.grd", package="raster"))

gfun <- function(v) {
    mg <- gstat(formula = as.formula(paste(v, "~1")), locations= ~x+y, data=meuse, nmax=7, set=list(idp = .5))
    z <- interpolate(r, mg)
    mask(z, r)
}

vars <- c('cadmium', 'copper', 'lead', 'zinc')
x <- sapply(vars, gfun)
s <- stack(x)
s
#class       : RasterStack 
#dimensions  : 115, 80, 9200, 4  (nrow, ncol, ncell, nlayers)
#resolution  : 40, 40  (x, y)
#extent      : 178400, 181600, 329400, 334000  (xmin, xmax, ymin, ymax)
#coord. ref. : +init=epsg:28992 +towgs84=565.237,50.0087,465.658,-0.406857,0.350733,-1.87035,4.0812 +proj=sterea +lat_0=52.15616055555555 +lon_0=5.38763888888889 +k=0.9999079 +x_0=155000 +y_0=463000 +ellps=bessel +units=m +no_defs 
#names       :     cadmium,      copper,        lead,        zinc 
#min values  :   0.3131223,  17.3430502,  42.3218862, 123.9369493 
#max values  :    11.91185,    97.28761,   470.70461,  1387.50615
...