UCLA Academic Technology Services HomeServicesClassesContactJobs
Help the Stat Consulting Group by giving a gift             
Loading

Stata FAQ
How can I access the random effects after xtmixed using _diparm?

Consider the following xtmixed model:
use http://www.ats.ucla.edu/stat/data/depression, clear

xtmixed deprisson i.group visit || sid: visit,  cov(un)

Performing EM optimization: 

Performing gradient-based optimization: 

Iteration 0:   log restricted-likelihood =  -825.5536  
Iteration 1:   log restricted-likelihood = -825.55355  

Computing standard errors:

Mixed-effects REML regression                   Number of obs      =       295
Group variable: subj                            Number of groups   =        61

                                                Obs per group: min =         1
                                                               avg =       4.8
                                                               max =         6


                                                Wald chi2(2)       =     63.98
Log restricted-likelihood = -825.55355          Prob > chi2        =    0.0000

------------------------------------------------------------------------------
  depression |      Coef.   Std. Err.      z    P>|z|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
     1.group |  -3.793401   1.201829    -3.16   0.002    -6.148942    -1.43786
       visit |  -1.202648   .1670519    -7.20   0.000    -1.530064   -.8752327
       _cons |    17.9614   1.009313    17.80   0.000     15.98319    19.93962
------------------------------------------------------------------------------

------------------------------------------------------------------------------
  Random-effects Parameters  |   Estimate   Std. Err.     [95% Conf. Interval]
-----------------------------+------------------------------------------------
subj: Unstructured           |
                   sd(visit) |   .9272044   .1567983      .6656254    1.291579
                   sd(_cons) |    5.15685   .6236461      4.068591    6.536194
           corr(visit,_cons) |   -.548247   .1364845      -.760894   -.2292597
-----------------------------+------------------------------------------------
                sd(Residual) |   2.891413   .1503662      2.611222    3.201668
------------------------------------------------------------------------------
LR test vs. linear regression:       chi2(3) =   176.17   Prob > chi2 = 0.0000

Note: LR test is conservative and provided only for reference.
Say that you want to use the random effect sd(visit) or cov(visit,_cons) for additional computations. You can access these values using the undocumented command _diparm (which stands for display parameter). When we say undocumented, we mean that it is not in the full manual but it is included in the online help system (help _diparm).

To use _diparm you have to understand how Stata computes the random effects. Stata computes the variances as the log of the standard deviation (ln_sigma) and computes covariances as the arc hyperbolic tangent of the correlation.

You also need to how stmixed names the random effects. The coeflegend option will not provide these names. The easiest way to get the names of the random effects is to list of the e(b) matrix, like this.

matrix list e(b)

e(b)[1,8]
     depression:  depression:  depression:  depression:    lns1_1_1:    lns1_1_2:  atr1_1_1_2:     lnsig_e:
             0b.           1.                                                                              
          group        group        visit        _cons        _cons        _cons        _cons        _cons
y1            0   -3.7934014   -1.2026483    17.961403   -.07558125    1.6403258   -.61587151    1.0617452
Thus, lns1_1_1 is the name for sd(visit), lns1_1_2 is the name for sd(_cons) and atr1_1_1_2 is the name for corr(visit,_cons).

Finally, you need to provide _diparm with three pieces of information: 1) the name of the parameter, 2) the inverse function (abbr f), and 3) the derivative (abbr d) of the inverse function.

For example, to get the estimate for sd(visit) the name is lns1_1_1. Since xtmixed estimates the ln_sigma, the inverse function is exp(). The derivative is easy because the derivative of exp() is just exp(). Here is _diparm in action.

_diparm lns1_1_1, f(exp(@)) d(exp(@))
 
   /lns1_1_1 |   .9272044   .1567983                      .6656254    1.291579l
Note the use of the @ symbol as a place holder for the argument of the function.

All of the information from the _diparm is stored in the return list.

return list

scalars:
                  r(z) =  .
                  r(p) =  .
                  r(i) =  0
                 r(ub) =  1.291579210210525
                 r(lb) =  .6656254461991183
                r(est) =  .9272043939164035
                 r(se) =  .156798306133924
For corr(visit,_cons) the name is atr1_1_1_2, the inverse function is tanh() and the derivative is 1-tanh()^2.
_diparm atr1_1_1_2, f(tanh(@)) d(1-tanh(@)^2)

 /atr1_1_1_2 |   -.548247   .1364845                      -.760894   -.2292597
Next we will rerun the xtmixed using the var option.
xtmixed, var

Mixed-effects REML regression                   Number of obs      =       295
Group variable: subj                            Number of groups   =        61

                                                Obs per group: min =         1
                                                               avg =       4.8
                                                               max =         6


                                                Wald chi2(2)       =     63.98
Log restricted-likelihood = -825.55355          Prob > chi2        =    0.0000

------------------------------------------------------------------------------
  depression |      Coef.   Std. Err.      z    P>|z|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
     1.group |  -3.793401   1.201829    -3.16   0.002    -6.148942    -1.43786
       visit |  -1.202648   .1670519    -7.20   0.000    -1.530064   -.8752327
       _cons |    17.9614   1.009313    17.80   0.000     15.98319    19.93962
------------------------------------------------------------------------------

------------------------------------------------------------------------------
  Random-effects Parameters  |   Estimate   Std. Err.     [95% Conf. Interval]
-----------------------------+------------------------------------------------
subj: Unstructured           |
                  var(visit) |    .859708   .2907682      .4430572    1.668177
                  var(_cons) |    26.5931   6.432099      16.55343    42.72183
            cov(visit,_cons) |  -2.621418   1.158134     -4.891318    -.351517
-----------------------------+------------------------------------------------
               var(Residual) |   8.360267   .8695415      6.818482    10.25068
------------------------------------------------------------------------------
LR test vs. linear regression:       chi2(3) =   176.17   Prob > chi2 = 0.0000

Note: LR test is conservative and provided only for reference.
We could always get the variance by squaring the standard deviation we computed earlier. But, that wouldn't give us the standard error or confidence intervals. To get all of these values, we will rerun _diparm with a slightly different inverse function and derivative.
_diparm lns1_1_1, f(exp(@)^2) d(2*exp(@)^2)

   /lns1_1_1 |    .859708   .2907682                      .4430572    1.668177
Getting cov(visit,_cons) is a bit more work. Computing a covariance from a correlation is just a matter of multiplying the correlation by the two standard deviations. Doing this with _diparmvmeans that we will need to use three arguments in the command. We will need to use three place holders for the argument; @1, @2 and @3.
_diparm atr1_1_1_2 lns1_1_1 lns1_1_2, f(tanh(@1)*exp(@2+@3))  ///
        d((1-tanh(@1)^2)*exp(@2+@3) tanh(@1)*exp(@2+@3) tanh(@1)*exp(@2+@3))

 /atr1_1_1_2 |  -2.621418   1.158134                     -4.891318    -.351517
This time the inverse function had one term with three arguments while the derivative had three terms, one for each argument. Also note that exp(@2)*exp(@3) was expressed as exp(@2+@3).

How to cite this page

Report an error on this page or leave a comment

UCLA Researchers are invited to our Statistical Consulting Services
We recommend others to our list of Other Resources for Statistical Computing Help
These pages are Copyrighted (c) by UCLA Academic Technology Services


The content of this web site should not be construed as an endorsement of any particular web site, book, or software product by the University of California.