Login

News

The goal of moderator analysis in metaSEM::wls()

Posted on
No user picture. sharmel Joined: 07/31/2023
Forums
metaSEM
Dear Mike and Colleagues,

My understanding is that moderator analysis at the 2nd stage of metaSEM approach using `metaSEM::wls()` involves multi-group analysis.

Thus, as a user I have to `subset` my data by categories of my moderator, and fit the **same model syntax** to each subset as shown below.

But at the end, how can I conclude whether I have a significant moderator or not? For example, should/can I do: `anova(wls1, wls2)` and if one model is significantly different from another conclude that the moderator is significant? (I appreciate a demonstration.)


library(metaSEM)

## Sample sizes
n1 <- 100
n2 <- 200

## Variable labels
vars <- c("y", "x1", "x2")

## Group 1 data and model
R1 <- matrix(c(1.00, 0.22, 0.24,
0.22, 1.00, 0.30,
0.24, 0.30, 1.00), ncol=3, nrow=3,
dimnames=list(vars, vars))
acov1 <- asyCov(R1, n1)

model1 <- "y ~ b1a*x1 + b2a*x2
## Variances of x1 and x2 are 1
x1 ~~ 1*x1
x2 ~~ 1*x2
## x1 and x2 are correlated
x1 ~~ ra*x2"

RAM1 <- lavaan2RAM(model1, obs.variables=vars)

wls1 <- wls(model="Group1", Cov=R1, aCov=acov1, n=n1, RAM=RAM1)

## Group 2 data and model
R2 <- matrix(c(1.00, 0.33, 0.41,
0.33, 1.00, 0.35,
0.41, 0.35, 1.00), ncol=3, nrow=3,
dimnames=list(vars, vars))
acov2 <- asyCov(R2, n2)

model2 <- "y ~ b1b*x1 + b2b*x2
## Variances of x1 and x2 are 1
x1 ~~ 1*x1
x2 ~~ 1*x2
## x1 and x2 are correlated
x1 ~~ rb*x2"

RAM2 <- lavaan2RAM(model2, obs.variables=vars)

wls2 <- wls(model="Group2", Cov=R2, aCov=acov2, n=n2, RAM=RAM2)

Replied on Wed, 10/11/2023 - 20:48
No user picture. sharmel Joined: 07/31/2023

Dear Mike,

