2
我正在使用R來執行我的迴歸。我成功地將r平方剩餘標準錯誤提取到了rasterbrick中。然後我需要運行另一個代碼來獲得p值(F stat)。我怎麼能把fun1和fun2結合起來,這樣我就可以製作一個包含這些信息的rasterbrick?如何將回歸摘要(例如p值和coeff)輸出到rasterbrick中?
這裏是我的代碼:
library(raster)
#1 create test data
r <- raster(nrow=10, ncol=10)
set.seed(0)
s <- stack(lapply(1:12, function(i) setValues(r, rnorm(ncell(r), i, 3))))
time <- 1:nlayers(s)
s[1:5] <- NA
#2 Run function1 to obtain r-squared and residual standard error
fun1 <- function(x) {
if (all(is.na(x))) {
return(cbind(NA,NA))
}
m = lm(x~time)
s <- summary(m)
r2 <- s$r.squared
resid.s.e <- s$sigma
cbind(r2, resid.s.e)
}
#3 Run function to obtaion p-value(from F stat)
fun2 <- function(x) {
if (all(is.na(x))) {
return(cbind(NA,NA))
}
m = lm(x~time)
s <- summary(m)
r2 <- s$r.squared
pf<- pf(s$fstatistic[1], s$fstatistic[2], s$fstatistic[3],lower.tail = FALSE)
cbind(r2, pf)
}
#Apply both functions with rasterstack and plot
r <- calc(s, fun)
plot(r)
r2 <- calc(s, fun2)
plot(r2)
在此先感謝。
你的文字說p值,你的標題說係數。這是什麼? –
謝謝@羅曼·盧斯特里克。編輯我的標題和帖子。我很匆忙,混了起來。道歉。 – Eddie
爲什麼不讓1個函數返回'resid.s.e'和'pf',然後選擇變量傳遞給'calc'?重新計算相同的'lm'似乎是多餘的。 – zx8754