Raphael Fraser wrote:
I am trying to reproduce the code below in SAS into Stata. proc mixed
by default uses restricted maximum likelihood. The model estimates the
fixed effects and standard error.
The options R and RCORR estimates the covariance and correlation
matrices respectively.
I can't seem to replicate the standard errors. Is there anyone that can help?
PROC MIXED;
CLASS id group time;
MODEL y = group time group*time/ S CHISQ;
REPEATED time/ TYPE = UN SUBJECT = id R RCORR;
--------------------------------------------------------------------------------
The PROC MIXED code calls for an unstructured patterned covariance matrix on
the errors. Stata's -xtmixed- can't do that. -xtmixed- can apply an
unstructured covariance matrix only to higher-level random effects.
Someone showed a couple of years ago how you can get what you want, though.
(I can't find the post in the Statalist archives at the moment, but I believe
it was Yulia Marchenko who posted the method.) You basically insert a fake
extra level so that all except the common variance of the residuals are now at
the second level, where -xtmixed- can handle them as higher-level random
effects.
You can compare the results of the do-file below to the corresponding
SAS PROC MIXED output shown in
http://users.uoa.gr/~fsiannis/SAS/longitudinal/lab/tlc_profileanalysis.pdf to
see that they're the same, including the standard errors of the fixed-effects
estimates. (The SAS PROC MIXED code shown there happens to be identical to
what you show above except for the name of the response variable.)
Joseph Coveney
clear *
set more off
infile id str1 group lead0 lead1 lead4 lead6 using ///
http://users.uoa.gr/~fsiannis/data/longitudinal/tlcdata.txt
quietly {
reshape long lead, i(id) j(time)
tabulate time, generate(time_)
compress
}
* SAS drops the high values:
char group[omit] "P"
char time[omit] 6
xi: xtmixed lead i.group*i.time || id: time_*, covariance(unstructured) ///
noconstant nostderr variance nolrtest nolog
* SAS displays the deviance, and not the log-likelihood:
display in smcl as result %13.8f -2 * e(ll)
* Covariances can be taken directly from the -xtmixed- output, but
* variances need to be re-assembled:
matrix define B = e(b)
tempname temporary_matrix
matrix define `temporary_matrix' = B[1,"lnsig_e:_cons"]
scalar define s2_residual = exp(`temporary_matrix'[1,1])^2
forvalues time_point = 1/4 {
matrix define `temporary_matrix' = B[1,"lns1_1_`time_point':_cons"]
scalar define s2_time`time_point' = exp(`temporary_matrix'[1,1])^2 + ///
s2_residual
display in smcl as result %7.4f s2_time`time_point'
}
exit
*
* 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/