自动克里金法和手动克里金法的不同结果
DIfferent outcomes from autokriging and manual kriging
有人可以帮助我理解为什么自动克里金法和手动克里金法会得到如此不同的结果吗?我看到这两种算法使用的是不同的变差函数模型,但这是否是造成所有差异的唯一原因?我也对 autokrige 在没有数据的位置预测显着更高的值这一事实感到不安,例如网格的左下角和右下角。它是否与通过取幂将(对数转换的)克里格输出转换回来有关?此外,这两种方法预测的值都比数据低得多。看来我对其中一个或两个算法都做错了,非常感谢任何帮助。
编辑:GitHub link to the data used here (I used this 共享 csv 文件的建议)。
library(sp)
library(gstat)
library(automap)
# get data
df <- read.csv("data.csv")
---------------------------
> head(df)
long lat n
1 -76.7116 39.2992 1
2 -76.7113 39.3553 1
3 -76.7113 39.3607 1
4 -76.7112 39.3126 2
5 -76.7110 39.2973 1
6 -76.7110 39.3364 1
---------------------------
# generate regular grid
coordinates(df) <- ~ long + lat
grd <- spsample(df, type = "regular", n = 10000)
colnames(grd@coords) <- c("long", "lat") # rename coords as in data
gridded(grd) <- TRUE
# manual kriging
form <- as.formula("log(n) ~ long + lat") # universal kriging formula
v <- variogram(object = form, locations = df)
vmfit <- fit.variogram(v,
model = vgm(model = c("Sph", "Exp", "Mat", "Gau", "Ste")))
-----------------------------------
> vmfit
model psill range kappa
1 Ste 0.2886451 0.001786237 0.5
-----------------------------------
krg <- krige(formula = form, locations = df, newdata = grd, model = vmfit)
krg_df <- data.frame(krg@coords, # extract kriged data
pred = exp(krg@data$var1.pred))
# auto kriging
krg2 <- autoKrige(formula = form, input_data = df, new_data = grd)
krg2_df <- data.frame(krg2$krige_output@coords,
pred = exp(krg2$krige_output@data$var1.pred))
-----------------------------------
> krg2$var_model
model psill range
1 Nug 0.26473893 0.00000000
2 Sph 0.02574092 0.01657061
-----------------------------------
我认为问题出在你的公式上。应该是n ~ 1
示例数据
library(gstat)
library(automap)
data(meuse)
coordinates(meuse) =~ x+y
data(meuse.grid)
gridded(meuse.grid) =~ x+y
解决方案:
form <- zinc ~ 1
# manual
v <- variogram(form, meuse)
m <- fit.variogram(v, model=vgm(model=c("Sph", "Exp", "Mat", "Gau", "Ste")))
km <- krige(form, meuse, meuse.grid, model = m)
# auto
ka <- autoKrige(form, meuse, meuse.grid)
看看:
g <- cbind(km[,1], ka[[1]][,1])
names(g) <- c("manual", "auto")
spplot(g)
上面罗伯特的回答给了我尝试公式 n ~ long + lat
的想法(而不是我之前使用的 log(n) ~ long + lat
)。现在,这两种方法似乎表现一致,差异合理。我想知道这是不是因为我的数据偏斜,70%的数据有n = 1
。在较早的公式中使用 log(n)
仅留下 30% non-zero 值,导致空间相关性非常小的稀疏场供 krige 使用。
上面例子中的meuse
数据没有这个问题:log(zinc) > 0
for all zinc
values.
有人可以帮助我理解为什么自动克里金法和手动克里金法会得到如此不同的结果吗?我看到这两种算法使用的是不同的变差函数模型,但这是否是造成所有差异的唯一原因?我也对 autokrige 在没有数据的位置预测显着更高的值这一事实感到不安,例如网格的左下角和右下角。它是否与通过取幂将(对数转换的)克里格输出转换回来有关?此外,这两种方法预测的值都比数据低得多。看来我对其中一个或两个算法都做错了,非常感谢任何帮助。
编辑:GitHub link to the data used here (I used this 共享 csv 文件的建议)。
library(sp)
library(gstat)
library(automap)
# get data
df <- read.csv("data.csv")
---------------------------
> head(df)
long lat n
1 -76.7116 39.2992 1
2 -76.7113 39.3553 1
3 -76.7113 39.3607 1
4 -76.7112 39.3126 2
5 -76.7110 39.2973 1
6 -76.7110 39.3364 1
---------------------------
# generate regular grid
coordinates(df) <- ~ long + lat
grd <- spsample(df, type = "regular", n = 10000)
colnames(grd@coords) <- c("long", "lat") # rename coords as in data
gridded(grd) <- TRUE
# manual kriging
form <- as.formula("log(n) ~ long + lat") # universal kriging formula
v <- variogram(object = form, locations = df)
vmfit <- fit.variogram(v,
model = vgm(model = c("Sph", "Exp", "Mat", "Gau", "Ste")))
-----------------------------------
> vmfit
model psill range kappa
1 Ste 0.2886451 0.001786237 0.5
-----------------------------------
krg <- krige(formula = form, locations = df, newdata = grd, model = vmfit)
krg_df <- data.frame(krg@coords, # extract kriged data
pred = exp(krg@data$var1.pred))
# auto kriging
krg2 <- autoKrige(formula = form, input_data = df, new_data = grd)
krg2_df <- data.frame(krg2$krige_output@coords,
pred = exp(krg2$krige_output@data$var1.pred))
-----------------------------------
> krg2$var_model
model psill range
1 Nug 0.26473893 0.00000000
2 Sph 0.02574092 0.01657061
-----------------------------------
我认为问题出在你的公式上。应该是n ~ 1
示例数据
library(gstat)
library(automap)
data(meuse)
coordinates(meuse) =~ x+y
data(meuse.grid)
gridded(meuse.grid) =~ x+y
解决方案:
form <- zinc ~ 1
# manual
v <- variogram(form, meuse)
m <- fit.variogram(v, model=vgm(model=c("Sph", "Exp", "Mat", "Gau", "Ste")))
km <- krige(form, meuse, meuse.grid, model = m)
# auto
ka <- autoKrige(form, meuse, meuse.grid)
看看:
g <- cbind(km[,1], ka[[1]][,1])
names(g) <- c("manual", "auto")
spplot(g)
上面罗伯特的回答给了我尝试公式 n ~ long + lat
的想法(而不是我之前使用的 log(n) ~ long + lat
)。现在,这两种方法似乎表现一致,差异合理。我想知道这是不是因为我的数据偏斜,70%的数据有n = 1
。在较早的公式中使用 log(n)
仅留下 30% non-zero 值,导致空间相关性非常小的稀疏场供 krige 使用。
meuse
数据没有这个问题:log(zinc) > 0
for all zinc
values.