You are here

multilevel path/structural equation models in OpenMx

33 posts / 0 new
Last post
Soyoung's picture
Offline
Joined: 05/08/2017 - 11:45
multilevel path/structural equation models in OpenMx
AttachmentSize
File Data1.26 MB
PDF icon Mplus code158.41 KB

Hi,

Is it possible to estimate multilevel path/structural equation models with OpenMx?
If it is possible, is there an example code for multilevel mediation (path) model?

Attached files are the example of multilevel path model using Mplus.
Source: Heck, R. H., & Thomas, S. L. (2015). An introduction to multilevel modeling techniques: MLM and SEM approaches using Mplus. Routledge.

Hope I can do the same analysis using OpenMx.

Thank you in advance.
Soyoung.

AdminHunter's picture
Offline
Joined: 03/01/2013 - 11:03
Yep

It is possible to estimate multilevel structural equations in OpenMx. This paper by Pritikin et al. discusses many of the details and provides some example scripts.

Many-Level Multilevel Structural Equation Modeling: An Efficient Evaluation Strategy

AdminJosh's picture
Offline
Joined: 12/12/2012 - 12:15
voila

Here's a similar model (see attachment). You probably need to add/remove a few paths to make it identical to what you posted. I find it tedious to read Mplus model diagrams or tables.

File attachments: 
Soyoung's picture
Offline
Joined: 05/08/2017 - 11:45
fabulous

Thank you all.

Two more questions.
Please check the attached file.

[1]
I revised the code, and i though I made the same model.
BUT I found these results are not identical in terms of estimates.
Why did these have different results?
Which Estimation method the OpenMx use?
Mplus uses Robust-ML.

[2]
I want to estimate indirect effects in Multilevel SEM.
How can I set a product term (indirect effect parameter) in OpenMx?

AdminJosh's picture
Offline
Joined: 12/12/2012 - 12:15
questions
  1. See imxRobustSE

  2. Create an mxAlgebra and use mxCI or mxSE to request confidence intervals on it. See attached code.

File attachments: 
Soyoung's picture
Offline
Joined: 05/08/2017 - 11:45
Thank you

Thank you so much for your kind response.
It was very helpful.

All the best,
Soyoung.

AdminNeale's picture
Offline
Joined: 03/01/2013 - 14:09
Agreement

Hi Soyoung

Could you confirm that you were able to get the same estimates and goodness-of-fit statistics with the multilevel model in OpenMx when compared to Mplus?

Thanks
Mike

Soyoung's picture
Offline
Joined: 05/08/2017 - 11:45
could not get the identical results

Hi, Mike.
I still could not get the identical results.
I am curious why these are different.

First, in terms of model fit statistics, two results were almost same.

  Mplus OpenMx
RMSEA 0 0
CFI 1 NA
TLI 1.001 1

BUT I don't think this example is suitable here, because these fit values are too perfect. It should be tested in the other data.

 

Second, in terms of estimates, I could not get the same estimates.
I'm trying to figure out why there are different results.

File attachments: 
AdminRobK's picture
Offline
Joined: 01/24/2014 - 12:15
imxRobustSE()

I am pretty sure that imxRobustSE() does not correctly handle multilevel data; the matter of whether or not it does is currently a matter of disagreement among the OpenMx development team. I suggest you not rely too much on the robust SEs you get out of OpenMx.

Soyoung, could you post the OpenMx script you're currently working from?

Soyoung's picture
Offline
Joined: 05/08/2017 - 11:45
sure

Hi Rob

Following code is based on what Dr.Pritikin provided.
You can find the data in post #1 .

library(OpenMx)
flatData <- read.table("ch4mv.dat")
colnames(flatData) <- c("orgid",'female','white',
                   'satpay','morale','org1','org2','benefit','cond',
                   'resour','zproduct','lev1wt','lev2wt')
 
level1 <- flatData[,c('benefit','cond','female','white','orgid')]
level2 <- flatData[!duplicated(flatData$orgid),c('orgid', 'org1','org2','zproduct')]
 
