These data are from a longitudinal growth model with two level-1 predictors, time and time^2. There are two level-2 predictors of the random intercepts, x (a continuous predictor) and grp (a 2 level categorical predictor). x is also used as a predictor of the random slopes for both time and time^2.
Multilevel models are analyzed in Stata as mixed models. Although it is not really very difficult to convert from multilevel notation to mixed notation some find the conversion a bit confusing. Therefore, we will make use of the ml2mixed program to assist in the process. If you wish to download ml2mixed you can get it by typing, search ml2mixed, in Stata’s command window and following the installation instructions. Here’s how we would use the program for our multilevel model.
ml2mixed, dep(dv) l1(time time#time) l2id(sid) l2i(grp x) l2s([time] x [time#time] x) notes Multilevel Model Level 1 Model dv = time time#time Level 2 Model -- id = sid [int] = grp x slope[time] = x slope[time#time] = x Stata Mixed Model -- Stata 11 notation xtmixed dv time time#time grp x time#x time#time#x /// || sid: time time#time , cov(unstr) Stata Notes 1) Categorical predictors need the i. prefix. 2) Continuous variables in interactions need the c. prefix. 3) Use var option to get variances instead of standard deviations. 4) If outcome variable is binary use -xtmelogit- command. 5) If outcome variable is a count use -xtmepoisson- command.
The options for the ml2mixed command above work as follows:
dep(dv) -- declare the dependent variable l1(time time#time) -- declare the level-1 predictors l2id(sid) -- declare the level-2 identifier l2i(grp x) -- declare the level-2 predictors of the random intercept l2s([time] x [time#time] x) -- declare the level-2 predictors for the random slopes notes -- request the Stata notes
Let’s begin by downloading a multilevel dataset named longitudinal.
The xtmixed command provided by ml2mixed is pretty close to what we want. We will just need to add an i. prefix to the categorical predictor, c. prefixes to each of the continuous predictors in interactions and replace the interaction term in the random slopes portion of the command with a variable t2, where t2 = time^2. We need to do this because the xtmixed command does not allow interaction notation in the random effects portion of the command.
So, here is the multilevel growth model run using xtmixed.
use https://stats.idre.ucla.edu/stat/data/longitudinal, clear xtmixed dv time c.time#c.time i.grp x c.time#c.x c.time#c.time#c.x /// || sid: time t2, cov(unstr) Performing EM optimization: Performing gradient-based optimization: Iteration 0: log restricted-likelihood = -589.52489 Iteration 1: log restricted-likelihood = -588.95211 Iteration 2: log restricted-likelihood = -588.95021 Iteration 3: log restricted-likelihood = -588.9502 Computing standard errors: Mixed-effects REML regression Number of obs = 205 Group variable: sid Number of groups = 61 Obs per group: min = 1 avg = 3.4 max = 4 Wald chi2(6) = 52.75 Log restricted-likelihood = -588.9502 Prob > chi2 = 0.0000 ------------------------------------------------------------------------------ dv | Coef. Std. Err. z P>|z| [95% Conf. Interval] -------------+---------------------------------------------------------------- time | -11.62801 5.121689 -2.27 0.023 -21.66633 -1.58968 | c.time#| c.time | 3.646695 1.450853 2.51 0.012 .8030761 6.490315 | 1.grp | -3.575066 1.133945 -3.15 0.002 -5.797559 -1.352574 x | .3374781 .1823358 1.85 0.064 -.0198936 .6948498 | c.time#c.x | .4888845 .2394994 2.04 0.041 .0194744 .9582947 | c.time#| c.time#c.x| -.1728845 .0678586 -2.55 0.011 -.3058849 -.0398841 | _cons | 9.693063 3.913496 2.48 0.013 2.022752 17.36337 ------------------------------------------------------------------------------ ------------------------------------------------------------------------------ Random-effects Parameters | Estimate Std. Err. [95% Conf. Interval] -----------------------------+------------------------------------------------ sid: Unstructured | sd(time) | 5.151503 .9204459 3.629487 7.311771 sd(t2) | 1.299955 .3063831 .8190507 2.063221 sd(_cons) | 4.614089 .5709101 3.620472 5.880399 corr(time,t2) | -.9803668 .0179769 -.9967681 -.8855039 corr(time,_cons) | -.3448718 .1728846 -.6316764 .0249725 corr(t2,_cons) | .2348463 .2098262 -.1934751 .587978 -----------------------------+------------------------------------------------ sd(Residual) | 2.556747 .2716945 2.076033 3.148772 ------------------------------------------------------------------------------ LR test vs. linear regression: chi2(6) = 85.97 Prob > chi2 = 0.0000 Note: LR test is conservative and provided only for reference.
margins grp Predictive margins Number of obs = 205 Expression : Linear prediction, fixed portion, predict() ------------------------------------------------------------------------------ | Delta-method | Margin Std. Err. z P>|z| [95% Conf. Interval] -------------+---------------------------------------------------------------- grp | 0 | 14.99115 .8771161 17.09 0.000 13.27204 16.71027 1 | 11.41609 .7612376 15.00 0.000 9.924087 12.90808 ------------------------------------------------------------------------------ margins grp, at(time=(0 1 2 3)) vsquish Predictive margins Number of obs = 205 Expression : Linear prediction, fixed portion, predict() 1._at : time = 0 2._at : time = 1 3._at : time = 2 4._at : time = 3 ------------------------------------------------------------------------------ | Delta-method | Margin Std. Err. z P>|z| [95% Conf. Interval] -------------+---------------------------------------------------------------- _at#grp | 1 0 | 16.79301 .92251 18.20 0.000 14.98492 18.6011 1 1 | 13.21794 .8392378 15.75 0.000 11.57307 14.86282 2 0 | 15.45978 .9504585 16.27 0.000 13.59692 17.32265 2 1 | 11.88472 .8481288 14.01 0.000 10.22242 13.54702 3 0 | 14.14558 1.019025 13.88 0.000 12.14833 16.14283 3 1 | 10.57051 .9087843 11.63 0.000 8.789328 12.3517 4 0 | 12.85039 .9767999 13.16 0.000 10.9359 14.76488 4 1 | 9.275325 .8497111 10.92 0.000 7.609922 10.94073 ------------------------------------------------------------------------------ margins, dydx(x) at(time=(0 1 2 3)) vsquish Average marginal effects Number of obs = 205 Expression : Linear prediction, fixed portion, predict() dy/dx w.r.t. : x 1._at : time = 0 2._at : time = 1 3._at : time = 2 4._at : time = 3 ------------------------------------------------------------------------------ | Delta-method | dy/dx Std. Err. z P>|z| [95% Conf. Interval] -------------+---------------------------------------------------------------- x | _at | 1 | .3374781 .1823358 1.85 0.064 -.0198936 .6948498 2 | .6534781 .186003 3.51 0.000 .2889189 1.018037 3 | .6237092 .2070978 3.01 0.003 .2178049 1.029614 4 | .2481713 .1914065 1.30 0.195 -.1269785 .6233211 ------------------------------------------------------------------------------ margins, dydx(time) at(x=(16(2)28)) vsquish Average marginal effects Number of obs = 205 Expression : Linear prediction, fixed portion, predict() dy/dx w.r.t. : time 1._at : x = 16 2._at : x = 18 3._at : x = 20 4._at : x = 22 5._at : x = 24 6._at : x = 26 7._at : x = 28 ------------------------------------------------------------------------------ | Delta-method | dy/dx Std. Err. z P>|z| [95% Conf. Interval] -------------+---------------------------------------------------------------- time | _at | 1 | -1.400466 .5039248 -2.78 0.005 -2.388141 -.4127917 2 | -1.367237 .386552 -3.54 0.000 -2.124865 -.6096086 3 | -1.334007 .3108527 -4.29 0.000 -1.943267 -.7247471 4 | -1.300778 .3091294 -4.21 0.000 -1.90666 -.6948954 5 | -1.267548 .3823835 -3.31 0.001 -2.017006 -.5180904 6 | -1.234319 .498596 -2.48 0.013 -2.211549 -.2570887 7 | -1.201089 .6345892 -1.89 0.058 -2.444861 .0426825 ------------------------------------------------------------------------------