0

我尝试在雅加达湾进行克里金法。我有一组具有适当坐标和属性(pH、盐度、...)的测量点

为了进行克里金法,我首先需要为我的变异函数找到一个模型。当我使用“变异函数”函数时,输出并不完美,但应该没问题,但是当我尝试拟合变异函数时,我收到一条警告消息:In fit.variogram(ph.vgm, model = vgm(0.12, " Sph", 0.1, 0.01)) :警告:变异函数拟合中的奇异模型,我有一个奇异模型。

在这里,我阅读了与变异函数计算相关的奇异模型。我可以做点什么让它变得更好吗?

我怎样才能更好地拟合我的变异函数?为什么我只能在测量点周围获得小圆圈?我想拥有带有预测值的完整地图。

我还尝试了更不灵活的“automap”库,我没有获得好的结果。

 library(sp)
library(gstat)
library(automap)

x = c(11878417.51,11882987.17,11887690.42,11892582.91,11897119.18,11902527.08,11879348.14,11884237.29,11888933.86,11893819.67,11898835.73,11903940.84,11908386.94,11885529.71,11889836.66,11900118.13,11905765.37,11896037.16,11901234.67,11906244.04,11892136.86,11900822.56,11904493.1,11907692.42,11910346.05,11888709,11887268.41,11885237.28,11883450.38,11880668.5)
y = c(-668537.7429,-667290.838,-666043.9586,-663943.1247,-663992.3709,-662612.3726,-672878.6036,-672014.4364,-671960.7062,-669604.4601,-668541.1009,-667203.5333,-666181.6289,-676933.1896,-676566.0044,-673095.7667,-671736.8309,-679340.0992,-677788.4711,-676606.3051,-682542.446,-680607.5158,-680131.1539,-679733.0503,-662307.2774,-680754.1755,-681408.3272,-680494.7783,-680491.4197,-679426.19)
ph = c(7.1,7.76,7.14,7.19,7.56,7.56,7.11,8.14,7.22,7.17,7.33,7.37,7.36,7.23,7.12,7.54,7.96,7.98,7.96,7.2,7.44,7.36,7.71,7.71,8.01,7.73,8.11,7.03,7.26,7.77)
TSS = c(13.7,21,17.7,18.8,4.7,12.4,17.3,18.8,20.2,18.3,5.6,NA,NA,NA,21.9,11.1,NA,NA,21.2,29.1,31.3,29.3,21.3,25.4,31.8,14.5,2.9,11.7,8.4,NA)

df = data.frame(x,y,ph,TSS)
coordinates(df) = ~x+y
proj4string(df) <- CRS("+init=epsg:3857")
spplot(df)

grid <- data.frame(x=c(11877328.43,11879828.43,11882328.43,11884828.43,11887328.43,11889828.43,11892328.43,11894828.43,11897328.43,11899828.43,11902328.43,11904828.43,11907328.43,11909828.43,11912328.43,11877328.43,11879828.43,11882328.43,11884828.43,11887328.43,11889828.43,11892328.43,11894828.43,11897328.43,11899828.43,11902328.43,11904828.43,11907328.43,11877328.43,11879828.43,11882328.43,11884828.43,11887328.43,11889828.43,11892328.43,11894828.43,11897328.43,11899828.43,11902328.43,11904828.43,11907328.43,11909828.43,11912328.43,11877328.43,11879828.43,11882328.43,11884828.43,11887328.43,11889828.43,11892328.43,11894828.43,11897328.43,11899828.43,11902328.43,11904828.43,11907328.43,11909828.43,11912328.43,11877328.43,11879828.43,11882328.43,11884828.43,11887328.43,11889828.43,11892328.43,11894828.43,11897328.43,11899828.43,11902328.43,11904828.43,11907328.43,11909828.43,11879828.43,11882328.43,11884828.43,11887328.43,11889828.43,11892328.43,11894828.43,11897328.43,11899828.43,11902328.43,11904828.43,11907328.43,11909828.43,11912328.43,11879828.43,11882328.43,11884828.43,11887328.43,11889828.43,11892328.43,11894828.43,11897328.43,11899828.43,11902328.43,11904828.43,11907328.43,11879828.43,11882328.43,11884828.43,11887328.43,11889828.43,11892328.43,11894828.43,11897328.43,11899828.43,11902328.43,11904828.43,11907328.43,11909828.43,11884828.43,11887328.43,11889828.43,11892328.43,11894828.43,11897328.43,11899828.43,11904828.43,11894828.43),y=c(-659719.4518,-659719.4518,-659719.4518,-659719.4518,-659719.4518,-659719.4518,-659719.4518,-659719.4518,-659719.4518,-659719.4518,-659719.4518,-659719.4518,-659719.4518,-659719.4518,-659719.4518,-662219.4518,-662219.4518,-662219.4518,-662219.4518,-662219.4518,-662219.4518,-662219.4518,-662219.4518,-662219.4518,-662219.4518,-662219.4518,-662219.4518,-662219.4518,-664719.4518,-664719.4518,-664719.4518,-664719.4518,-664719.4518,-664719.4518,-664719.4518,-664719.4518,-664719.4518,-664719.4518,-664719.4518,-664719.4518,-664719.4518,-664719.4518,-664719.4518,-667219.4518,-667219.4518,-667219.4518,-667219.4518,-667219.4518,-667219.4518,-667219.4518,-667219.4518,-667219.4518,-667219.4518,-667219.4518,-667219.4518,-667219.4518,-667219.4518,-667219.4518,-669719.4518,-669719.4518,-669719.4518,-669719.4518,-669719.4518,-669719.4518,-669719.4518,-669719.4518,-669719.4518,-669719.4518,-669719.4518,-669719.4518,-669719.4518,-669719.4518,-672219.4518,-672219.4518,-672219.4518,-672219.4518,-672219.4518,-672219.4518,-672219.4518,-672219.4518,-672219.4518,-672219.4518,-672219.4518,-672219.4518,-672219.4518,-672219.4518,-674719.4518,-674719.4518,-674719.4518,-674719.4518,-674719.4518,-674719.4518,-674719.4518,-674719.4518,-674719.4518,-674719.4518,-674719.4518,-674719.4518,-677219.4518,-677219.4518,-677219.4518,-677219.4518,-677219.4518,-677219.4518,-677219.4518,-677219.4518,-677219.4518,-677219.4518,-677219.4518,-677219.4518,-677219.4518,-679719.4518,-679719.4518,-679719.4518,-679719.4518,-679719.4518,-679719.4518,-679719.4518,-679719.4518,-682219.4518))

coordinates(grid)=~x+y
proj4string(grid) <- CRS("+init=epsg:3857")
gridded(grid)=T
spplot(grid)

ph.vgm <- variogram(ph~1, df[!is.na(df@data$ph),]); plot(ph.vgm)
ph.fit = fit.variogram(ph.vgm, model = vgm(0.12, "Sph", 4000, 0.01), warn.if.neg = T); ph.fit
plot(ph.vgm, ph.fit)
ph.kriged = krige(ph~1, df[!is.na(df@data$ph),], grid, model = ph.fit)
spplot(ph.kriged["var1.pred"])
4

1 回答 1

0

和小样本。

如果您想查看每个点对的 $0.5(Z(x)-Z(x+h))^2$,针对 $h$ 绘制,然后使用

plot(variogram(ph~1, df[!is.na(df@data$ph),], cloud=TRUE))

但这也可能不会让你开心。最重要的是,您的观察结果非常少 (30),而且观察结果如此之少,您基本上永远不会得到看起来不错的变异函数。

于 2015-09-04T17:08:18.377 回答