[R-sig-ME] [R] different results from lme and lmer function

Pantelis Z. Hadjipantelis kalakouentin at gmail.com
Wed May 27 03:45:24 CEST 2015


Ben,

I fully accept that the two procedures should converge at the same 
solution but as the OP does not give the versions used. Therefore I am 
not sure that the fitting is done using the same procedure ("bobyqa" vs. 
"Nelder_Mead"). Different optim. algorithms with different initial 
values might converge to different local minima. The OP has only 39 
subjects so for a model of this size over-fitting will not be unheard 
off given the number of parameters.

I am curious though for your later comment: I understand that in the 
case of an unstable model you might expect lme4 to be slightly better, 
but wouldn't the singular fit in the correlation (-1) suggest that 
lmer's fit is sub-optimal?

All best,
Pantelis


On 05/26/2015 05:50 PM, Ben Bolker wrote:
> I agree that the difference is trivially small/practically
> unimportant.  The point here is that -- just for those of us who are
> interested in the details of the methods -- lme4 and lme are fitting
> *exactly* the same model, by similar methods, so in general they
> should converge to the same answer (to a somewhat closer tolerance
> than this).  Generally when they don't it's because the model is
> slightly unstable, and I have generally found that lme4 does slightly
> better (but I wouldn't rule out the opposite case).
>
>   cheers
>      Ben Bolker
>
>
> On Tue, May 26, 2015 at 8:47 PM, John Sorkin
> <jsorkin at grecc.umaryland.edu> wrote:
>> Ben,
>> I doubt the very small difference in log likelihood gives much, if any
>> information about which model is a better fit. Even if we overlook the
>> limited precision of the estimate of the REML criterion, the difference is
>> so small as to me of minimal importance.
>> John
>>
>> John David Sorkin M.D., Ph.D.
>>
>> Professor of Medicine
>>
>> Chief, Biostatistics and Informatics
>>
>> University of Maryland School of Medicine Division of Gerontology and
>> Geriatric Medicine
>>
>> Baltimore VA Medical Center
>>
>> 10 North Greene Street
>>
>> GRECC (BT/18/GR)
>>
>> Baltimore, MD 21201-1524
>>
>> (Phone) 410-605-7119
>>
>> (Fax) 410-605-7913 (Please call phone number above prior to faxing)
>>
>>
>> On May 26, 2015, at 8:03 PM, Ben Bolker <bbolker at gmail.com> wrote:
>>
>>   These actually aren't terribly different from each other.  I suspect
>> that lmer is slightly closer to the correct answer, because lme
>> reports a "log-likelihood" (really -1/2 times the REML criterion) of
>> 49.30021, while lmer reports a REML criterion of  -98.8 -> slightly
>> better fit at -R/2 = 49.4.  The residual sds are 0.0447 (lme) vs.
>> 0.0442 (lmer); the intercept sd estimate is 0.016 vs 0.0089,
>> admittedly a bit low, and both month sds are very small.  lmer
>> indicates a singular fit (correlation of -1).    If you look at the
>> confidence intervals on these estimates (confint(fitted_model) in
>> lme4; intervals(fitted_model) in lme) I think you'll find that the
>> confidence intervals are much wider than these differences (you may
>> even find that lme reports that it can't give you the intervals
>> because the Hessian [curvature] matrix is not positive definite).
>>
>>   Both should be comparable to SAS PROC MIXED results, I think, if
>> you get the syntax right ...
>>
>> On Tue, May 26, 2015 at 7:09 PM, li li <hannah.hlx at gmail.com> wrote:
>>
>> Hi all,
>>
>>   I am fitting a random slope and random intercept model using R. I
>>
>> used both lme and lmer funciton for the same model. However I got
>>
>> different results as shown below (different variance component
>>
>> estimates and so on). I think that is really confusing. They should
>>
>> produce close results. Anyone has any thoughts or suggestions. Also,
>>
>> which one should be comparable to sas results?
>>
>> Thanks!
>>
>>   Hanna
>>
>>
>> ## using lme function
>>
>> mod_lme <- lme(ti  ~ type*months, random=~ 1+months|lot, na.action=na.omit,
>>
>> + data=one, control = lmeControl(opt = "optim"))
>>
>> summary(mod_lme)
>>
>> Linear mixed-effects model fit by REML
>>
>> Data: one
>>
>>         AIC       BIC   logLik
>>
>>   -82.60042 -70.15763 49.30021
>>
>>
>> Random effects:
>>
>> Formula: ~1 + months | lot
>>
>> Structure: General positive-definite, Log-Cholesky parametrization
>>
>>             StdDev       Corr
>>
>> (Intercept) 8.907584e-03 (Intr)
>>
>> months      6.039781e-05 -0.096
>>
>> Residual    4.471243e-02
>>
>>
>> Fixed effects: ti ~ type * months
>>
>>                      Value   Std.Error DF   t-value p-value
>>
>> (Intercept)     0.25831245 0.016891587 31 15.292373  0.0000
>>
>> type            0.13502089 0.026676101  4  5.061493  0.0072
>>
>> months          0.00804790 0.001218941 31  6.602368  0.0000
>>
>> type:months -0.00693679 0.002981859 31 -2.326329  0.0267
>>
>> Correlation:
>>
>>                (Intr) typPPQ months
>>
>> type           -0.633
>>
>> months         -0.785  0.497
>>
>> type:months  0.321 -0.762 -0.409
>>
>>
>> Standardized Within-Group Residuals:
>>
>>           Min            Q1           Med            Q3           Max
>>
>> -2.162856e+00 -1.962972e-01 -2.771184e-05  3.749035e-01  2.088392e+00
>>
>>
>> Number of Observations: 39
>>
>> Number of Groups: 6
>>
>>
>>
>>
>>
>> ###Using lmer function
>>
>> mod_lmer <-lmer(ti  ~ type*months+(1+months|lot), na.action=na.omit,
>> data=one)
>>
>> summary(mod_lmer)
>>
>> Linear mixed model fit by REML t-tests use Satterthwaite approximations to
>>
>>   degrees of freedom [merModLmerTest]
>>
>> Formula: ti ~ type * months + (1 + months | lot)
>>
>>    Data: one
>>
>>
>> REML criterion at convergence: -98.8
>>
>>
>> Scaled residuals:
>>
>>     Min      1Q  Median      3Q     Max
>>
>> -2.1347 -0.2156 -0.0067  0.3615  2.0840
>>
>>
>> Random effects:
>>
>> Groups   Name        Variance  Std.Dev.  Corr
>>
>> lot      (Intercept) 2.870e-04 0.0169424
>>
>>           months      4.135e-07 0.0006431 -1.00
>>
>> Residual             1.950e-03 0.0441644
>>
>> Number of obs: 39, groups:  lot, 6
>>
>>
>> Fixed effects:
>>
>>                 Estimate Std. Error        df t value Pr(>|t|)
>>
>> (Intercept)     0.258312   0.018661  4.820000  13.842 4.59e-05 ***
>>
>> type         0.135021   0.028880  6.802000   4.675  0.00245 **
>>
>> months          0.008048   0.001259 11.943000   6.390 3.53e-05 ***
>>
>> type:months -0.006937   0.002991 28.910000  -2.319  0.02767 *
>>
>> ---
>>
>> Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
>>
>>
>> Correlation of Fixed Effects:
>>
>>             (Intr) typPPQ months
>>
>> type     -0.646
>>
>> months      -0.825  0.533
>>
>> type:month  0.347 -0.768 -0.421
>>
>>
>> _______________________________________________
>>
>> R-sig-mixed-models at r-project.org mailing list
>>
>> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>>
>>
>> ______________________________________________
>> R-help at r-project.org mailing list -- To UNSUBSCRIBE and more, see
>> https://stat.ethz.ch/mailman/listinfo/r-help
>> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
>> and provide commented, minimal, self-contained, reproducible code.
>>
>>
>> Confidentiality Statement:
>>
>> This email message, including any attachments, is for ...{{dropped:7}}
> _______________________________________________
> R-sig-mixed-models at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models



More information about the R-sig-mixed-models mailing list