#### Approach 2 - Estimate ONE BIG (with all effect sizes) Random-Effects Pooled Correlation Matrix. #### Use the BIG Pooled correlation matrices and estimate models 1, 2, and 3. #### Use the select variable approach where correlation matrices are constructed #### from the BIG Pooled correlation matrices #### Adapted from https://openmx.ssri.psu.edu/sites/default/files/select_variables.pdf ### Step 1: Import Libraries library("OpenMx") library("lavaan") library("metaSEM") library("semPlot") ### Step 2: Load Data. ### SAME DATA FILE USED IN APPROACH 1. ### Using the processed file "Data - Final - Positive Definite Removed.csv" data3 <- read.csv("/Users/srikanthparameswaran/Desktop/TSSEM Project 2/Approach 2/Approach 2 Data.csv") ### Step 3: Setup correlation matrices nvar <- 9 varnames <- c("X","Y","A","B","C","M","D","E","F") labels <- list(varnames,varnames) cordat <- list() for (i in 1:nrow(data3)){ cordat[[i]] <- vec2symMat(as.matrix(data3[i,2:37]),diag = FALSE) dimnames(cordat[[i]]) <- labels} data3\$data<-cordat ### Step 4: Put NA on diagonal if variable is missing for (i in 1:length(data3\$data)){ for (j in 1:nrow(data3\$data[[i]])){ if (sum(is.na(data3\$data[[i]][j,]))==nvar-1) {data3\$data[[i]][j,j] <- NA} }} ### Step 5: Put NA on diagonal for variable with least correlations for (i in 1:length(data3\$data)){ for (j in 1:nrow(data3\$data[[i]])){ for (k in 1:nvar){ if (is.na(data3\$data[[i]][j,k])==TRUE &is.na(data3\$data[[i]][j,j])!=TRUE &is.na(data3\$data[[i]][k,k])!=TRUE){ if(sum(is.na(data3\$data[[i]])[j,])>sum(is.na(data3\$data[[i]])[k,])) {data3\$data[[i]][k,k] = NA} if(sum(is.na(data3\$data[[i]])[j,])<=sum(is.na(data3\$data[[i]])[k,])) {data3\$data[[i]][j,j] = NA} }}}} # Display processed data data3\$data # Assess missingness pattern.na(data3\$data, show.na = FALSE) pattern.n(data3\$data, data3\$Corr_Sample) ### Step 6: Check Positive-Definiteness. Everything should be TRUE. is.pd(data3\$data) length(is.pd(data3\$data)) ### Step 7: Random-effects Stage 1 estimation random_stage1<- tssem1(data3\$data, data3\$Corr_Sample, method="REM", RE.type="Diag") summary(random_stage1) vec2symMat(coef(random_stage1, select="fixed"), diag = FALSE) ### Step 8: Creating variables for Random-effects Stage 2 estimation pooled <- vec2symMat(coef(random_stage1, select="fixed"), diag = FALSE) dimnames(pooled)<- list(varnames,varnames) pooled aCov <- vcov(random_stage1, select="fixed") aCov n2 <- random_stage1\$total.n n2 ### Step 9: Model 1 ## Step 9a: Variables to keep ## Let's exclude "Y", "M", and "C" var.to.keep <- c(TRUE, FALSE, TRUE, TRUE, FALSE, FALSE,TRUE, TRUE, TRUE) names(var.to.keep) <- varnames var.to.keep ## Correlation coefficients to keep cor.to.keep <- vechs(outer(var.to.keep, var.to.keep, function(y, z) y&z)) names(cor.to.keep) <- colnames(aCov) cor.to.keep ## New correlation matrix Cov_new <- pooled[var.to.keep, var.to.keep] Cov_new ## New sampling covariance matrix aCov_new <- aCov[cor.to.keep, cor.to.keep] aCov_new ## Step 9b: Model setup in Lavaan. Converting Lavaan to RAM model1 <- 'X ~ A2X*A + B2X*B + D2X*D + E2X*E + F2X*F' plot(model1, col="yellow") RAM1 <- lavaan2RAM(model1, obs.variables = varnames[var.to.keep]) RAM1 RAM1\$S # ## Step 9c: All variances except the DVs were set to 1. RAM1\$S[2,2] <- 1 RAM1\$S[3,3] <- 1 RAM1\$S[4,4] <- 1 RAM1\$S[5,5] <- 1 RAM1\$S[6,6] <- 1 RAM1\$S RAM1 ## Step 9d: Model estimation fit1 <- wls(Cov=Cov_new , aCov=aCov_new, n= n2, RAM=RAM1, diag.constraints = TRUE) summary(fit1) my.plot1 <- meta2semPlot(fit1) semPaths(my.plot1, whatLabels="est", layout = "tree2", sizeMan=12, edge.label.cex=2, color="yellow", edge.color = "black", weighted=FALSE) (Smatrix <- diag(mxEval(Smatrix, fit1\$mx.fit))) 1 - Smatrix ### Step 10: Model 2 ## Step 10a: Variables to keep ## Let's exclude "Y" var.to.keep <- c(TRUE, FALSE, TRUE, TRUE, TRUE, TRUE,TRUE, TRUE, TRUE) names(var.to.keep) <- varnames var.to.keep ## Correlation coefficients to keep cor.to.keep <- vechs(outer(var.to.keep, var.to.keep, function(y, z) y&z)) names(cor.to.keep) <- colnames(aCov) cor.to.keep ## New correlation matrix Cov_new1 <- pooled[var.to.keep, var.to.keep] Cov_new1 ## New sampling covariance matrix aCov_new1 <- aCov[cor.to.keep, cor.to.keep] aCov_new1 ## Step 10b: Model setup in Lavaan. Converting Lavaan to RAM model2 <- 'X ~ A2X*A + M2X*M + B2X*B + C2X*C + D2X*D + E2X*E + F2X*F M ~ A2M*A + B2M*B' plot(model2, col="yellow") RAM2 <- lavaan2RAM(model2, obs.variables = varnames[var.to.keep]) RAM2 RAM2\$S ## Step 10c: All variances except the DVs were set to 1. RAM2\$S[2,2] <- 1 RAM2\$S[3,3] <- 1 RAM2\$S[4,4] <- 1 RAM2\$S[6,6] <- 1 RAM2\$S[7,7] <- 1 RAM2\$S[8,8] <- 1 RAM2\$S RAM2 ## Step 10d: Model estimation fit2 <- wls(Cov=Cov_new1 , aCov=aCov_new1, n= n2, RAM=RAM2, diag.constraints = TRUE) summary(fit2) fit2 <-rerun(fit2) summary(fit2) my.plot2 <- meta2semPlot(fit2) semPaths(my.plot2, whatLabels="est", layout = "tree2", sizeMan=12, edge.label.cex=2, color="yellow", edge.color = "black", weighted=FALSE) (Smatrix <- diag(mxEval(Smatrix, fit2\$mx.fit))) 1 - Smatrix ### Step 11: Model 3 ## Step 11a: Variables to keep - NOT APPLICABLE ## Step 11b: Model setup in Lavaan. Converting Lavaan to RAM model3 <- 'X ~ A2X*A + M2X*M + B2X*B + C2X*C + D2X*D + E2X*E + F2X*F M ~ A2M*A + B2M*B Y ~ A2Y*A + M2Y*M + B2Y*B + C2Y*C + D2Y*D + E2Y*E + F2Y*F + X2Y*X' plot(model3, col="yellow") RAM3 <- lavaan2RAM(model3, obs.variables=c("X","Y","A","B","C","M","D","E","F")) RAM3 ## Step 11c: All variances except the DVs were set to 1. RAM3\$S RAM3\$S[3,3] <- 1 RAM3\$S[4,4] <- 1 RAM3\$S[5,5] <- 1 RAM3\$S[7,7] <- 1 RAM3\$S[8,8] <- 1 RAM3\$S[9,9] <- 1 RAM3\$S RAM3 ## Step 11d: Model estimation fit3 <- wls(Cov=pooled, aCov=aCov, n= n2, RAM=RAM3, diag.constraints = TRUE) summary(fit3) fit3 <-rerun(fit3) summary(fit3) my.plot3 <- meta2semPlot(fit3) semPaths(my.plot3, whatLabels="est", layout = "tree2", sizeMan=12, edge.label.cex=2, color="yellow", edge.color = "black", weighted=FALSE) (Smatrix <- diag(mxEval(Smatrix, fit3\$mx.fit))) 1 - Smatrix