Focal regression 3x3 group of cells

I have two raster with same resolution and dimension: water-occurrence index and rainfall (mm). I want to create a new raster of SLOPE (a) and INTERCEPT (b) by performing linear regression between every 3×3 pixels of the two rasters (rainfall and water-occurrence), such that each pixel of the SLOPE and INTERCEPT will hold the regression slope and intercept value obtained from linear regression of the corresponding 3×3 pixels in water-occurrence and rainfall that surround that pixel.

I adapt the code provide by @jeffrey-evans in this thread https://gis.stackexchange.com/a/278985

# Linear regression between every 3×3 pixels of two rasters
# Based on thread from SE: https://gis.stackexchange.com/questions/278979/linear-regression-between-every-3×3-pixels-between-two-rasters-using-r

library(raster)
library(gstat)
library(sp)

# Set working directory
setwd("/Users/mac/data/")


# Raster import and stacking
# Change the path variable
fs <- list.files(path="/Users/mac/data", pattern = ".tif$", full.names = TRUE)
r <- raster::stack(fs)

# Focal linear regression
w <- c(3,3) # 3x3 neighbour/group of cells

# Slope
lapse_1 <- r[[1]]
lapse_1[] <- NA
for( rl in 1:nrow(r) ) { 
  v <- getValuesFocal(r[[1:2]], row=rl, nrows=1, ngb = w, array = FALSE)
  fit <- rep(NA,nrow(v[[1]]))
  for(i in 1:nrow(v[[1]]) ) {
    xy <- na.omit( data.frame(x=v[[1]][i,],y=v[[2]][i,]) )    
    if( nrow(xy) > 4 ) {
      fit[i] <- coefficients(lm(as.numeric(xy$y) ~ as.numeric(xy$x)))[2]
      if(is.null(fit)) fit = 1
    } else {
      fit[i] <- NA 
    }
  }
  lapse_1[rl,] <- fit
}

lapse_1
plot(lapse_1)

# Saving SLOPE into geotiff
writeRaster(lapse_1, filename=file.path("/Users/mac/output/Day1_01_Jan_Slope_3x3.tif"), format="GTiff", overwrite=TRUE)

# Intercept
lapse_2 <- r[[1]]
lapse_2[] <- NA
for( rl in 1:nrow(r) ) { 
  v <- getValuesFocal(r[[1:2]], row=rl, nrows=1, ngb = w, array = FALSE)
  fit <- rep(NA,nrow(v[[1]]))
  for(i in 1:nrow(v[[1]]) ) {
    xy <- na.omit( data.frame(x=v[[1]][i,],y=v[[2]][i,]) )    
    if( nrow(xy) > 4 ) {
      fit[i] <- coefficients(lm(as.numeric(xy$y) ~ as.numeric(xy$x)))[1]
      if(is.null(fit)) fit = 1
    } else {
      fit[i] <- NA 
    }
  }
  lapse_2[rl,] <- fit
}
lapse_2
plot(lapse_2)

# Saving INTERCEPT into geotiff
writeRaster(lapse_2, filename=file.path("/Users/mac/output/Day1_01_Jan_Intercept_3x3.tif"), format="GTiff", overwrite=TRUE)

I am not familiar with R and programming, I am wondering whether above script is correct or not for calculating slope and intercept?