# How to do regression analysis out-of-memory on a set of large rasters in R?

I am trying to do a regression analysis on a set of large rasters (254,004,000 cells each). Ultimately, I want to run something like the following (or a bit more complex, but let's start simple!):

model<-lm(dv ~ iv1+iv2+iv3... data=df,na.action=na.exclude)

where "dv" is the values from one raster and "iv1", "iv2" "iv3" ... are values from other rasters (up to 10 variables) with the same extent and resolution. It seems I should be able to do this out-of-memory using the Raster package, but I am confused how. Whether I create a brick, stack, or set of individual Raster objects, I cannot figure out how to send the variables to the lm function without using getValues and thus calling everything into memory (mine cannot even handle two variables).

A point in the right direction would be much appreciated!

• I would take a moment and consider this in statistical terms. 1) you are effectively using the population and not a sample thus, negating the need for a regression. 2) using all the cells in a rasters is going to certainty add an unnecessary autocorrelation issue to a linear model. 3) In classical statistical terms, you will have a psuedoreplication (lack of independence) issue. 4) I highly doubt that you would meet iid assumptions. I would recommend taking a sample of the raster(s), use the sample data to build your regression model then estimate the model to your raster(s). – Jeffrey Evans Sep 26 '13 at 23:26
• Thanks, Jeffery, I appreciate the note. This would not be my final statistical product, but in conjunction with autocorrelation plots for each variable, I find it helps me with diagnostics. The answer below seems like it might be a fruitful path. – Kendra Walker Sep 27 '13 at 13:30
• I beg to differ (slightly) with some of @Jeffrey Evans' points. First, regression for an entire population is meaningful: it describes relationships among variables. Second, autocorrelation is not necessarily a problem, but the advice to worry about it is excellent. There is a direct solution: tile your rasters. For each tile compute the mean, the count, and the [SSP matrix]. You can combine these statistics and proceed with the solution. There's no limit to the raster size this applies to. Another approach (using 2 rasters at a time) is given at stats.stackexchange.com/a/71257. – whuber Sep 27 '13 at 15:56
• I should be more specific. I do believe that regression approaches on rasters are useful in the context of "exploratory" analysis. One thought, have you considered an OLS rather than a straight linear model? The resulting residual error in OLS is a bit more robust to autocorrelation issues. – Jeffrey Evans Sep 27 '13 at 19:22
• Sorry, I should have been more specific as well. I do not necessarily want to pin myself to a linear model, I just thought that if I could get something to run with lm, I could carry it over to other similar packages (certainly not the most direct approach; a direct route to an OLS or other more robust method would be very welcome!) I have seen several examples of lm run in memory with the Raster package and get the impression that it can manage problems like these out of memory as well using a brick/stack object – Kendra Walker Sep 27 '13 at 20:18

The help for `lm` references `biglm`:
`biglm` in package biglm for an alternative way to fit linear models to large datasets.
The help pages for `biglm` indicate this package was developed for precisely such problems. The algorithm it references, AS274, is an updating procedure, allowing a solution based on a subset of the cases (cells) to be modified as additional cases are given.