0

可能有更好的方法，但是由于我是 R 新手并且已经设置了 IDW 代码，我一直在尝试通过调整 IDW 代码来获得 2000 米内所有点的中值，设置权重（idp ) 接近零，因此较近的点的权重与较远的点相同。

我猜当我用 maxdist=2000 运行下面的代码时它会说 NA 因为有些点在 2000 米内没有任何邻居。即使我将 nmin 设置为零，我可以让它使用的最小 maxdist 约为 40,000。

有没有办法告诉它忽略 2000 米内没有邻居的点，或者有人知道更好的方法吗？

这是我的代码：

library(gstat)
clean3145 = read.csv("clean3145.csv")

#Set up the k-fold validation
set.seed(88)
groups <- sample(1:5, nrow(clean3145), replace=TRUE)

#res=result=R=Pearson's correlation between predicted and actual arsenic concentration
MEDres<- rep(NA, 5)

r <- list()
for (k in 1:5) {
  print(k)
  flush.console()
  train <- clean3145[groups!=k, ]
  test <- clean3145[groups==k, ]

  med <- gstat(formula = As1~1, locations = ~UTMNM+UTMEM, data=train, nmin=0, maxdist=40000, set=list(idp = .01))
  medpred <- predict(med, test)$var1.pred
  MEDres[k] <- cor(test$As1, medpred)

  }

#Show the mean correlation for the 5 different training-test dataset pairs in K-fold validation
mean(MEDres)

谢谢你的帮助！

4

2 回答 2

0

谢谢埃泽！

我会保存它以供将来参考。我们让它以这种方式工作，也有深度标准（我正在尝试估计地下水中的砷）：

#Load required packages and data
library(raster)
depth = read.csv("depth.csv")

设置 k 折验证，确保每次选择相同的随机样本以进行可比性

set.seed(88)
groups <- sample(1:5, nrow(depth), replace=TRUE)

使用以米为单位的 UTM 东和北坐标 (UTMEM、UTMNM) 计算测试 (tst) 井的某个点距离 (pd) 内所有训练 (trn) 井的中值砷浓度。忽略或“移除”在 148 米内没有邻居的测试井 (pd>148=NA, na.rm=TRUE)

  computeMed <- function(trn, tst) {
  pd <- pointDistance(trn[ , c('UTMEM', 'UTMNM')], tst[ , c('UTMEM','UTMNM')], lonlat=FALSE)  
  pd[pd > 148] <- NA

  as <- trn$As1
  as <- matrix(rep(as, ncol(pd)), ncol=ncol(pd))
  aspd <- as * (pd >= 0)
  apply(aspd, 2, median, na.rm=TRUE)

    }

再次计算中位数，这次使用深度标准（例如，如果测试井靠近 Fallon（Tcan2car=1=从 Truckee Canal 到 Carson Basin 的井和下坡度）并且深度超过 40 m，则仅给出也 > 40 的邻居的中位数米深）

r <- rd <- list()
Fallon <- FALSE
for (k in 1:5) {
  print(k)
  flush.console()
  depth$deep <- TRUE
  depth$deep[depth$Depth_m < 40] <- FALSE
  if (Fallon) {
    d  <- depth[depth$Tcan2car==1]
   } else {
     d <- depth
  }
  train <- d[groups!=k, ]
  test <- d[groups==k, ]

   p <- computeMed(train,test)
   r[[k]] <- cbind(k=k, prd=p, obs=test$As1)


   pdeep <- computeMed(train[train$deep,],test[test$deep,])
   pshallow <- computeMed(train[!train$deep,],test[!test$deep,])


  rd[[k]] <- cbind(k=k, prd=c(pdeep, pshallow),    obs=test$As1[c(which(test$deep), which(!test$deep))])

 }

显示 K 折验证中 5 个不同训练-测试数据集对的平均 Pearson 相关性。cr 和 r 仅指基于距离的相关性。crd 和 rd refer 还包括深度标准

cr <- sapply(r, function(x) {x <- na.omit(x); cor(x[,2:3])[2]})
cr
mean(cr)

crd <- sapply(rd, function(x) {x <- na.omit(x); cor(x[,2:3])[2]})
crd
mean(crd)
于 2015-01-19T17:06:23.353 回答
0

我看不出您的代码如何帮助回答您的原始问题，但对于当地的中位数，我会尝试

library(sp)
demo(meuse, ask = FALSE)
library(gstat)
x = krige(zinc~1, meuse, meuse.grid, maxdist = 1000, set = list(method = "med"))

如果邻域不包含数据，您可以通过最近点的数量来定义它nmax，在这种情况下，当然不再控制距离。

于 2015-01-18T22:42:36.660 回答