无法有趣地计算 R 中光栅堆栈上的像素级回归

我正在处理栅格，我有一个带有 7n 层的 RasterStack。我想使用下面的公式计算像素级回归。我试图这样做raster::calc，但我的功能失败并显示消息：

' lm.fit(x, y, 偏移量 = 偏移量, 奇异.ok = 奇异.ok, 中的错误， ...) : 0 (非 NA) 案例。'

但所有栅格都可以，并且包含数字（不仅是 NA），我可以绘制它， 我可以用公式计算一般线性回归

 cr.sig=lm (raster::as.array(MK_trend.EVI.sig_Only) ~ raster::as.array(stack.pet)+raster::as.array(stack.tmp)+raster::as.array(stack.vap)+raster::as.array(stack.pre)+raster::as.array(stack.wet)+raster::as.array(stack.dtr))

但是当我堆叠图层时

allData = stack(MK_trend.EVI.sig_Only,stack.dtr,stack.wet,stack.pre,stack.vap,stack.tmp,stack.pet)

并尝试 calc 函数

    # Regression Function, R2
lmFun=function(x){
    x1=as.vector(x);
    if (is.na(x1[1])){
        NA 
    } else {
        m = lm(x1[1] ~ x1[2]+x1[3]+x1[4]+x1[5]+x1[6]+x1[7])
        return(summary(m)$r.squared)
    }
}

我看到错误消息。
我对 R 和编程还很陌生，所以，也许有一些愚蠢的错误？ 为了使处理工作顺利进行，我将不胜感激。

您可以使用calc对于像素级（局部）回归，但您的公式似乎表明您想要其他东西（全局模型）。

如果回归是按像素进行的，则每个单元格的 x 和 y 值数量相同，您可以使用calc. See ?calc举些例子。

相反，每个单元格有 1 个 y（独立）和 6 个 x（因）变量值。这表明您想要一个全球模型。为此，您可以执行以下操作：

library(raster)
# example data
r <- raster(nrow=10, ncol=10)
set.seed(0) 
s <- stack(lapply(1:7, function(i) setValues(r, rnorm(ncell(r), i, 3))))  
x <- values(s)

# model
m <- lm(layer.1 ~ ., data=data.frame(x))

# prediction
p <- predict(s, m)

这需要将所有值加载到内存中。如果你做不到这一点，你可以采取大量的常规样本。看sampleRegular

并说明为什么你的方法不起作用：

testFun=function(x1){
    lm(x1[1] ~ x1[2]+x1[3]+x1[4]+x1[5]+x1[6]+x1[7])
}

# first cell
v <- s[1]
v
#      layer.1  layer.2   layer.3  layer.4  layer.5  layer.6  layer.7
#[1,] 4.788863 4.345578 -0.137153 3.626695 3.829971 4.120895 1.936597

m <- testFun(v)
m
#Call:
#lm(formula = x1[1] ~ x1[2] + x1[3] + x1[4] + x1[5] + x1[6] + x1[7])

#Coefficients:
#(Intercept)        x1[2]        x1[3]        x1[4]        x1[5]        x1[6]        x1[7]  
#      4.789           NA           NA           NA           NA           NA           NA  


summary(m)$r.squared
# 0

即使我没有收到您报告的错误消息（但所有 R^2 值均为零）。