Raster linear and conditional regression using raster stacks by month in R

Question

I have two raster stacks and I want to carry out a refression analysis. If each raster in each stack was a month in the year (6 data points would be three months in two years i.e. January, February and March for two different years), how do I calculate the slope using the indices such that the result generates 3 slope rasters (one for each month) please?

#First raster track
r <- raster(ncol=10, nrow=10)
r[]=1:ncell(r)
S <- stack(r,r,r,r,r,r)

#Second raster stack
r1 <- raster(ncol=10, nrow=10)
r1[]=1:ncell(r1)
N <- stack(r1,r1,r1,r1,r1,r1)

#combine both raster stacks
s <- stack(S,N)

#function to calculate slope
fun=function(x) { if (is.na(x[1])){ NA } else { lm(x[7:12] ~ x[1:6] )$coefficients [2]}}

#apply function
slope <- calc(s, fun)

Result should be 3 rasters.

A second question: If I wanted to do a conditional regression using a third raster stack, what would the codes be?

Answer 1

Try fun with 1:12

fun(1:12)

# Error in model.frame.default(formula = x[6:12] ~ x[1:6], drop.unused.levels = TRUE) : 
#  variable lengths differ (found for 'x[1:6]')

it should be

fun=function(x) { if (is.na(x[1])){ NA } else { lm(x[7:12] ~ x[1:6] )$coefficients [2]}}

Working example

library(raster)
r <- raster(ncol=10, nrow=10)
set.seed(99)
s <- stack(sapply(1:12, function(i) setValues(r, runif(ncell(r)))))

fun <- function(x) { if (is.na(x[1])){ NA } else { lm(x[7:12] ~ x[1:6] )$coefficients [2]}}

slope <- calc(s, fun)

For the three slopes:

fun3 <- function(x) {
    r <- rep(NA, 3)
    if (!is.na(x[1])) {
       r[1] <- lm(x[3:4] ~ x[1:2] )$coefficients[2]
       r[2] <- lm(x[7:8] ~ x[5:6] )$coefficients[2]
       r[3] <- lm(x[11:12] ~ x[9:10] )$coefficients[2]
    }
    r
 }

slope3 <- calc(s, fun3)