orgModel <- mxModel(
  "org", type="RAM",
  mxData(level2, 'raw', primaryKey = "orgid"),
  manifestVars = c('zproduct'),
  latentVars = c('Benefits','Conditions', 'org1','org2'),
  mxPath(c('zproduct','Benefits','Conditions'), arrows=2, values=1),
  mxPath('Benefits','Conditions', arrows=2),
  mxPath('one', 'org1', labels='data.org1', free=FALSE),
  mxPath('one', 'org2', labels='data.org2', free=FALSE),
  mxPath('org1', c('Benefits', 'Conditions')),
  mxPath('org2', 'zproduct', labels="a"),
  mxPath('one', 'zproduct'),
  mxPath('zproduct', c('Benefits','Conditions'), labels=paste0('b',1:2)),
  mxMatrix(nrow=2, ncol=1, free=TRUE, labels=paste0('b',1:2), name='B'),
  mxAlgebra(B * a, name="indirect")
)
 
empModel <- mxModel(
  "emp", type="RAM", orgModel,
  mxData(level1, 'raw'),
  manifestVars = c('benefit','cond'),
  latentVars = c('female','white'),
  mxPath('one','female', labels='data.female', free=FALSE),
  mxPath('one','white', labels='data.white', free=FALSE),
  mxPath(c('female','white'), c('benefit','cond'), connect = "all.bivariate"),
  mxPath(c('benefit','cond'), arrows=2, connect = "unique.pairs", values=c(1,0,1)),
  mxPath('org.Benefits', 'benefit', values=1, free=FALSE, joinKey = "orgid"),
  mxPath('org.Conditions', 'cond', values=1, free=FALSE, joinKey = "orgid"),
  mxPath('one', c('benefit','cond')),
  mxCI("org.indirect")
)
 
ch4M3 <- mxRun(empModel, intervals=TRUE)
summary(ch4M3)
imxRobustSE(ch4M3)
AdminRobK's picture
Offline
Joined: 01/24/2014 - 12:15
Under x86_64 Linux/GNU, I

Under x86_64 Linux/GNU, I substantively reproduce the OpenMx results you report with all 3 gradient-descent optimizers (script and output are attached).

Could you post the full output from MPlus? In particular, I'm curious if the -2logL matches.

Soyoung's picture
Offline
Joined: 05/08/2017 - 11:45
full output

Yes. Attachment is full output.

Mplus reports LL -40521.828.
(-2) * (-40521.828) = 81043.656

Thank you.
Soyoung.

File attachments: 
AdminNeale's picture
Offline
Joined: 03/01/2013 - 14:09
ML vs MLR?

Hi Soyoung, I notice that the estimator in MPlus is MLR, whereas in OpenMx it is FIML. Could you try it with ML (FIML) in Mplus to see if that changes the results?

Cheers
Mike

Soyoung's picture
Offline
Joined: 05/08/2017 - 11:45
different

Hi Mike.
Attachment is the result of FIML,
but this is almost same as the result of MLR of Mplus,
and the estimates does not match the result of OpenMx.

I wonder why the two results are different...Where does the difference come from?
Is there a way to check local-maxima in OpenMx, by any chance?

Thank you.
Soyoung.

File attachments: 
jpritikin's picture
Offline
Joined: 05/24/2012 - 00:35
fit value

One way to debug this is to load the parameter vectors into the OpenMx model using omxSetParameters and then evaluate the fit at that point using mxRun(mxModel(model, mxComputeOnce('fitfunction','fit'))). If the Mplus fit matches at a specific parameter vector then the models are likely the same. If the fit doesn't match then the models are different.

AdminRobK's picture
Offline
Joined: 01/24/2014 - 12:15
already tried it
One way to debug this is to load the parameter vectors into the OpenMx model using omxSetParameters and then evaluate the fit at that point using mxRun(mxModel(model, mxComputeOnce('fitfunction','fit'))).

