通过calc函数在R上从Lineair回归中提取p值在R中

问题描述

对于几年前回答的问题,我基本上有一个后续问题,涉及在栅格堆栈上进行线性回归。 (请参见UPDATE statement (Microsoft Access SQL)) 我进行了线性回归,并使用$ coefficients中的值计算了趋势,但现在想知道相关的p值(每个栅格像素一个)。

但是,calc抱怨无法为签名“” integer”,“ function””找到函数“ calc”的继承方法

您可以使用以下代码重现此错误

library(raster)
names = c('...','...','...')
s <- stack(names)
y <- values(s)
x <- log(c(10,20,30,40,50))    
funa <- function(y) { 
        if(all(is.na(y))) {
            c(NA,NA)
        } else {
            summary(lm(y ~ x))$coefficients
        }
    }
    r <- calc(s,funa)

我可以理解calc不知道如何将摘要输出转换为新的栅格堆栈。因此,我尝试将lm的输出重塑为其他形式,例如使用“ broom”包,但无济于事。现在calc抱怨:is.infinite(v)中的错误:即使试图将输出强制为data.frame或as.numeric,也无法为“列表”类型实现方法。例如这样

library(broom)
funlm <- function(y) {
    if(all(is.na(y))) {
        c(NA,NA)
    } else {
        as.data.frame(glance(lm(y ~ x)))
    }
}
r <- calc(s,funlm)

请帮助

解决方法

这里是一个示例,基于?calc

中的示例
# create data
library(raster)
r <- raster(nrow=10,ncol=10)
s1 <- lapply(1:12,function(i) setValues(r,rnorm(ncell(r),i,3)))
s2 <- lapply(1:12,3)))
s1 <- stack(s1)
s2 <- stack(s2)

# regression of values in one brick (or stack) with another,extract p-values
s <- stack(s1,s2)
fun <- function(x) {
    if (all(is.na(x))) {
       return(c(NA,NA))
    } 
    m <- lm(x[1:12] ~ x[13:24])
    summary(m)$coefficients[,4]
}
x1 <- calc(s,fun)