[R-sig-ME] Fwd: lmer2 slower?

[Douglas Bates]

I am forwarding (with permission) a conversation with Bill Auty who
discovered that in one of his examples lmer2 is considerably slower
than lmer, although they get to the same estimates.  This can happen.

As I mentioned in one reply I hope that forcing a supernodal Cholesky
decomposition for models with crossed or partially crossed grouping
factors for the random effects will help but I haven't had time to
explore this yet.

I'll post some other timing results in another message.

---------- Forwarded message ----------
From: Douglas Bates <bates at stat.wisc.edu>
Date: Jan 26, 2007 3:38 PM
Subject: Re: lmer2 slower?
To: Bill Auty <bill at edmeasure.com>


On 1/26/07, Bill Auty <bill at edmeasure.com> wrote:
>  > class(xq2 at L)
>
> [1] "dCHMsimpl"
>
> attr(,"package")
>
> [1] "Matrix"

Ah, good.  That 's the answer that I wanted.

Briefly, the CHOLMOD sparse matrix code allows for two types of
Cholesky decompositions of sparse, symmetric matrices.  These are
called "simplicial" and "supernodal".  In the lmer code I forced the
supernodal decomposition.  In the lmer2 code I allow the sparse matrix
manipulation code to choose one or the other based on some criteria
which I can "tune".  I think that the default values that are
currently being used are too conservative about switching to the
supernodal decomposition so I may be able to gain back some of the
loss of speed.

Thanks again for checking this.

May I send a copy of this reply to the R-SIG-mixed-models list?

>
>  >
>
> I don't have the math chops to stay with you on the theory of this, but
> I'm glad to help by trying things on my data.
>
> Douglas Bates wrote:
> > Thanks.  That's valuable information.
> >
> > Could you check something for me?  If you still have the object fit by
> > lmer2 could you send me the result of
> >
> > class(xq2 at L)
> >
> > The response should be either "dCHMsimpl" or "dCHMsuper"
> >
> > On 1/26/07, Bill Auty <bill at edmeasure.com> wrote:
> >> I've attached the output from running lmer and lmer2 on the same model.
> >> The variables are:
> >>
> >> rit = scale score on a math test
> >> gr = grade (range 0 - 2)
> >> stu = student ID
> >> sch = school ID
> >> dist = district ID
> >>
> >> The students are partially crossed in school and district.
> >>
> >>
> >> >
> >> system.time(xq2<-lmer2(rit~gr+I(gr^2)+(gr|stu)+(gr|sch)+(gr|dist),e2math,
> >>
> >> + control = list(niterEM =0, gradient = FALSE),method="ML"))
> >> [1] 428.223  30.569 458.986   0.000   0.000
> >> > xq2
> >> Linear mixed-effects model fit by maximum likelihood
> >> Formula: rit ~ gr + I(gr^2) + (gr | stu) + (gr | sch) + (gr | dist)
> >>    Data: e2math
> >>     AIC    BIC  logLik MLdeviance REMLdeviance
> >>  977133 977251 -488554     977109       977117
> >> Random effects:
> >>  Groups   Name        Variance Std.Dev. Corr
> >>  stu      (Intercept) 78.60373 8.86587
> >>           gr           1.88712 1.37373  -0.227
> >>  sch      (Intercept)  7.24479 2.69161
> >>           gr           0.75162 0.86696  -0.512
> >>  dist     (Intercept)  4.30657 2.07523
> >>           gr           0.37241 0.61026  0.139
> >>  Residual             26.61534 5.15901
> >> Number of obs: 136484, groups: stu, 80224; sch, 182; dist, 63
> >>
> >> Fixed effects:
> >>              Estimate Std. Error t value
> >> (Intercept) 209.65220    0.38078   550.6
> >> gr            7.67870    0.13085    58.7
> >> I(gr^2)      -0.76376    0.01931   -39.6
> >>
> >> Correlation of Fixed Effects:
> >>         (Intr) gr
> >> gr      -0.188
> >> I(gr^2)  0.049 -0.411
> >> >
> >> system.time(q2<-lmer(rit~gr+I(gr^2)+(gr|stu)+(gr|sch)+(gr|dist),e2math,
> >> + control = list(niterEM =0, gradient = FALSE),method="ML"))
> >> [1] 241.323  24.070 265.412   0.000   0.000
> >> > q2
> >> Linear mixed-effects model fit by maximum likelihood
> >> Formula: rit ~ gr + I(gr^2) + (gr | stu) + (gr | sch) + (gr | dist)
> >>    Data: e2math
> >>     AIC    BIC  logLik MLdeviance REMLdeviance
> >>  977133 977251 -488554     977109       977117
> >> Random effects:
> >>  Groups   Name        Variance Std.Dev. Corr
> >>  stu      (Intercept) 78.60255 8.86581
> >>           gr           1.88767 1.37393  -0.227
> >>  sch      (Intercept)  7.24380 2.69143
> >>           gr           0.75110 0.86666  -0.512
> >>  dist     (Intercept)  4.30710 2.07536
> >>           gr           0.37265 0.61045  0.139
> >>  Residual             26.61352 5.15883
> >> number of obs: 136484, groups: stu, 80224; sch, 182; dist, 63
> >>
> >> Fixed effects:
> >>              Estimate Std. Error t value
> >> (Intercept) 209.65218    0.38079   550.6
> >> gr            7.67867    0.13086    58.7
> >> I(gr^2)      -0.76376    0.01931   -39.6
> >>
> >> Correlation of Fixed Effects:
> >>         (Intr) gr
> >> gr      -0.188
> >> I(gr^2)  0.049 -0.411
> >> >
> >>
> >>
> >>
> >>
> >
>