I already tried that (see my post, presently #17, below). The fit values don't match.

AdminNeale's picture
Offline
Joined: 03/01/2013 - 14:09
Local maxima possibility

Hi Soyoung

OpenMx does have facilities for testing for local minima, such as mxTryHard(). However, in this case it seems to me that Mplus is the software that might need it, because the fit function from OpenMx's estimation is indicating a much higher likelihood (41000.78 vs 40522.071) than is Mplus. You could perhaps try putting the OpenMx estimates into Mplus and fixing all the parameters and seeing what it thinks the answer is with those values. If the answers still don't agree, it indicates that the models are not the same, and we should try to figure out why they aren't.

Cheers
Mike

jpritikin's picture
Offline
Joined: 05/24/2012 - 00:35
robust

Good point Rob. I am not sure what Robust-ML means for multilevel models.

AdminRobK's picture
Offline
Joined: 01/24/2014 - 12:15
It appears that imxRobustSE()

It appears that imxRobustSE() does not calculate robust SEs in a manner theoretically correct for multilevel data (e.g., Freedman 2006).

It also appears that what OpenMx and MPlus are doing in this case are different enough that their results shouldn't be compared to each other. See the attached script and output. When OpenMx calculates the fit at MPlus' solution, OpenMx returns a different fit from MPlus. Further, I can start OpenMx at MPlus' solution, and OpenMx still returns the same solution it has before (with a fit value of 82001.55).

I don't know OpenMx's multilevel feature very well, and I don't know MPlus at all, so I'm not sure what might explain the difference (or even if Soyoung's syntax is even fitting equivalent models in the two software packages). Perhaps a difference in how the two programs internally handle zero-variance exogenous regressors? Does MPlus drop normalizing constants from its loglikelihood?

Soyoung's picture
Offline
Joined: 05/08/2017 - 11:45
I have been told that the

I have been told that the Mplus software excludes any variable with no variance from analysis.
and Mplus does not drop normalizing constants.

In my eyes, OpenMx and Mplus are different in the manner of providing results.
In OpenMx, when there are level-1 manifest variables causes level-2 latent variables,
intercepts are showed at the lower level results.
In Mplus, when there are level-1 manifest variables causes level-2 latent variables,
intercepts are showed at the between level (upper level) results.
When I look at the mathematical formulas of the model, I feel easy to understand to treat the intercept at the upper level.
(probably because I got used to it)

Will the outcome differ depending on the level of handling the intercept?

Thank you again,
Soyoung.

jpritikin's picture
Offline
Joined: 05/24/2012 - 00:35
means

> Will the outcome differ depending on the level of handling the intercept?

No and you can prove it by trying it both ways.

jpritikin's picture
Offline
Joined: 05/24/2012 - 00:35
bug

OK, there is an OpenMx bug. Try the model again with,

empModel$expectation$.useSufficientSets <- FALSE
jpritikin's picture
Offline
Joined: 05/24/2012 - 00:35
looks better now
Summary of emp 
 
free parameters:
         name matrix        row        col    Estimate  Std.Error A
1  emp.A[1,3]      A    benefit     female  0.69700745 0.02567046  
2  emp.A[2,3]      A       cond     female  0.70677248 0.02514620  
3  emp.A[1,4]      A    benefit      white  0.22406279 0.02694068  
4  emp.A[2,4]      A       cond      white  0.24780137 0.02630581  
5  emp.S[1,1]      S    benefit    benefit  2.02929959 0.02590178  
6  emp.S[1,2]      S    benefit       cond  1.34397122 0.02165825  
7  emp.S[2,2]      S       cond       cond  1.94820623 0.02485886  
8  emp.M[1,1]      M          1    benefit  4.41229957 0.04570310  
9  emp.M[1,2]      M          1       cond  4.67685717 0.04025743  
10         b1  org.A   Benefits   zproduct  0.31132159 0.03144856  
11         b2  org.A Conditions   zproduct  0.29100892 0.02721743  
12 org.A[2,4]  org.A   Benefits       org1  0.28727342 0.07941363  
13 org.A[3,4]  org.A Conditions       org1  0.31525677 0.06869425  
14          a  org.A   zproduct       org2 -0.45125514 0.20589852  
15 org.A[2,5]  org.A   Benefits       org2  0.03564603 0.08489777  
16 org.A[3,5]  org.A Conditions       org2  0.02533898 0.07302397  
17 org.S[1,1]  org.S   zproduct   zproduct  1.07059845 0.11934389 !
18 org.S[2,2]  org.S   Benefits   Benefits  0.13261944 0.01900351  
19 org.S[2,3]  org.S   Benefits Conditions  0.09517323 0.01495373  
20 org.S[3,3]  org.S Conditions Conditions  0.09256918 0.01401026  
21 org.M[1,1]  org.M          1   zproduct  0.08692317 0.08979952  
 
Model Statistics: 
               |  Parameters  |  Degrees of Freedom  |  Fit (-2lnL units)
       Model:             21                  25029              81043.97
   Saturated:             NA                     NA                    NA
Independence:             NA                     NA                    NA
Number of observations/statistics: 12605/25050
Soyoung's picture
Offline
Joined: 05/08/2017 - 11:45
Estimates are almost same/Where is the result of indirect effect

Hi Joshua

I am so pleased to see these estimates are almost same as the results of Mplus :)

Where are the results of indirect effects such as a*b1 and a*b2?
I am curious whether they are also same as the results of Mplus.

Mplus (FIML)
 
Loglikelihood 
H0 Value   -40521.984
H1 Value   -40521.828

Thank you so much,
Soyoung.

Soyoung's picture
Offline
Joined: 05/08/2017 - 11:45
Please never mind.

Please never mind.
I am going to try! :)

