Linear regression between one global data(raster) and one variable(a fixed value for all grid cells)

Question

I have time series of a global dataset(raster) and a variable(a fixed value for all grid cells) related to this dataset for the corresponding year.

Global dataset example(not sure why cannot assign the order as a value in this example):

library(terra)
# Set dimensions of the raster
nrow <- 360
ncol <- 720

# Set number of rasters to create
n <- 4

r_final <- rast(nrows=nrow, ncols=ncol)
# Create a loop to generate each raster
for(i in 1:n) {
  # Create empty raster with specified dimensions
  r <- rast(nrow=nrow, ncol=ncol, nlyr=1,vals=i)
  
  # Set the extent and resolution of the raster
  ext(r) <- c(-180, 180, -90, 90)
  res(r) <- 0.5
 
  raster_final <- c(raster_final,r)
}
raster_final = raster_final[[2:nlyr(r_final)]]

variable example

variable =data.table(Year=c(1:4),value=0:3) 

I want to bulid a linear regression between each grid cell value of the raster and value in variable [something like lm(raster_final ~ variable$value, data= ?)] and create a raster of the R2 global map.

How could I do that?

Answer 1

Example data

library(terra)
s <- rast(system.file("ex/logo.tif", package="terra"))   
time <- 1:nlyr(s)

Using lm, here only extracting the slope

fun <- function(x) { lm(x ~ time)$coefficients[2] }
x <- app(s, fun)

A much faster approach with linear algebra

## add 1 for a model with an intercept
X <- cbind(1, time)

## pre-compute constant part of least squares
invXtX <- solve(t(X) %*% X) %*% t(X)

## regression model; [2] is to get the slope
quickfun <- function(y) (invXtX %*% y)[2]

y <- app(s, quickfun) 

I do not know what your R code is for, but to create four layers where all cell values correspond to the layer number:

r <- rast(nrow=360, ncol=720, res=0.5, nlyr=4, vals=rbind(1:4))