如何使用gstat predict进行普通克里金插值

Question

我正在尝试用R语言编写代码，使用gstat库进行插值。 我已经阅读了gstat手册，并基于互联网上的一些示例成功编写了以下代码（这仅是部分内容）：

 g <- gstat(id="tec", formula=TEC ~ 1, data=data)  ##I create an object
 v <- variogram(g) # plot the empirical variogram
 plot(v)
 mod<-vgm(sill=var(data$TEC),model="Sph",range=200,nugget=200) #create the variogram model

v.fit <- fit.variogram(v, model=mod,fit.method=1)  #fit the empirical variogram 
Theor_variogram=plot(variogram(g),v.fit,main="WLS Model") #plot the theoretical variogram
plot(Theor_variogram)
 ## Kriging interpolation
 p <- predict.gstat(g, model=v.fit, newdata=predGrid)

我的问题是，当我运行最后一个命令（predict）时，我得到的结果不是普通克里金插值，而是反距离加权（IDW）插值。 我在gstat手册中读到：“当没有指定变差函数时，使用反距离加权插值作为默认操作。当指定了变差函数时，默认预测方法是普通克里金法。” 但是，正如您在我的代码中所看到的，我指定了经验和理论变差函数。 你知道我为什么一直得到IDW而不是普通克里金吗？这是否与我拥有的坐标类型有关？例如，如果我有彼此接近的坐标或者感兴趣的区域太大？ 任何帮助都将非常有用。 提前谢谢 Dimitris

- user1919

Answer 1

在预测阶段之前，您需要包括gstat对象的创建：

g <- gstat(id="tec", formula=TEC ~ 1, data=data, model = v.fit)

然而，我建议使用标准接口krige来使用gstat。这将把构建gstat对象和预测结合在一个函数中。很少需要自己构建gstat对象。例如：

data(meuse)
coordinates(meuse) = ~x+y
data(meuse.grid)
gridded(meuse.grid) = ~x+y
m <- vgm(.59, "Sph", 874, .04) 
# OK:
x <- krige(log(zinc)~1, meuse, meuse.grid, model = m)

你也可以使用我编写的automap包，让变差函数模型自动适应数据。例如，使用meuse数据集：

library(automap)
kr = autoKrige(log(zinc)~1, meuse, meuse.grid)

这将自动构建一个样本半方差函数，并对该样本半方差函数拟合一个半方差函数模型。

- Paul Hiemstra

---

谢谢你的回答。Autokrige函数听起来非常有趣。在使用它时，您可以指定拟合方法（例如WLS、OLS）吗？此外，它使用什么初始值来计算 sill/range/nugget？无论如何，我会尝试找到一些相关文档。谢谢。 - user1919