[R-sig-ME] lmer fails when too many observations
Asaf Weinstein
asafw.at.wharton at gmail.com
Wed Mar 11 04:24:47 CET 2015
Dear Thierry and Douglas,
Thank you very much for your responses; It seems like, as Douglas remarked,
the error message I got was not a result of too many observations; I tried
a similar example on a different computer, and it worked fine (as a matter
of fact, I tried it on the same computer again and I it also worked, only
gave a warning message now rather than an error).
I regret sending the original email and sorry for the false alarm!
Here is the example, in any case:
## generate data from a two-way random-effects model with balanced design
##
R <- 100
C <- 100
sig <- 5 # sd, not variance
mu <- 0
theta.A <- 1/sqrt(4*C) # variance order of sigma^2/(# of col's) for
shrinkage to take place
theta.B <- 1/sqrt(4*R)
K <- matrix( 10, nrow=R, ncol=C )
level.A.all <- as.factor( rep(t(row(K)),t(K)) ) #same as factor( rep(1:R,
apply(K,1,sum)) )
level.B.all <- as.factor( rep( t(col(K)),t(K) ) ) #same as factor( unlist(
apply(K,1,function(t) rep( 1:C, t )) ) )
#cbind(level.A.all, level.B.all)
#table(level.A.all,level.B.all) # same as K
effects.A <- rnorm(R, mean=0, sd= sig * theta.A)
effects.B <- rnorm(C, mean=0, sd= sig * theta.B)
eta.all <- mu + effects.A[level.A.all] + effects.B[level.B.all] # cell
means for individual obs model
eps.all <- rnorm( sum(K), mean= 0, sd = sig )
y.all <- eta.all + eps.all
#head(cbind(level.A.all,level.B.all,eta.all, y.all))
## estimation ##
# LS estimates
fm.LS <- lm(y.all~1+level.A.all+level.B.all) # Least Squares
etahat.LS <- c( t( tapply(fitted(fm.LS),list(level.A.all,level.B.all),mean)
) )
# EB Maximum-Likelihood
fm.EBML <- blmer( y.all ~ 1+(1|level.A.all)+(1|level.B.all), cov.prior =
NULL, resid.prior = point(value = sig, posterior.scale='sd'), REML=FALSE )
theta.EBML <- getME(fm.EBML,'theta')
etahat.EBML <- c( t( tapply(fitted(fm.EBML),list(level.A.all,level.B.all),
mean) ) )
Thank you again,
Asaf
On 10 March 2015 at 12:40, Douglas Bates <bates at stat.wisc.edu> wrote:
> It is not clear that the error message you quote is a result of too many
> observations.
>
> As Thierry stated, it would help to have a small reproducible example.
>
> On Tue, Mar 10, 2015 at 9:24 AM Asaf Weinstein <asafw.at.wharton at gmail.com>
> wrote:
>
>> Dear lmer community,
>>
>> I am trying to run a simulation for a two-way random-effects model with
>> unbalanced design (ie, unequal number of observations per cell) and no
>> interaction.
>> It's especially important for me to be able to run the lmer/blmer
>> functions
>> when the number of (column and row) random effects is large, say 100, and
>> with possible replicates in each cell.
>> The problem is that lmer() works with the full vector of observations, as
>> opposed to working with the cell averages (which is a sufficient
>> statistic), and the methods fails pretty quickly when there are replicates
>> (because the response vector is too big, I suppose). I get the following
>> error:
>>
>> *Error in get("checkConv", lme4Env)(attr(opt, "derivs"), opt$par, ctrl =
>> control$checkConv, : *
>> * (converted from warning) Model failed to converge with max|grad| =
>> 0.00244385 (tol = 0.002)*
>>
>> Just to give an example: suppose there are R=100 row effects, C=100 column
>> effects, and 5 replicates in each cell. The vector of individual
>> observations is of length 100^5 (lmer fails), while the vector of cell
>> averages is of length 100^2 (a size which causes no problem for lmer).
>> My question is whether there is a way to tell lmer() to work with the
>> sufficient statistic (of course, the conditional covariance is no longer
>> c*Identity, a fact which is used in the implementation of lmer (according
>> to documentation) ).
>>
>> Thank you very much and I hope I was clear!
>>
>> Asaf
>>
>> [[alternative HTML version deleted]]
>>
>> _______________________________________________
>> R-sig-mixed-models at r-project.org mailing list
>> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>>
>
[[alternative HTML version deleted]]
More information about the R-sig-mixed-models
mailing list