2

Is there any way to estimate random effect panel linear regression in R using weights option?

For fixed effects, I can simply use

lm(y ~ x + factors(id_variable), weights=wgt)

However I am looking for RE as the number of individuals is too large and I hit memory limits.

Weighting is needed due to survey data type

If not available in R, is it available in stata or SAS?

ewe
  • 21
  • 2

1 Answers1

2

In my (limited) experience survey weighting for multilevel models is a can of worms (I mostly follow Andrew Gelman's opinion on this, e.g this paper or this blog post about it; here's an (unanswered) CrossValidated question on the topic). In other words, you might have to solve a statistical problem (what is the right way to handle this?) before you worry about the computational problem.

I do think Stata and/or SAS offer such capabilities (see CrossValidated link above), but I can't tell you much about them.

However, if you don't want/need random effects, i.e. you just want to handle a problem with a large number of fixed effects, you should be able to use a sparse model matrix as in this question, as in the example below.

It will still be up to you to figure out how to deal with the weights. I don't believe that the weights= argument in lm, glm, lme4, etc. handle survey weights. I believe the survey package is very sophisticated, but I don't know if its surveyglm() function can handle/has an option for sparse model matrices ...

library("Matrix")
library("MatrixModels")
set.seed(101)
ngrps <- 1e3
nobs <- 10
ntot <- nobs*ngrps
d <- data.frame(y=rnorm(ntot),
                x=rnorm(ntot),
                f=gl(ngrps,nobs))
object.size(X1 <- model.matrix(~x+f,data=d))  ## 80 MB
object.size(X2 <- sparse.model.matrix(~x+f,data=d))  ## 709K
system.time(m1 <- MatrixModels:::lm.fit.sparse(X2,d$y)) ## 0.004 seconds
system.time(m2 <- lm.fit(X1,d$y))  ## 9 seconds
Community
  • 1
  • 1
Ben Bolker
  • 211,554
  • 25
  • 370
  • 453