Soyoung's picture
Offline
Joined: 05/08/2017 - 11:45
almost same

Hello,

l'd like to confirm that I can get the same estimates.
and I also want to ask a question if I can ignore small difference of lbound/ubound.

When I ran the following code,

library(OpenMx)
flatData <- read.table("ch4mv.dat")
colnames(flatData) <- c("orgid",'female','white',
                   'satpay','morale','org1','org2','benefit','cond',
                   'resour','zproduct','lev1wt','lev2wt')
 
level1 <- flatData[,c('benefit','cond','female','white','orgid')]
level2 <- flatData[!duplicated(flatData$orgid),c('orgid', 'org1','org2','zproduct')]
 
orgModel <- mxModel(
  "org", type="RAM",
  mxData(level2, 'raw', primaryKey = "orgid"),
  manifestVars = c('zproduct'),
  latentVars = c('Benefits','Conditions', 'org1','org2'),
  mxPath(c('zproduct','Benefits','Conditions'), arrows=2, values=1),
  mxPath('Benefits','Conditions', arrows=2),
  mxPath('one', 'org1', labels='data.org1', free=FALSE),
  mxPath('one', 'org2', labels='data.org2', free=FALSE),
  mxPath('org1', c('Benefits', 'Conditions')),
  mxPath('org2', 'zproduct', labels="a"),
  mxPath('one', 'zproduct'),
  mxPath('zproduct', c('Benefits','Conditions'), labels=paste0('b',1:2)),
  mxMatrix(nrow=2, ncol=1, free=TRUE, labels=paste0('b',1:2), name='B'),
  mxAlgebra(B * a, name="indirect")
)
 
empModel <- mxModel(
  "emp", type="RAM", orgModel,
  mxData(level1, 'raw'),
  manifestVars = c('benefit','cond'),
  latentVars = c('female','white'),
  mxPath('one','female', labels='data.female', free=FALSE),
  mxPath('one','white', labels='data.white', free=FALSE),
  mxPath(c('female','white'), c('benefit','cond'), connect = "all.bivariate"),
  mxPath(c('benefit','cond'), arrows=2, connect = "unique.pairs", values=c(1,0,1)),
  mxPath('org.Benefits', 'benefit', values=1, free=FALSE, joinKey = "orgid"),
  mxPath('org.Conditions', 'cond', values=1, free=FALSE, joinKey = "orgid"),
  mxPath('one', c('benefit','cond')),
  mxCI("org.indirect")
)
empModel$expectation$.useSufficientSets <- FALSE
ch4M3 <- mxRun(empModel, intervals=TRUE)
summary(ch4M3)

I got the following result.

Summary of emp 
 
free parameters:
         name matrix        row        col    Estimate  Std.Error A
