Login

News

Converting a Stata -sem- path model

Posted on
Picture of user. Robert Long Joined: 02/22/2013
Forums
Converting to and from OpenMx
Hi all

I am trying to convert a simple path model from Stata. In Stata I use

clear all
ssd init EE DR RO SE
ssd set observations 321

ssd set sd 15.23 10.15 16.27 1.01

ssd set correlations ///
1 \ ///
0.25 1 \ ///
0.78 0.21 1 \ ///
-0.47 -0.27 -0.35 1

sem (SE -> EE) (EE -> DR) (EE -> RO) (DR -> RO), nocapslatent

This generates the following output


------------------------------------------------------------------------------
| OIM
| Coef. Std. Err. z P>|z| [95% Conf. Interval]
-------------+----------------------------------------------------------------
Structural |
EE <- |
SE | -7.087228 .7428868 -9.54 0.000 -8.543259 -5.631196
_cons | 59.18695 3.629161 16.31 0.000 52.07392 66.29997
-----------+----------------------------------------------------------------
DR <- |
EE | .166612 .0360163 4.63 0.000 .0960213 .2372026
_cons | 16.53305 1.063443 15.55 0.000 14.44874 18.61736
-----------+----------------------------------------------------------------
RO <- |
EE | .8289902 .0385246 21.52 0.000 .7534834 .9044969
DR | .0256473 .0578058 0.44 0.657 -.08765 .1389446
_cons | 51.46608 1.458225 35.29 0.000 48.60801 54.32415
-------------+----------------------------------------------------------------
Variance |
e.EE | 180.1515 14.22003 154.3297 210.2938
e.DR | 96.28271 7.599951 82.48212 112.3924
e.RO | 103.2753 8.151902 88.47244 120.5549
------------------------------------------------------------------------------
LR test of model vs. saturated: chi2(2) = 10.78, Prob > chi2 = 0.0046

My attempt to reproduce this in OpenMX is


vars <- c("SE","EE","DR","RO")

cormat <- matrix(c(1,0.25, 0.78, -0.47,
0.25, 1, 0.21, -0.27,
0.78, 0.21 , 1, -0.35,
-0.47, -0.27, -0.35, 1
), nrow=4)

sds <- c(15.23, 10.15, 16.27, 1.01)

rownames(cormat) <- colnames(cormat) <- vars
covmat <- round(cor2cov(cormat,sds),2)

m.cleary <- mxModel("mcleary",
manifestVars = vars,
type="RAM",
mxPath(
from=vars[c(1,2,3,2)], to=vars[c(2,3,4,4)],
values=1,
labels=c("p1","p2","p3","p4")
),
mxPath(
from=vars, # “to=vars” is implicit
arrows=2, # the default is 1
values=.9,
labels=paste("var", vars, sep="")
),
mxData(
observed=covmat,
type="cov", # data is a covariance matrix
numObs=321
)
)

f.cleary <- mxRun(m.cleary)
summary(f.cleary)

and this produced:


data:
$mcleary.data
$mcleary.data$cov
SE EE DR RO
SE 231.95 38.65 193.28 -7.23
EE 38.65 103.02 34.68 -2.77
DR 193.28 34.68 264.71 -5.75
RO -7.23 -2.77 -5.75 1.02

free parameters:
name matrix row col Estimate Std.Error Std.Estimate Std.SE lbound ubound
1 p1 A EE SE 0.16663073 0.036072013 0.2500296 0.05412611
2 p2 A DR EE 0.33663366 0.087610155 0.2100067 0.05465502
3 p4 A RO EE -0.02047884 0.005205191 -0.2058098 0.05231152
4 p3 A RO DR -0.01903892 0.003247224 -0.3067097 0.05231152
5 varSE S SE SE 231.94999676 18.337149808 1.0000000 0.07905648
6 varEE S EE EE 96.57972313 7.635160662 0.9374852 0.07411338
7 varDR S DR DR 253.03554502 20.004099673 0.9558972 0.07556987
8 varRO S RO RO 0.85379980 0.067498669 0.8370586 0.06617517

observed statistics: 10
estimated parameters: 8
degrees of freedom: 2
-2 log likelihood: 6205.556
saturated -2 log likelihood: 5886.007
number of observations: 321
chi-square: 319.5493
p: 4.080808e-70

Could anyone advise where I am going wrong ?

Thanks
RL

Replied on Mon, 02/25/2013 - 17:18
Picture of user. mhunter Joined: 07/31/2009

At first glance everything looks like it should be correct. The only thing I can see that might be odd is the "_cons" term that goes with all the coefficients in the Stata output you provided. I'm guessing it's a kind of constant/intercept term. This was not in the output for the previous model you asked about. http://openmx.psyc.virginia.edu/thread/1934.

The "_cons" term makes me think the Stata output was generated off raw data, and/or observed means. If so, then the corresponding OpenMx model would need the same information. Could that be the case?

Log in or register to post comments
Replied on Mon, 02/25/2013 - 18:39
Picture of user. Robert Long Joined: 02/22/2013

I think I may have inadvertently posted the output for the model which ran from raw data instead of the covariance matrix, although the path coefficients and variances are the same, so that isn't part of the problem...

However, I have now solved it ! I made the mistake of changing the order of the variables (and standard deviations), but not changing the covariance matrix !! So the correct model is


vars <- c("EE","DR","RO","SE")

cormat <- matrix(c(1,0.25, 0.78, -0.47,
0.25, 1, 0.21, -0.27,
0.78, 0.21 , 1, -0.35,
-0.47, -0.27, -0.35, 1
), nrow=4)

sds <- c(15.23, 10.15, 16.27, 1.01)

rownames(cormat) <- colnames(cormat) <- vars
covmat <- round(cor2cov(cormat,sds),2)

m.cleary <- mxModel("mcleary",
manifestVars = vars,
type="RAM",
mxPath(
from=vars[c(4,1,1,2)], to=vars[c(1,2,3,3)],
values=1,
labels=c("p1","p2","p3","p4")
),
mxPath(
from=vars,
arrows=2,
values=.9,
labels=paste("var", vars, sep="")
),

mxData(
observed=covmat,
type="cov",
numObs=321
)
)

f.cleary <- mxRun(m.cleary)
summary(f.cleary)

And this now agrees with Stata.

Thanks again for your help
RL

Log in or register to post comments
  Login

News