如何在 R 中使用 gstat 拟合具有半变异函数的模型？

Question

我有一个包含 3 月 1 日下午 12 点的大气 PM10 浓度数据的 csv 文件。please, download.我想在 R 中使用 gstat 包绘制半变异函数。我试图在 R 中编写这些代码。但是有了这些数据，我无法拟合模型。

    library(sp)
    library(gstat)

    seoul3112<-read.csv("seoul3112.csv")
    seoul3112<-na.omit(seoul3112)

    g<-gstat(id="PM10",formula=PM10~LON+LAT,location=~LON+LAT,
             data=seoul3112)
    seoul3112.var<-variogram(g,width=0.04,cutoff=0.6)
    seoul3112.var
    plot(seoul3112.var, col="black", pch=16,cex=1.3,
         xlab="Distance",ylab="Semivariance",
         main="Omnidirectional Variogram for seoul 3112")

model.3112<- fit.variogram(seoul3112.var,vgm(700,"Gau",0.5,200), fit.method = 2)
    plot(seoul3112.var,model=model.3112, col="black", pch=16,cex=1.3,
         xlab="Distance",ylab="Semivariance",
         main="Omnidirectional Variogram for seoul 3112")

其实我是R和统计学的初学者。所以，我什至对变差函数一无所知。我有一些疑问：

a) 当我将数据绘制为半变异函数时，它看起来与典型的半变异函数不同！为什么会这样？我应该对我的数据做任何其他事情吗？

b)我怎样才能用这些数据拟合模型？我尝试过不同的模型，如 "Sph"、"Exp"，但它们看起来像线性的！为什么？

c)我如何理解我应该在 vgm() 函数中使用什么 sill,range,nugget 的初始值？

d)我如何理解模型与数据的正确拟合？

e)为了使用克里金法，我应该绘制什么样的半变异函数？只有全向半变异函数？或者我应该绘制方向半变异函数？

f)我该如何解释半变异函数？我的意思是我对半变异函数的数据实际上能理解什么？

提前致谢。

Answer 1

您的坐标是纬度和经度，但您没有告知 gstat 它们是。因此，gstat 将假设它可以根据这些数字计算欧几里德距离，这是没有意义的。

建议是在使用包 sp 将您的点转换为 SpatialPointsDataFrame 之后学习如何使用 gstat，然后学习如何投影您的数据以使欧几里得距离有意义。

Answer 2

我会回答您的代码相关问题。您的其余问题（d、e 和 f）更多地与理论相关。

首先，在您的评论中，当您更改 proj4string 时，绘图上的距离单位应该已经更改。他们做了吗？根据您的评论，这听起来好像并没有发生。

a) 除了玩转 cutoff 距离外，还要注意支持半变异函数上每个 bin 的 np（点对）。例如，使用您更新的 proj4string 信息，我尝试了 cutoff=80 和 width=80/10（10 个 bin 而不是 15 个）来查看半变异函数的形状如何变化。从 15 个减少到 10 个 bin 不会改变点对存在的位置，只会增加每个 bin 代表的距离。此外，这种方法不一定是您应该使用的方法，但它是一个示例，说明如何更改 bin 以获得更平滑的样本半变异函数（但更平滑并不意味着更好）。

b) 使用您的代码，"Sph" 和 "Exp" 模型 return Warning: singular model in variogram fit。此警告表明没有足够的数据来拟合球形和指数经验模型的某些参数。有关每个经验方程及其参数的指导，请参阅 gstat user manual。

c) 例如，vgm() 函数可用于眼睛拟合样本半变异函数。如果您对如何使用样本数据绘制 vgm() 模型感到困惑，请尝试

eye_vgm = vgm(psill=1200,model="Gau",range=60,nugget=350)
plot(seoul3112.var,model=eye_vgm, col="black", pch=16,cex=1.3)

您在 fit.variogram() 的调用中使用了 vgm()，因此只要您提供给 vgm() 的参数是合理的（例如，基于示例数据）并且经验模型可以有参数拟合，fit.variogram() 会根据 fit.method.

找到一个拟合

如何在 R 中使用 gstat 拟合具有半变异函数的模型？

How to fit model with semivariogram using gstat in R?

r

spatial

gstat

kriging

covariogram