1  emp.A[1,3]      A    benefit     female  0.69705464 0.02567670  
2  emp.A[2,3]      A       cond     female  0.70680474 0.02515014  
3  emp.A[1,4]      A    benefit      white  0.22419999 0.02695238  
4  emp.A[2,4]      A       cond      white  0.24789392 0.02630992  
5  emp.S[1,1]      S    benefit    benefit  2.02928474 0.02590308  
6  emp.S[1,2]      S    benefit       cond  1.34396875 0.02165891  
7  emp.S[2,2]      S       cond       cond  1.94820919 0.02486024  
8  emp.M[1,1]      M          1    benefit  4.42133171 0.04035614  
9  emp.M[1,2]      M          1       cond  4.68326778 0.03568958  
10         b1  org.A   Benefits   zproduct  0.30887125 0.03092515  
11         b2  org.A Conditions   zproduct  0.28927449 0.02673664  
12 org.A[2,4]  org.A   Benefits       org1  0.27821530 0.07658434  
13 org.A[3,4]  org.A Conditions       org1  0.30882772 0.06618441  
14          a  org.A   zproduct       org2 -0.45125100 0.21303763  
15 org.S[1,1]  org.S   zproduct   zproduct  1.07059364 0.11968246  
16 org.S[2,2]  org.S   Benefits   Benefits  0.13288795 0.01902547  
17 org.S[2,3]  org.S   Benefits Conditions  0.09533848 0.01496872  
18 org.S[3,3]  org.S Conditions Conditions  0.09267003 0.01402068  
19 org.M[1,1]  org.M          1   zproduct  0.08692014 0.09164268  
 
confidence intervals:
                      lbound   estimate      ubound note
org.indirect[1,1] -0.2768346 -0.1393785 -0.01158206     
org.indirect[2,1] -0.2580154 -0.1305354 -0.01085831     
 
Model Statistics: 
               |  Parameters  |  Degrees of Freedom  |  Fit (-2lnL units)
       Model:             19                  25031              81044.14
   Saturated:             NA                     NA                    NA
Independence:             NA                     NA                    NA
Number of observations/statistics: 12605/25050
 
Information Criteria: 
      |  df Penalty  |  Parameters Penalty  |  Sample-Size Adjusted
AIC:       30982.14               81082.14                       NA
BIC:     -155294.78               81223.54                 81163.16
CFI: NA 
TLI: 1   (also known as NNFI) 
RMSEA:  0  [95% CI (NA, NA)]
Prob(RMSEA <= 0.05): NA
To get additional fit indices, see help(mxRefModels)
timestamp: 2017-06-11 16:38:45 
Wall clock time (HH:MM:SS.hh): 00:02:07.38 
optimizer:  CSOLNP 
OpenMx version number: 2.7.11 
Need help?  See help(mxSummary) 

|

Comparing the results of Mplus(Estimator=ML),

(attachment is the full output of Mplus).
Fit valus are same,
path coefficents are same,
BUT, lower bound / upper bound of indirect effects are not perfectly identical.

followings are from OpenMx,

confidence intervals:
                      lbound   estimate      ubound
org.indirect[1,1] -0.2768346 -0.1393785 -0.01158206     
org.indirect[2,1] -0.2580154 -0.1305354 -0.01085831     

and followings are from Mplus.

CONFIDENCE INTERVALS OF MODEL RESULTS
           Lower 2.5%      Estimate       Upper 2.5%  
IND1       -0.269          -0.139          -0.010  
IND2       -0.252          -0.131          -0.009  
 

Can I ignore this difference?

Thank you,
Soyoung.

File attachments: 
jpritikin's picture
Offline
Joined: 05/24/2012 - 00:35
CI details

Use summary(model, verbose=TRUE) to get detailed information about the CIs. I'm not sure if Mplus makes this information available for comparison, but that's not a fault of OpenMx.

Soyoung's picture
Offline
Joined: 05/08/2017 - 11:45
this is another example (2-1-1 model)

Hi,
I would like to ask about the other example of Multilevel SEM.

Now, I am trying 2-1-1 mediation model (using the attachment ex9.5.dat)
I checked the values of -2LnL are almost same; 4303.159 (OpenMx); 4303.158 (Mplus),
but some estimates are not identical.