Sure, but I asked this to clarify your comment [**HERE**](https://openmx.ssri.psu.edu/comment/9874#comment-9874). In [that comment](https://openmx.ssri.psu.edu/comment/9874#comment-9874), what is the purpose behind creating a combined model like: `mxModel(model="combined", wls1, wls2, mxFitFunctionMultigroup(c("Group1", "Group2")))`? What do the statistics of this combined model tell us about the moderatory role of `"Group"`?

In fact, as I ask in this thread, if the path models are exactly the same across two subsets of data, and we fit `wls(...,run=F)` for each model, then is it possible to compare the two models using `mxModel(..., mxFitFunctionMultigroup(...))` and then `mxRun()`? I really appreciate a demonstration.

Log in or register to post comments
Replied on Wed, 10/11/2023 - 22:06
Picture of user. Mike Cheung Joined: 10/08/2009

The purpose is to conduct a subgroup analysis in meta-analysis, which is equivalent to a multiple-group analysis in SEM.

The models with equality constraints and without equality constraints are nested. Thus, a chi-square difference statistic can be used to compare them.

There are plenty of resources for subgroup analysis in meta-analysis and multiple-group SEM. Here are some materials for a three-day workshop in MASEM. https://github.com/mikewlcheung/masemWorkshop2023

Log in or register to post comments
Replied on Sun, 10/22/2023 - 17:33
No user picture. sharmel Joined: 07/31/2023

Dear Mike,

Is it possible to conduct *pairwise comparisons* between the elements of A.matrix and S.matrix from the two groups in the object `output` below?


## Sample sizes
n1 <- 100
n2 <- 200

## Variable labels
vars <- c("y", "x1", "x2")

## Group 1 data and model
R1 <- matrix(c(1.00, 0.22, 0.24,
0.22, 1.00, 0.30,
0.24, 0.30, 1.00), ncol=3, nrow=3,
dimnames=list(vars, vars))
acov1 <- asyCov(R1, n1)

model1 <- "y ~ b1a*x1 + b2a*x2
## Variances of x1 and x2 are 1
x1 ~~ 1*x1
x2 ~~ 1*x2
## x1 and x2 are correlated
x1 ~~ ra*x2"

## Group 2 data but model is the same as model1
R2 <- matrix(c(1.00, 0.33, 0.41,
0.33, 1.00, 0.35,
0.41, 0.35, 1.00), ncol=3, nrow=3,
dimnames=list(vars, vars))
acov2 <- asyCov(R2, n2)

model2 <- "y ~ b1b*x1 + b2b*x2
## Variances of x1 and x2 are 1
x1 ~~ 1*x1
x2 ~~ 1*x2
## x1 and x2 are correlated
x1 ~~ rb*x2"

RAM1 <- lavaan2RAM(model1, obs.variables=vars)
RAM2 <- lavaan2RAM(model2, obs.variables=vars)

wls1 <- wls(model="Group1", Cov=R1, aCov=acov1, n=n1, RAM=RAM1, run=F)
wls2 <- wls(model="Group2", Cov=R2, aCov=acov2, n=n2, RAM=RAM2, run=F)

## Combine both groups
wls.model <- mxModel(model="combined", wls1, wls2, mxFitFunctionMultigroup(c("Group1", "Group2")))

wls.fit <- mxRun(wls.model, intervals=TRUE)

zz = summary(wls.fit)

(output = cbind(zz$parameters[-c(3:4,7:10)], zz$CI[c("lbound","ubound")]))

name matrix Estimate Std.Error lbound ubound
1 b1a Group1.Amatrix 0.1626374 0.09928458 -0.034465926 0.3598924
2 b2a Group1.Amatrix 0.1912088 0.09881592 -0.004306087 0.3874504
3 ra Group1.S1 0.3000000 0.09100000 0.120857208 0.4785835
4 b1b Group2.Amatrix 0.2125356 0.06601287 0.082020597 0.3433026
5 b2b Group2.Amatrix 0.3356125 0.06395664 0.209506797 0.4626121
6 rb Group2.S1 0.3500000 0.06204862 0.228300669 0.4720198

Log in or register to post comments
Replied on Mon, 10/23/2023 - 04:40
Picture of user. Mike Cheung Joined: 10/08/2009

Yes, but it is not automatic.
You may test the constraints one by one, e.g., b1a=b1b, b2a=b2b. Please note that it does not adjust the issue of multiple comparisons.
Log in or register to post comments
Replied on Tue, 10/24/2023 - 05:53
Picture of user. Mike Cheung Joined: 10/08/2009

Here is an example.


model1 <- "y ~ b1*x1 + b2a*x2
## Variances of x1 and x2 are 1
x1 ~~ 1*x1
x2 ~~ 1*x2
## x1 and x2 are correlated
x1 ~~ ra*x2"

model2 <- "y ~ b1*x1 + b2b*x2
## Variances of x1 and x2 are 1
x1 ~~ 1*x1
x2 ~~ 1*x2
## x1 and x2 are correlated
x1 ~~ rb*x2"

RAM1 <- lavaan2RAM(model1, obs.variables=vars)
RAM2 <- lavaan2RAM(model2, obs.variables=vars)

wls1 <- wls(model="Group1", Cov=R1, aCov=acov1, n=n1, RAM=RAM1, run=F)
wls2 <- wls(model="Group2", Cov=R2, aCov=acov2, n=n2, RAM=RAM2, run=F)

## Combine both groups
wls.b1 <- mxModel(model="b1", wls1, wls2, mxFitFunctionMultigroup(c("Group1", "Group2")))
fit.b1 <- mxRun(wls.b1, intervals=TRUE)

summary(fit.b1)

anova(wls.fit, fit.b1)

Log in or register to post comments
Replied on Tue, 10/24/2023 - 13:03
No user picture. sharmel Joined: 07/31/2023

Truly appreciated, Mike. Just curious, given that estimates from each group use independent pieces of data, can't we compute the standard error of their coefficients' differences e.g., SE(b1a - b1b = 0) by doing `SE_dif = SQRT(SE_b1a^2 + SE_b1b^2)` assuming near-normally distributed sampling distribution for b1a - b1b?

If yes, then H0: b1a - b1b = 0, may be tested with standard wald-type, Z-test.

Log in or register to post comments
Replied on Tue, 10/24/2023 - 22:53
No user picture. sharmel Joined: 07/31/2023

Final confirmation, can Wald tests be applied both to the elements of A.matrix and S.matrix from the two groups?

Also, are there any references (papers, book chapters etc.) discussing the use of Wald tests in the manner that I suggested in my previous post for comparing across the estimates from independent groups?

Log in or register to post comments
Replied on Thu, 11/23/2023 - 11:14
No user picture. sharmel Joined: 07/31/2023

Dear Mike,

Is there any way to convert a `wls()` object with `run=TRUE` or `run=FALSE` to a corresponding lavaan object?

Log in or register to post comments
Replied on Wed, 10/09/2024 - 12:55
No user picture. JuanJMV Joined: 07/20/2016

Hi,

 

I am new with MetaSEM. I am trying to fit a simple meta-analysis (I just want to analyse 2 correlations; see attached file).

 

Here you can find the code that I am using with some simulated data but if there are no regressions I cannot make it work. It should be an easy model but I am stuck. I just want two covariances and I have fixed all the variances to one. 

Correlations between S and D, S and I and D and I are stored here:

rSD <- c(0.85, 0.8, 0.90, 0.78, 0.70,NA,NA,NA,NA,NA)
rSI <- c(NA,NA,NA,NA,NA,0.85,NA,NA,NA,NA)
rDI <- c(NA,NA,NA,NA,NA,0.60, 0.50, 0.55, 0.60, 0.55)

 

rSI could be al NAs since I am not interested in that correlation. 

Any help and guidance would be really appreciated.

 

library(metaSEM)
source("http://www.suzannejak.nl/MASEM_functions.R")

rSD <- c(0.85, 0.8, 0.90, 0.78, 0.70,NA,NA,NA,NA,NA)
rSI <- c(NA,NA,NA,NA,NA,0.85,NA,NA,NA,NA)
rDI <- c(NA,NA,NA,NA,NA,0.60, 0.50, 0.55, 0.60, 0.55)

 


N <- c(10000,40000, 30000, 10000, 10000,20000,10000, 50000, 10000, 10000  )

data <- as.data.frame(cbind(rSD,rSI,rSI, N))


nvar <- 3
varnames <- c("S","D", "I")
labels <- list(varnames,varnames)

cormatrices <- readstack(data[,c(2,1,3)], no.var = nvar, var.names = varnames, diag = FALSE)

n <-  data$N

pattern.na(cormatrices, show.na=F)
pattern.n(cormatrices, n=n)


my.df <- Cor2DataFrame(cormatrices, n, acov = "weighted")


## Specify model using lavaan syntax
model <-
 '

# Covariances
D ~~ I
D ~~ S
# Variances

D ~~ 1*D
S ~~ 1*S
I ~~ 1*I
'

RAM1 <- lavaan2RAM(model, obs.variables=varnames)
RAM1

 


## Create the model implied correlation structure with implicit diagonal constraints
M0 <- create.vechsR(A0=RAM1$A, S0=RAM1$S)

## Create the heterogeneity variance-covariance matrix
T0 <- create.Tau2(RAM=RAM1, RE.type="Diag", Transform="expLog", RE.startvalues=0.05)

## Fit the model
mx.fit0 <- osmasem(model.name="No moderator", Mmatrix=M0, Tmatrix=T0, data=my.df)


## View the results
summary(mx.fit0, fitIndices = TRUE)
VarCorr(mx.fit0)

 

 

 

File attachments
Model_2.pdf
Log in or register to post comments
  Login

News