Bookmark and Share

Notice: On April 23, 2014, Statalist moved from an email list to a forum, based at statalist.org.


[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]

st: Separate intercept mixed model


From   Michael Mitchell <Michael.Norman.Mitchell@gmail.com>
To   statalist@hsphsun2.harvard.edu
Subject   st: Separate intercept mixed model
Date   Tue, 8 Feb 2011 00:58:41 -0800

Greetings

  I am trying to estimate a mixed (multilevel) model using a separate
intercept approach, and am getting an error saying that the likelihood
evaluates to missing. Here is the model I am estimating with the
error...

. webuse pig, clear
(Longitudinal analysis of pig weights)
. replace week = week - 1
(432 real changes made)
. * create 2 groups, trt and control
. generate trt = (id > 24)
. * make weight 10 pounds more for treatment
. replace weight = weight + 10     if trt == 1
(216 real changes made)
. * Model 0:
. * Enter trt as factor variable
. xtmixed weight ibn.trt c.week , nocons || id:

Performing EM optimization:
likelihood evaluates to missing
r(430);

  I try to model this using a conventional coding for the variable
-trt-, and all works well...

. * Model 1:
. * Predict weight from week and treatment
. * express treatment as a factor variable
. xtmixed weight i.trt c.week || id:

Performing EM optimization:

Performing gradient-based optimization:

Iteration 0:   log restricted-likelihood = -1015.4632
Iteration 1:   log restricted-likelihood = -1015.4632

Computing standard errors:

Mixed-effects REML regression                   Number of obs      =       432
Group variable: id                              Number of groups   =        48

                                                Obs per group: min =         9
                                                               avg =       9.0
                                                               max =         9


                                                Wald chi2(2)       =  25363.93
Log restricted-likelihood = -1015.4632          Prob > chi2        =    0.0000

------------------------------------------------------------------------------
      weight |      Coef.   Std. Err.      z    P>|z|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
       1.trt |         11   1.144154     9.61   0.000     8.757499     13.2425
        week |   6.209896   .0390633   158.97   0.000     6.133333    6.286458
       _cons |   25.06551     .82399    30.42   0.000     23.45052     26.6805
------------------------------------------------------------------------------

------------------------------------------------------------------------------
  Random-effects Parameters  |   Estimate   Std. Err.     [95% Conf. Interval]
-----------------------------+------------------------------------------------
id: Identity                 |
                   sd(_cons) |    3.90138   .4198202      3.159527    4.817419
-----------------------------+------------------------------------------------
                sd(Residual) |   2.096356   .0757444      1.953034    2.250195
------------------------------------------------------------------------------
LR test vs. linear regression: chibar2(01) =   470.28 Prob >= chibar2 = 0.0000

And I am able to get the fitted mean by treatment group for week 0, shown below.

. * show means by treatment
. margins trt , at(week=0)

Adjusted predictions                              Number of obs   =        432

Expression   : Linear prediction, fixed portion, predict()
at           : week            =           0

------------------------------------------------------------------------------
             |            Delta-method
             |     Margin   Std. Err.      z    P>|z|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
         trt |
          0  |   25.06551     .82399    30.42   0.000     23.45052     26.6805
          1  |   36.06551     .82399    43.77   0.000     34.45052     37.6805
------------------------------------------------------------------------------

Now, I try modeling this using a separate intercept model, but fitting
it using manually constructed dummy variables and it fits nicely. I
feel like this model is identical to model 0, but model 0 fails and
model 2 succeeeds. Further estimates of the parameters -trt0- and
-trt1- from model 2 match the estimates from the previous margins
command from model 1.

. * Model 2:
. * Now express trt as a dummy variable
. * Use a separate intercept model, getting
. * estimates of the mean weight for week
. * 0 by treatment group
. gen trt0 = (trt==0)
. gen trt1 = (trt==1)
. xtmixed weight trt0 trt1 week, nocons || id:

Performing EM optimization:

Performing gradient-based optimization:

Iteration 0:   log restricted-likelihood = -1015.4632
Iteration 1:   log restricted-likelihood = -1015.4632

Computing standard errors:

Mixed-effects REML regression                   Number of obs      =       432
Group variable: id                              Number of groups   =        48

                                                Obs per group: min =         9
                                                               avg =       9.0
                                                               max =         9


                                                Wald chi2(3)       =  34743.66
Log restricted-likelihood = -1015.4632          Prob > chi2        =    0.0000

------------------------------------------------------------------------------
      weight |      Coef.   Std. Err.      z    P>|z|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
        trt0 |   25.06551     .82399    30.42   0.000     23.45052     26.6805
        trt1 |   36.06551     .82399    43.77   0.000     34.45052     37.6805
        week |   6.209896   .0390633   158.97   0.000     6.133333    6.286458
------------------------------------------------------------------------------

------------------------------------------------------------------------------
  Random-effects Parameters  |   Estimate   Std. Err.     [95% Conf. Interval]
-----------------------------+------------------------------------------------
id: Identity                 |
                   sd(_cons) |    3.90138   .4198202      3.159527    4.817419
-----------------------------+------------------------------------------------
                sd(Residual) |   2.096356   .0757444      1.953034    2.250195
------------------------------------------------------------------------------
LR test vs. linear regression: chibar2(01) =   470.28 Prob >= chibar2 = 0.0000

. * the estimate of the difference in means matches trt from model 1
. lincom trt1 - trt0

 ( 1)  - [weight]trt0 + [weight]trt1 = 0

------------------------------------------------------------------------------
      weight |      Coef.   Std. Err.      z    P>|z|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
         (1) |         11   1.144154     9.61   0.000     8.757499     13.2425
------------------------------------------------------------------------------

  So, I cannot figure out how to fit model 2 while expressing -trt-
using factor variable notation, i.e. ibn.trt, rather than manually
fitting dummy variables, i.e. -trt0- and -trt1-.

  When I fit model 0, omitting the random effects portion of the
model, it seems that the estimate for treatment group 0 is being
omitted even though I specified -ibn.trt- and -nocons-, see below.

. xtmixed weight ibn.trt c.week , nocons

Mixed-effects REML regression                   Number of obs      =       432

                                                Wald chi2(2)       =         .
Log restricted-likelihood =          .          Prob > chi2        =         .

------------------------------------------------------------------------------
      weight |      Coef.   Std. Err.      z    P>|z|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
         trt |
          0  |  (omitted)
          1  |   22.39341   1.078251    20.77   0.000     20.28008    24.50675
             |
        week |    9.62792   .1601441    60.12   0.000     9.314043    9.941796
------------------------------------------------------------------------------

------------------------------------------------------------------------------
  Random-effects Parameters  |   Estimate   Std. Err.     [95% Conf. Interval]
-----------------------------+------------------------------------------------
                sd(Residual) |   12.74731   .4336723      11.92505    13.62627
------------------------------------------------------------------------------

  Thank you for your patience in reading this far! I would be most
grateful for suggestions.

Many thanks,

Michael Mitchell
*
*   For searches and help try:
*   http://www.stata.com/help.cgi?search
*   http://www.stata.com/support/statalist/faq
*   http://www.ats.ucla.edu/stat/stata/


© Copyright 1996–2018 StataCorp LLC   |   Terms of use   |   Privacy   |   Contact us   |   Site index