If you look at the diagram file (attachment), a path and a intercept are different between two programs,
others are almost identical.

Could you please check if there is a way to get the same result?

I want to make this 2-1-1 model based on the article:
Preacher, K. J., Zhang, Z., & Zyphur, M. J. (2011). Alternative methods for assessing mediation in multilevel data: The advantages of multilevel SEM. Structural Equation Modeling, 18, 161-182.
Appendix of this article is here.
 
my 2-1-1 model is based on E. 2-1-1 model (MSEM) in the Appendix above.

Thank you in advance.
Soyoung.

P.S. following code is what I have run in OpenMx.

library(OpenMx)
exdata <- read.table("ex9.5.dat")
colnames(exdata) <- c('Y', 'M', 'V1', 'V2', 'X', 'CLUS')
Data.L1 <- exdata[,c('CLUS','M','Y')]
Data.L2 <- exdata[!duplicated(exdata$CLUS),c('CLUS', 'X')]
upperModel <- mxModel(
  "upper", type="RAM",
  mxData(Data.L2, 'raw', primaryKey = "CLUS"),
  manifestVars = c('X'),
  latentVars = c('Mb','Yb'),
  mxPath(c('X','Mb'), 'Yb', values=1),
  mxPath('X', 'Mb', values=1), 
  mxPath('one', 'X'), # means
  mxPath(c('X','Mb','Yb'), arrows=2, values=1)  )
lowerModel <- mxModel(
  "lower", type="RAM", upperModel,
  mxData(Data.L1, 'raw'),
  manifestVars = c('M', 'Y'),
  mxPath('M', 'Y', values=1),
  mxPath('one', c('M', 'Y')), # means
  mxPath(c('M', 'Y'), arrows=2), # variances
  mxPath('upper.Mb', 'M', values=1, free=FALSE, joinKey = "CLUS"),
  mxPath('upper.Yb', 'Y', values=1, free=FALSE, joinKey = "CLUS") )
lowerModel$expectation$.useSufficientSets <- FALSE
model211 <- mxRun(lowerModel, intervals=TRUE)
summary(model211)

and i got the following result.

free parameters:
           name  matrix row col    Estimate  Std.Error A
1  lower.A[2,1]       A   Y   M  0.33571331 0.02560867  
2  lower.S[1,1]       S   M   M  5.39126176 0.37565626  
3  lower.S[2,2]       S   Y   Y  1.45223197 0.10133471  
4  lower.M[1,1]       M   1   M  0.70880568 0.20585133  
5  lower.M[1,2]       M   1   Y -0.25604332 0.09926723  
6  upper.A[2,1] upper.A  Mb   X  1.59487423 0.22672706  
7  upper.A[3,1] upper.A  Yb   X  0.82937245 0.16255263  
8  upper.A[3,2] upper.A  Yb  Mb -0.08507563 0.08169814  
9  upper.S[1,1] upper.S   X   X  0.81028050 0.12078958  
10 upper.S[2,2] upper.S  Mb  Mb  2.53495826 0.56401216  
11 upper.S[3,3] upper.S  Yb  Yb  0.50668743 0.12456607  
12 upper.M[1,1] upper.M   1   X -0.13019981 0.09488538  
Model Statistics: 
               |  Parameters  |  Degrees of Freedom  |  Fit (-2lnL units)
       Model:             12                   1078              4303.159
   Saturated:             NA                     NA                    NA
Independence:             NA                     NA                    NA
Number of observations/statistics: 590/1090
Information Criteria: 
      |  df Penalty  |  Parameters Penalty  |  Sample-Size Adjusted
AIC:       2147.159               4327.159                       NA
BIC:      -2574.613               4379.720                 4341.624
CFI: NA 
TLI: 1   (also known as NNFI) 
RMSEA:  0  [95% CI (NA, NA)]
Prob(RMSEA <= 0.05): NA
To get additional fit indices, see help(mxRefModels)
timestamp: 2017-06-12 15:34:32 
Wall clock time (HH:MM:SS.hh): 00:00:01.09 
optimizer:  CSOLNP 
OpenMx version number: 2.7.11 
Need help?  See help(mxSummary)
AdminNeale's picture
Offline
Joined: 03/01/2013 - 14:09
Thank you for clarity

