移动 window 回归

Question

我想对代表 Landsat 数据 Band3 和 Band4 的两个栅格堆栈的每个像素执行移动 window 回归。结果应该是两个额外的堆栈，一个代表截距，另一个代表回归的斜率。因此堆栈的第 1 层 "B3" 和堆栈 "B4" 导致堆栈的第 1 层 "intercept" 和堆栈 "slope"。堆栈 B3 的第 2 层和堆栈 B4 导致第 2 层，...等等。

我已经有了 gwr 功能，但想留在光栅包中。我不知何故知道必须包含 focal 才能设置我的移动 window（应该是 3x3 像素）和某种线性模型，例如：lm(as.matrix(b3)~as.matrix(b4)) 虽然我不认为这获取像素值...

也可以采用逐层方法代替栅格堆栈。（所以不一定是Band3的rasterstack。

有人知道如何在 R 中对此进行编程吗？

Answer 1

这是一种方法，改编自 ?raster::localFun

set.seed(0)
b <- stack(system.file("external/rlogo.grd", package="raster"))
x <- flip(b[[2]], 'y') + runif(ncell(b))
y <- b[[1]] + runif(ncell(b))

# local regression:
rfun <- function(x, y, ...) {
    d <- na.omit(data.frame(x, y))
    if (nrow(d) < 3) return(NA)
    m <- lm(y~x, data=d)
    # return slope
    coefficients(m)[2]
}

ff <- localFun(x, y, fun=rfun)
plot(ff)

不幸的是，您必须运行这两次才能获得斜率和截距 (coefficients(m)[1])。

移动 window 回归

Moving window regression

r

raster

linear-regression