在R中插值栅格数据的时间序列的高效方法

Question

在R中插值栅格数据的时间序列的高效方法

3

我有一些大型栅格数据，代表时间序列（每个图层是一个时间点）。我希望通过在层之间进行插值来增加时间分辨率（在下面的示例中，使用zoo::na.spline进行插值，尽管如果有助于解决问题，我也可以使用其他插值方法）。我不能直接在栅格数据上使用na.spline - 它不是为此设计的，并且只会给出Error in as.array.default(X) : 'dimnames' applied to non-array。

目前，我正在使用一系列循环来完成此操作（第1个循环用于生成具有更多时间层的栅格数据，第2个循环用于将原始栅格数据插入到其中的正确层中，第3个循环用于对每个单元格使用na.spline）。然而，在大型栅格数据上这种方法非常缓慢，似乎是一种相当低效的方法。

一个最小的可重现示例。

首先让我们制作一个表示低时间分辨率数据的初始栅格数据集。请注意，左上角的单元格始终为NA，因为任何解决方案都必须对仅包含NA的单元格进行鲁棒性检查：

library(raster); library(rasterVis); library(zoo)
r1 = raster(matrix(c(NA, rep(1,8)), 3, 3))
r2 = raster(matrix(c(NA, rep(2,8)), 3, 3))
r3 = raster(matrix(c(NA, rep(3,8)), 3, 3))    
b1 = brick(list(r1,r2,r3))
levelplot(b1)

现在，让我们创建一个空的栅格砖块，以便插入值：

b2 = brick(lapply(1:9, function(x) raster(matrix(NA, 3, 3))))

接下来，我们将b1的图层插入到b2的相应图层中。请注意，我甚至在这里使用了一个循环，因为将选择的图层子分配到光栅砖中不起作用。

old.layers = c(1,4,7)
for (i in 1:nlayers(b1)) {
  b2[[old.layers[i]]] = b1[[i]]
} 
levelplot(b2, layout=c(9,1))

最后，我们循环遍历每个单元格进行时间序列插值，并将结果插入回b2中：

for (cell in 1:ncell(b2)) {
  if (all(is.na(as.vector(b2[cell])))) {
    # na.spline wil fail if all are NA
    new.values = as.numeric(rep(NA, dim(b2)[3])) 
  } else {
    new.values = na.spline(as.vector(b2[cell]), na.rm = F) 
  }
  b2[cell][] = new.values
}
levelplot(b2, layout=c(9,1))

这个方法可以实现，但似乎非常低效和不优雅。有没有更快（最好也更优雅）的方法来做到这一点？

- dww

1个回答

网页内容由stack overflow 提供, 点击上面的

可以查看英文原文，
原文链接

- Robert Hijmans · Accepted Answer

你可以使用像这样的函数与raster::calc一起使用：

yy <- rep(NA, 9)
fun <- function(y) { 
  if (all(is.na(y)))  yy else na.spline((yy[c(1,4,7)] <- y), xout=1:9)
}
b2 <- calc(b1,  fun)

供参考，raster::approxNA 接近于您所需的功能（但它不会添加图层，并且使用线性插值）。