First, thank you for providing an example that runs nicely, and for clearly documenting the differences! This helps tremendously. Second, I have been able to reproduce the OpenMx side of the analysis and recover your estimates. The Mplus output is a bit brief and it doesn't seem to contain expected means or covariances - though these would seem likely to equal those of OpenMx because the model fit is the same. Sometimes expectations for the model implied means differ, and this may be the case here. I note the warning in the Mplus output:

*** WARNING in MODEL command
In the MODEL command, the following variable is a y-variable on the BETWEEN
level and an x-variable on the WITHIN level. This variable will be treated
as a y-variable on both levels: M

and I'm wondering if some of the components of the expected mean are omitted because of this 'treatment.' The difference between the two expected means of Y under OpenMx vs Mplus is .251*1.595*-.13-(-.085*1.595*-.13) = -0.0696696, which is the same as the expected mean 'due' to the path from X to Mb to M to Y, .336*(1.595*-.13) = -0.0696696. So it seems that the 'treated as y-variable' translates into 'the expected mean of Y lacks the effect of the mean of X transmitted through Mb to M to Y, but includes the effect of the mean of X transmitted through Yb.' Just how reasonable this is I will defer to others - maybe the Muthens will chime in here ;).

Hope this helps!

Soyoung's picture
Offline
Joined: 05/08/2017 - 11:45
WARNING in Mplus is not error

Hi, Mike.
Thank you so much for your comments. :)
I have emailed the Mplus support team, and I have got the mail that the WARNING in Mplus is not error message. I think we can ignore the warning.

I have one more question.
I would like to ask about the way to handle level-1 variables.
Does OpenMx do latent variable decomposition?
If I want to see the latent variable decomposition results using OpenMx, what should I do?
 
If you look at the diagram file (attachment) variable M and Y are decomposed into two parts.
Mplus does a latent variable decomposition into a Within and Between part of M and Y.

There is some explanation about latent variable decomposition:
Please click here to view page 274 of Mplus User's Guide. 

I would be happy, if I can find the same solution with two programs.

gbufton's picture
Offline
Joined: 03/23/2019 - 06:26
Computing Fit Statistics for Multilevel Path Analysis

Hello,

This forum has been immensely helpful for me to set up my code properly for a multilevel path analysis I'm working on for my dissertation. I am running into issues when I compute my fit statistics with the following code:

library(OpenMx)
library(tidyverse)
library(magrittr)
library(RCurl)
 
x <- getURL("https://raw.githubusercontent.com/gbufton/dissertation/master/path_data.csv")
data <- read.csv(text = x)
 
data <- read.csv("path_data.csv")
names(data)
 
manifests2 <- c( "Tech_Rel_Avg","Tech_LoA_Avg2","Tech_PEU_Avg", "Tech_PU_Avg","Tech_Fam_Avg", "Tech_Skill_Avg",
                "AffTech", "Age", "TechSE")
 
all_tech <- c("Tech_Rel_Avg","Tech_LoA_Avg2","Tech_PEU_Avg", "Tech_PU_Avg","Tech_Fam_Avg", "Tech_Skill_Avg")
 
Person_data <- data %>%
  mutate(Tech_LoA_Avg2 = 10 - Tech_LoA_Avg) %>%
  select(ID,manifests2) %>%
  distinct(ID,.keep_all = TRUE) %>%
  ungroup()
 
personModel <- mxModel(
  model = "person", type="RAM",
  mxData(Person_data, 'raw', primaryKey = "ID"),
  manifestVars = manifests2,
  latentVars = c('Agency_IO_PEAvg', 'Agency_Control_PEAvg'), 
  mxPath('AffTech', c('Tech_PEU_Avg', 'Tech_PU_Avg', 'Tech_Rel_Avg'), values = .10),
  mxPath('Age', c('Tech_PEU_Avg', 'Tech_PU_Avg', 'Tech_Rel_Avg'), values = .10),
  mxPath('TechSE', c('Tech_PEU_Avg', 'Tech_PU_Avg', 'Tech_Fam_Avg', 'Tech_Skill_Avg'), values = .10),
  mxPath(c(all_tech, 'Agency_IO_PEAvg'), "Agency_Control_PEAvg"), 
  mxPath(all_tech,"Agency_IO_PEAvg"), 
  mxPath('one',manifests2), # means
  mxPath(c(manifests2,"Agency_IO_PEAvg",'Agency_Control_PEAvg'), arrows=2)) # variances
 
