large dataset - lmer2 vector size specified is too large
On 3/2/07, Gregor Gorjanc <gregor.gorjanc at bfro.uni-lj.si> wrote:
Douglas Bates <bates at ...> writes: ...
but even that formula means that you are estimating 448 fixed effects for which the model matrix is of size 448 * 43680 * 8 bytes (about 150 MB). In addition you are attempting to estimate
Not directly related. Is it really necessary to build model matrix i.e. X.
Necessary? For linear mixed models no - for generalized linear mixed models yes (think of how the IRLS algorithm works). Easy to change? Regretably, no.
Building X'X directly would be more efficient i.e. only 448 * 448 * 8 bytes.
That wouldn't be sufficient. As described in the "Implementation" vignette in the lme4 package, the calculations in lmer2 are based on the matrix A which is [Z X -y]'[Z X -y]. You can't use just X'X - you need the crossproducts of all the columns of Z, X and y. The matrix A is potentially huge but sparse. Assuming that Florian really does want to try to estimate a model with (affyID|cephID) then Z has 224 * 195 = 43680 columns. In lmer Z'Z is stored as a sparse matrix. In lmer2 A is stored as a sparse matrix. Generating A as a sparse matrix directly from the data or directly from a model.frame object would result in considerable savings in memory usage but doing so is far from trivial. Most useRs are not aware (nor need they be aware) of the complexity of the operations involved in taking a linear model formula and a data specification (and all the arguments like subset, na.action, offset, weights, ... that are expected to be available in model-fitting functions) and creating first a model.frame then a model.matrix. It is one of those things that "just works" but is incredibly complicated beneath the surface. Even the first stage (producing a model frame) is complicated - to the extent that Thomas Lumley has a report on "The Standard Non-standard Evaluation in R" dealing with how arguments are evaluated in the model.frame function. There are various approaches that could be taken to building A as a sparse matrix directly from the model specification and a model frame. Probably the easiest to implement would be to use the existing model.matrix code to produce horizontal sections of [Z X -y] and accumulate A from these sections. This is still far from trivial because you either need two passes (one to determine the positions of the nonzeros and one to evaluate them) or you have to somehow update the structure of A on the fly. The other thing that gets tricky is that the columns of X may change from section to section because not all levels of a factor may be used in a given section. It would be necessary to ensure that the columns from all sections are aligned.