I am fitting the mixed model below:
. mixed y trt || clst:trt, nocons reml dfmethod(sat)
Performing EM optimization:
Performing gradient-based optimization:
Iteration 0: log restricted-likelihood = -1295.3123
Iteration 1: log restricted-likelihood = -1295.3098
Iteration 2: log restricted-likelihood = -1295.3098
Computing standard errors:
Computing degrees of freedom:
Mixed-effects REML regression Number of obs = 919
Group variable: clst Number of groups = 49
Obs per group:
min = 1
avg = 18.8
max = 30
DF method: Satterthwaite DF: min = 888.00
avg = 900.91
max = 913.83
F(1, 913.83) = 0.40
Log restricted-likelihood = -1295.3098 Prob > F = 0.5251
------------------------------------------------------------------------------
y | Coef. Std. Err. t P>|t| [95% Conf. Interval]
-------------+----------------------------------------------------------------
trt | .1455914 .2290005 0.64 0.525 -.3038366 .5950193
_cons | .3951269 .2241477 1.76 0.078 -.0447941 .835048
------------------------------------------------------------------------------
------------------------------------------------------------------------------
Random-effects Parameters | Estimate Std. Err. [95% Conf. Interval]
-----------------------------+------------------------------------------------
clst: Identity |
var(trt) | .0341507 .0173905 .0125877 .092652
-----------------------------+------------------------------------------------
var(Residual) | .9546016 .0453034 .8698131 1.047655
------------------------------------------------------------------------------
LR test vs. linear model: chibar2(01) = 9.46 Prob >= chibar2 = 0.0010
. return list
scalars:
r(level) = 95
matrices:
r(table) : 9 x 4
Next, I calculate the ICC as follows:
. nlcom (icc_est: (exp(_b[lns1_1_1:_cons])^2)/((exp(_b[lns1_1_1:_cons])^2)+(exp(_b[lnsig_e:_cons])^2)))
icc_est: (exp(_b[lns1_1_1:_cons])^2)/((exp(_b[lns1_1_1:_cons])^2)+(exp(_b[lnsig_e:_cons])^2))
------------------------------------------------------------------------------
y | Coef. Std. Err. z P>|z| [95% Conf. Interval]
-------------+----------------------------------------------------------------
icc_est | .0345392 .0171907 2.01 0.045 .0008461 .0682323
------------------------------------------------------------------------------
How can I save the results in the dataset?
I want to keep all the three tables shown: fixed effects, random effects and the ICC results.
Consider the following reproducible example using Stata's pig toy dataset:
webuse pig, clear
mixed weight week || id:week, nocons reml dfmethod(sat)
nlcom (icc_est: (exp(_b[lns1_1_1:_cons])^2)/((exp(_b[lns1_1_1:_cons])^2)+(exp(_b[lnsig_e:_cons])^2))), post
------------------------------------------------------------------------------
weight | Coef. Std. Err. z P>|z| [95% Conf. Interval]
-------------+----------------------------------------------------------------
icc_est | .1380299 .0265754 5.19 0.000 .0859431 .1901167
------------------------------------------------------------------------------
The following works for me:
generate double coef = _b[icc_est]
generate double se = _se[icc_est]
generate p = string(2 * (normal(-(_b[icc_est] / _se[icc_est]))), "%9.3f")
generate double upper = _b[icc_est] + _se[icc_est] * invnormal(0.025)
generate double lower = _b[icc_est] + _se[icc_est] * invnormal(0.975)
list coef se p upper lower in 1
+-------------------------------------------------------+
| coef se p upper lower |
|-------------------------------------------------------|
1. | .13802987 .02657538 0.000 .08594308 .19011667 |
+-------------------------------------------------------+
save mydata.dta
The process is similar for the results of the main model.