manifests1 <- c("PE_Agency_IO", "PE_Agency_Control", "Valence", "Engagement", "Authenticity")
 
Within_data <- data %>%
  select(ID, manifests1)
 
withinModel <- mxModel(
  model = "within", type="RAM", personModel,
  mxData(Within_data, 'raw'),
  manifestVars = manifests1, 
  mxPath("PE_Agency_IO","PE_Agency_Control", values = .5), ## m --> y
  mxPath("PE_Agency_Control", c("Valence", "Engagement", "Authenticity"), values = .25),
  mxPath("one", manifests1), # means
  mxPath(manifests1, arrows = 2), # variances
  mxPath('person.Agency_IO_PEAvg', 'PE_Agency_IO', values=1, free=FALSE, joinKey = "ID"),
  mxPath('person.Agency_Control_PEAvg', 'PE_Agency_Control', values=1, free=FALSE, joinKey = "ID")  
)
 
 
model <- mxRun(withinModel, intervals=TRUE)
OpenMx:::summary.MxModel(model, refModels = mxRefModels(model, run = TRUE))

The error message I get is the following:

Warning message:
In computeFitStatistics(likelihood, DoF, chi, chiDoF, retval[["numObs"]],  :
Your model may be mis-specified (and fit worse than an independence model), or you may be
using the wrong independence model, see ?mxRefModels

I know this is a very complicated model, although I get the same error message when I break it down to only three variables. I am wondering if I am using the wrong reference model given that it's a two-level model or something that could be a similar easy fix. Thank you in advance for any help you can offer.

Best,
Gina

AdminRobK's picture
Offline
Joined: 01/24/2014 - 12:15
compatibility with multilevel

I have a feeling that mxRefModels() wasn't designed to be used with multilevel models.

AnaUrmeneta's picture
Offline
Joined: 10/01/2019 - 17:39
Can I include latent variables in higher levels? Model 2-2-1

Hello, I am new in this forum. I am very interested in OpenMx capabilities to analyze multilevel SEM models.

I have read this node to see examples of how I could implement my own model. The problem is that I don't know if the cases presented here would be appropriate for my study.

I want to evaluate a 2-2-1 mediation model. The independent variable X is "leadership", a latent variable with 16 indicators, grouped into 3 factors (second-order factorial model). The mediating variable M is "job satisfaction", a latent variable measured by 20 indicators. Both variables X and M were measured in a sample of 80 middle-level workers. Finally, the response variable Y is "customer satisfaction", latent variable with 6 indicators, measured in a large sample of approx. 9000 customers.
These data are nested. These are customers and workers of 8 hotels. Clients evaluate the hotel in general and five departments, from where we extract the sample of 80 middle-level workers.
The examples of multilevel mediation that I have seen so far have no latent variables in the upper levels. In fact, I have seen few examples 2-2-1. I don't know if I can build a model with these characteristics with OpenMx. If someone could show me an example, I would greatly appreciate it! :-)

I also wanted to ask if it would be better in my case, first to analyze the responses of the workers separately. Build a model with two levels for this group. In the first one, I would include the latent variables and in the second the variance explained by the structure of the data (workers nested in departments and hotels) would be collected. Once adjusted this model I could treat it as a submodel, as in the example of Pritikin (message 10 of this node), and join it to a new model with customers. Another option would be to calculate the factor scores of the latent variables X and M, and treat them as variables observed in a 2-2-1 model, in which the only latent variable would be in the first level and correspond to the customers. In the second level would be the observed variables X and M, and the grouping structure.

Please, I would be very grateful if you could guide me, since I have been looking for models similar to what I want to study for several weeks and I have not been lucky.