带有 autoKrige() 的 R 通用克里金法
R universal kriging with autoKrige()
我正在尝试使用 automap 包中的 autoKrige() 函数来简单应用通用克里金法。我有一个不规则间隔的测量网格,我想在它们之间进行精细空间尺度的插值。示例代码:
library('automap')
# create an irregularly spaced grid
y <-x <-c(-5,-4,-2,-1,-0.5,0,0.5,1,2,4,5)
grid <-expand.grid(x,y)
names(grid) <-c('x', 'y')
# create some measurements, greatest in the centre, with some noise
vals <-apply(grid,1, function(x) {12/(0.1+sqrt(x[1]^2 + x[2]^2))+rnorm(1,2,1.5)})
# get data into sp format
s <-SpatialPointsDataFrame(grid, data.frame(vals))
# make some prediction locations and get them into sp format
pred <-expand.grid(seq(-5,5,by=0.5), seq(-5,5,by=0.5))
pred <-cbind(pred[,1], pred[,2]) # this seems to be needed, not sure why
pred <-SpatialPoints(pred)
# try universal kriging
surf <-autoKrige(vals~x+y, s, new_data=pred)
这会导致错误:
Error in gstat.formula.predict(d$formula, newdata, na.action = na.action, :
NROW(locs) != NROW(X): this should not occur
我试过让new_data和原始数据有相同的行数,甚至试过让new_data中的坐标和原始数据完全一样,但是我仍然出现此错误。我是地统计学技术的新手,如果我犯了一个基本错误,我深表歉意。谁能告诉我哪里出错了?谢谢
问题是 autoKrige
函数的语法错误。 autoKrige
的公式输入指定您要使用的线性模型,例如:
log(zinc) ~ dist
来自 meuse
数据集。在这种情况下,您使用线性模型对 log(zinc)
与 dist
进行建模,并使用变差函数对该模型的残差进行插值。本质上,通用克里金法是具有空间相关残差的线性回归。
在您的情况下,您指定:
val ~ x+y
so autoKrige(gstat
实际上)将尝试首先对 vals
与 x
和 y
(多元回归)的线性模型进行建模,然后对残差进行插值使用变异函数模型。但是,x
和 y
变量不在 SpatialPointsDataFrame
中。
我认为您想做的是仅使用变差函数模型进行空间插值。那样的话,线性模型就很简单了,其实只是拟合一个平均值:
vals ~ 1
其中确定了 vals
的平均值,并使用变差函数模型对残差进行插值。这实际上称为普通克里金法。您对 autoKrige
的调用类似于:
surf <-autoKrige(vals ~ 1, s, new_data=pred)
我正在尝试使用 automap 包中的 autoKrige() 函数来简单应用通用克里金法。我有一个不规则间隔的测量网格,我想在它们之间进行精细空间尺度的插值。示例代码:
library('automap')
# create an irregularly spaced grid
y <-x <-c(-5,-4,-2,-1,-0.5,0,0.5,1,2,4,5)
grid <-expand.grid(x,y)
names(grid) <-c('x', 'y')
# create some measurements, greatest in the centre, with some noise
vals <-apply(grid,1, function(x) {12/(0.1+sqrt(x[1]^2 + x[2]^2))+rnorm(1,2,1.5)})
# get data into sp format
s <-SpatialPointsDataFrame(grid, data.frame(vals))
# make some prediction locations and get them into sp format
pred <-expand.grid(seq(-5,5,by=0.5), seq(-5,5,by=0.5))
pred <-cbind(pred[,1], pred[,2]) # this seems to be needed, not sure why
pred <-SpatialPoints(pred)
# try universal kriging
surf <-autoKrige(vals~x+y, s, new_data=pred)
这会导致错误:
Error in gstat.formula.predict(d$formula, newdata, na.action = na.action, :
NROW(locs) != NROW(X): this should not occur
我试过让new_data和原始数据有相同的行数,甚至试过让new_data中的坐标和原始数据完全一样,但是我仍然出现此错误。我是地统计学技术的新手,如果我犯了一个基本错误,我深表歉意。谁能告诉我哪里出错了?谢谢
问题是 autoKrige
函数的语法错误。 autoKrige
的公式输入指定您要使用的线性模型,例如:
log(zinc) ~ dist
来自 meuse
数据集。在这种情况下,您使用线性模型对 log(zinc)
与 dist
进行建模,并使用变差函数对该模型的残差进行插值。本质上,通用克里金法是具有空间相关残差的线性回归。
在您的情况下,您指定:
val ~ x+y
so autoKrige(gstat
实际上)将尝试首先对 vals
与 x
和 y
(多元回归)的线性模型进行建模,然后对残差进行插值使用变异函数模型。但是,x
和 y
变量不在 SpatialPointsDataFrame
中。
我认为您想做的是仅使用变差函数模型进行空间插值。那样的话,线性模型就很简单了,其实只是拟合一个平均值:
vals ~ 1
其中确定了 vals
的平均值,并使用变差函数模型对残差进行插值。这实际上称为普通克里金法。您对 autoKrige
的调用类似于:
surf <-autoKrige(vals ~ 1, s, new_data=pred)