use*_*703 6 interpolation r geospatial
我有一个网格，其中包含我要使用插值法填充的间隙（NA）。我的网格在x和y维度上显示自相关，因此我想尝试双线性插值。我发现的大多数解决方案都集中在“上采样”（为了增加样本数量/网格大小而进行插值）上，但我不想/不需要更改网格大小。我只想使用插值填充NA。其他潜在解决方案似乎无法处理输入值网格（“ z矩阵”）的NA，或者是基于邻域的解决方案，而不是双线性插值法，或者根本没有答案。
我发现，使用栅格数据包，我可以输入一个包含NA的栅格（作为栅格），并使用“ resample”命令输出相同大小的栅格。但是，结果看起来像最近邻插值，而不是双线性插值。
我是否缺少某种方法，使得可以使用栅格数据包进行双线性插值？还是有一种更好的方法来简单地填充NA来进行双线性插值？
library(raster)
# raster containing gap
r <- raster(nrow=10, ncol=10)
r[] <- 1:ncell(r)
r[25] <- NA
# The s raster is the same size as the r raster
s <- raster(nrow=10, ncol=10)
s <- resample(r, s, method='bilinear')
plot(r)
plot(s)
s[25]
s[35]
# s[25] appears to have been filled with neighbor s[35]
更新
Akima软件包似乎是上述栅格方法的一个有前途的替代方法，但是如果值的输入网格（Z矩阵）中存在NA，我将遇到麻烦。这是一个与上述示例相似的示例。（再次，我将插入到与原始尺寸相同的网格中）。
library(akima)
# Use bilinear interpolation (no NAs in input)
rmat<-matrix(seq(1,100,1), nrow = 10, ncol = 10, byrow = T)
x <- seq(1,10,1)
y <- seq(1,10,1)
smat <- bilinear.grid(x, y, rmat, nx = 10, ny = 10) # works
plot(raster(rmat), main = "original")
plot(raster(smat$z), main = "interpolated")
# Try using bilinear interpolation but with an NA
rmat<-matrix(seq(1,100,1), nrow = 10, ncol = 10, byrow = T)
rmat[3,5] <- NA
x <- seq(1,10,1)
y <- seq(1,10,1)
smat <- bilinear.grid(x, y, rmat, nx = 10, ny = 10) # Error about NAs
更新2
@Robert Hijmans提出了一个很大的问题，即为什么不在栅格数据包中的focus（）命令中使用移动窗口平均值。原因是我想尝试双线性插值，而且我不认为移动窗口平均值总是能提供与双线性插值相同的答案。但是，在我发布的示例中并不清楚（在该示例中，移动窗口和双线性插值确实给出了相同的答案），因此我将在下面的新示例中进行演示。请注意，以下示例的双线性插值解应为8（这是用于测试的便捷计算器）。
library(raster)
r <- raster(nrow=10, ncol=10)
# Different grid values than earlier examples
values(r) <- c(rep(1:5, 4), rep(4:8, 4), rep(1:5, 4), rep(4:8, 4), rep(1:5, 4))
r[25] <- NA
plot(r)
# See what the mean of the moving window produces
f <- focal(r, w=matrix(1,nrow=3, ncol=3), fun=mean, NAonly=TRUE, na.rm=TRUE)
f[25] # Moving window gives 5 but bilinear interp gives 8
# Note that this seems to be how the moving window works with equal weights
window_test <- c(r[14:16], r[24:26], r[34:36])
mean(window_test, na.rm = T)
我在这里想念什么吗？也许focus（）的权重参数有些聪明，可以产生双线性插值解？
让我们使用等距像元来避免因像元大小随 lon/lat 数据而变化的差异
library(raster)
r <- raster(nrow=10, ncol=10, crs='+proj=utm +zone=1 +datum=WGS84', xmn=0, xmx=1, ymn=0, ymx=1)
对于此示例，您可能会使用 focal
values(r) <- 1:ncell(r)
r[25] <- NA
f <- focal(r, w=matrix(1,nrow=3, ncol=3), fun=mean, NAonly=TRUE, na.rm=TRUE)
我看到您忽略了“基于邻域的解决方案而不是双线性插值”。但问题是为什么。在这种情况下，您可能需要基于邻域的解决方案。
更新。再说一次，在单元格不是近似正方形的情况下，双线性将是可取的。
values(r) <- c(rep(1:5, 4), rep(4:8, 4), rep(1:5, 4), rep(4:8, 4), rep(1:5, 4))
r[25] <- NA
双线性插值的问题通常使用 4 个连续的单元格，但在这种情况下，如果您想要单元格中心的值，适当的单元格将是单元格本身的值，因为到该单元格的距离为零，并且因此这就是插值结束的地方。例如，对于单元格 23
extract(r, xyFromCell(r, 23))
#6
extract(r, xyFromCell(r, 23), method='bilinear')
#[1] 6
在这种情况下，焦点单元格为 NA，因此您将获得焦点单元格和另外 3 个单元格的平均值。问题是哪三个？它是任意的，但要使其工作，NA 单元格必须获得一个值。该raster算法将 NA 单元格下方的值分配给该单元格（此处也是 8）。我认为，这可以很好地处理边缘（例如陆地/海洋）的 NA 值，但在这种情况下可能不是。`extract(r, xyFromCell(r, 25)) #NA extract(r, xyFromCell(r, 25), method='bilinear') #[1] 8
这也是什么resample给
resample(r, r)[25]
# 8
这也是在线计算器所建议的吗？
这对小的变化非常敏感
extract(r, xyFromCell(r, 25)+0.0001, method='bilinear')
#[1] 4.998997
在这种情况下，我真正想要的是车邻居的平均值
mean(r[adjacent(r, 25, pairs=FALSE)])
[1] 6
或者，更一般地，局部反距离加权平均值。您可以通过设置带有焦点的权重矩阵来计算
# compute weights matrix
a <- sort(adjacent(r, 25, 8, pairs=F, include=TRUE))
axy <- xyFromCell(r, a)
d <- pointDistance(axy, xyFromCell(r, 25), lonlat=F)
w <- matrix(d, 3, 3)
w[2,2] <- 0
w <- w / sum(w)
# A simpler approach could be:
# w <- matrix(c(0,.25,0,.25,0,.25,0,.25,0), 3, 3)
foc <- focal(r, w, na.rm=TRUE, NAonly=TRUE)
foc[25]
在这个例子中，这很好；但如果在焦点区域有多个 NA 值（因为权重之和将不再是 1），这将是不正确的。我们可以通过计算权重之和来纠正
x <- as.integer(r/r)
sum_weights <- focal(x, w, na.rm=TRUE, NAonly=TRUE)
fw <- foc/sum_weights
done <- cover(r, fw)
done[25]