Prediction using spatial layers (rasters) with brms.fit

predict.brmsfit expects a data.frame. Would it be possible to use a brmsfit object to predict using a raster stack?
If so, can you share a guide?

Can you elaborate a bit more?

Thanks for the quick reply Paul (as usual)!

In my workflow I extract values from a raster stack and then use those values as the independent variables in a mixed effects model.
When attempting to use the brmsfit object to predict values using the same raster stack I get an error stating that the input values need to be a dataframe.
If I run a simple lm, I can use it to predict a raster stack.

Here’s some rough code to demonstrate my workflow:

datXY <- SpatialPoints(datXY, proj4string=CRS('+proj=longlat +datum=WGS84'))
xVar <- rnorm(length(datXY))
filesraster <- list.files('rasters', full.names = T)
layers <- stack(filesraster)
dat <- data.frame(extract(layers, datXY))
dat <- cbind(dat, xVar)
mod <- brm(xVar ~ layer1 + layer2 + (1 | layer3), data = dat, prior = c(set_prior("normal(0.1, 1)", class = "b")), family = gaussian(), chains = 3, iter = 8000, warmup = 1000)
p <- predict(mod, layers, allow.new.levels = F)

For the brm call, you seemed to have transformed the stack into a data.frame. Wouldn’t it be possible to do the same for the newdata in predict()?

I did try converting to a DF, but I reach the max memory allocation. I can’t export directly as CSV and then import it with a clean workspace either.
I’m already allocating 16 gb of ram to R …

I see. In any case, your transformed DF should contain the variable layer1, layer2 and layer3 or otherwise brms won’t be able to make sense of it.

Anyway, a conversion to DF is necessary because brms has no idea how to deal with rasters. It is not designed to know. The newdata should have the same structure as the original data.

Ok. Good to know.
The model is run on values extracted from the rasters at specific points.

I’ll work to figure it out and then will post if I get a workaround.

Curious if you ever found a solution to this issue, @sfrav? The raster package says “Any type of model (e.g. glm, gam, randomForest) for which a predict method has been implemented (or can be implemented) can be used”; however, that does not appear to be the case for brms’ predict function. In my case, I fit a simple brm multilevel model with a quadratic effect of elevation and latitude on animal population density, and am trying to predict the density values for existing elevation and latitude values across the entire study area as a rasterized map (e.g., by adapting the lm/glm methods: HKU Movement Eco in R Workshop - The homepage for the Movement Ecology in R Workshop hosted at Hong Kong University, January 3-12, 2018).

Continuing the discussion from Prediction using spatial layers (rasters) with brms.fit:

Sounds like an interesting analysis. What kind of scale are you working on? If local to regional, you might find the R-INLA package useful. See here for tutorials. I’m not sure how well this would scale globally though.

1 Like

It’s a pseudo-global-scale analysis (North and South America). The species’ extant range spans from about 50N to 60S, and we have density estimates from throughout that latitude gradient. I’ll checkout R-INLA. Thanks for the tip and links to tutorials.

Someone was using Stan for kriging here Geostatistical modelling with R and Stan - The Academic Health Economists' Blog that might be a starting point.

1 Like

raster::as.data.frame will convert a raster * object to a data.frame. However, the raster package allows R to interact with rasters that are stored on disk without loading them into RAM. Therefore, converting to data.frame can exhaust available memory if the underlying files are too large. To circumvent this issue, it should work as expected to first crop the raster into chunks small enough to fit in RAM, and then to predict to these chunks one-at-a-time.

2 Likes