Login

News

five-factor model

Posted on
No user picture. loveyu3317 Joined: 07/25/2024
Forums
OpenMx General Help
Hi,

I'm using OpenMx to run a five-factor CFA model for 14 variables. "phenfile0" is the dataset and I used its 3rd to 16th columns as a input with the colnames of "y1" to "y14".
"y12" "y13" "y14" belong to "F1"; "y9" "y10" "y11" belong to "F2"; "y6" "y7" "y8" belong to F3; "y3" "y4" "y5" belong to F4; and "y1" "y2" belong to F5.
**Issue**: the output gave a larger SE for some loadings. I'm not sure if the code is correct. Thanks for your help!


phenofacmodel <- mxModel(
  name="PhenFactor",
  
  # data 
  mxData(observed=phenfile0[,3:16],type="raw"),
  
  mxMatrix(type="Full",nrow=1,ncol=14,free=T,values=colMeans(phenfile0[,3:16],na.rm=T),
           labels=c("m1","m2","m3","m4","m5","m6","m7","m8","m9","m10","m11","m12","m13","m14"),name="Mu",
           dimnames=list(NULL,c('y1', 'y2', 'y3', 'y4', 'y5', 'y6', 'y7', 'y8', 'y9', 'y10', 'y11', 'y12', 'y13', 'y14'))),
  
  mxMatrix(type="Full",nrow=14,ncol=1,name="Lambda1",
           free=  c(F,F,F,F,F,F,F,F,F,F,F,T,T,T),
           values=c(0,0,0,0,0,0,0,0,0,0,0,1,1,1),
           labels=
             c(NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,"l1_12","l1_13","l1_14"),
           dimnames=list(
             c(NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,'y12','y13','y14'),"F1")),
  mxMatrix(type="Full",nrow=14,ncol=1,name="Lambda2",
           free=  c(F,F,F,F,F,F,F,F,T,T,T,F,F,F),
           values=c(0,0,0,0,0,0,0,0,1,1,1,0,0,0),
           labels=
             c(NA,NA,NA,NA,NA,NA,NA,NA,"l2_9","l2_10","l2_11",NA,NA,NA),
           dimnames=list(
             c(NA,NA,NA,NA,NA,NA,NA,NA,'y9','y10','y11',NA,NA,NA),"F2")),
  mxMatrix(type="Full",nrow=14,ncol=1,name="Lambda3",
           free=  c(F,F,F,F,F,T,T,T,F,F,F,F,F,F),
           values=c(0,0,0,0,0,1,1,1,0,0,0,0,0,0),
           labels=
             c(NA,NA,NA,NA,NA,"l3_6","l3_7","l3_8",NA,NA,NA,NA,NA,NA),
           dimnames=list(
             c(NA,NA,NA,NA,NA,'y6','y7','y8',NA,NA,NA,NA,NA,NA),"F3")),
  mxMatrix(type="Full",nrow=14,ncol=1,name="Lambda4",
           free=  c(F,F,T,T,T,F,F,F,F,F,F,F,F,F),
           values=c(0,0,1,1,1,0,0,0,0,0,0,0,0,0),
           labels=
             c(NA,NA,"l4_3","l4_4","l4_5",NA,NA,NA,NA,NA,NA,NA,NA,NA),
           dimnames=list(
             c(NA,NA,'y3','y4','y5',NA,NA,NA,NA,NA,NA,NA,NA,NA),"F4")),
  mxMatrix(type="Full",nrow=14,ncol=1,name="Lambda5",
           free=  c(T,T,F,F,F,F,F,F,F,F,F,F,F,F),
           values=c(1,1,0,0,0,0,0,0,0,0,0,0,0,0),
           labels=
             c("l5_1","l5_2",NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA),
           dimnames=list(
             c('y1','y2',NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA),"F5")),
  
  #The common factor's variance, fixed to 1.0 to identify the scale:
  mxMatrix(type="Symm",nrow=1,free=F,values=1,labels="fvar1",lbound=0,name="Fvar1"),
  mxMatrix(type="Symm",nrow=1,free=F,values=1,labels="fvar2",lbound=0,name="Fvar2"),
  mxMatrix(type="Symm",nrow=1,free=F,values=1,labels="fvar3",lbound=0,name="Fvar3"),
  mxMatrix(type="Symm",nrow=1,free=F,values=1,labels="fvar4",lbound=0,name="Fvar4"),
  mxMatrix(type="Symm",nrow=1,free=F,values=1,labels="fvar5",lbound=0,name="Fvar5"),
  
  #The unique variances:
  mxMatrix(
    type="Diag",nrow=14,free=T,
    labels=c("u1","u2","u3","u4","u5","u6","u7","u8","u9","u10","u11","u12","u13","u14"),lbound=0.0001,name="U",
    values=diag(cov(phenfile0[,3:16],use="pair"))),
  
  #Compute
  mxAlgebra(Lambda1%*%Fvar1%*%t(Lambda1),name="Comm1"),
  mxAlgebra(Lambda2%*%Fvar2%*%t(Lambda2),name="Comm2"),
  mxAlgebra(Lambda3%*%Fvar3%*%t(Lambda3),name="Comm3"),
  mxAlgebra(Lambda4%*%Fvar4%*%t(Lambda4),name="Comm4"),
  mxAlgebra(Lambda5%*%Fvar5%*%t(Lambda5),name="Comm5"),
  
  mxAlgebra(Comm1+Comm2+Comm3+Comm4+Comm5+U,name="Sigmay",dimnames=list(
    c('y1', 'y2', 'y3', 'y4', 'y5', 'y6', 'y7', 'y8', 'y9', 'y10', 'y11', 'y12', 'y13', 'y14'),
    c('y1', 'y2', 'y3', 'y4', 'y5', 'y6', 'y7', 'y8', 'y9', 'y10', 'y11', 'y12', 'y13', 'y14'))),
  mxExpectationNormal(covariance="Sigmay",means="Mu"),
  mxFitFunctionML()
)

phenofacrun <- mxRun(phenofacmodel)
summary(phenofacrun, refModels=mxRefModels(phenofacrun, run = TRUE))
Replied on Mon, 08/26/2024 - 21:52
Picture of user. AdminNeale Joined: 03/01/2013

Hi

It's quite possible for factor loadings to have different standard errors. This will occur if the variables have different variances, or different response probabilities if the data are binary or ordinal. BTW are they continuous, normally distributed variables or something else?

To simplify model specification for those who get path diagrams, the Onyx-sem software is really handy - it will spit out OpenMx code from a diagram.

I note that there are also some possible simplifications of your code, where the lambda vectors are cbind() (column-wise glued together) into a matrix Lambda (or just set up that way initially), and the variances specified as a diagonal matrix FCov- actually an identity matrix since all the variances are 1 and no covariances between factors are specified. Then the covariance algebra Lambda %*% FCov %*% t(Lambda) stays the same however many factors you have vary the number of columns, and you also have the option of letting them correlate.

Log in or register to post comments
Replied on Wed, 08/28/2024 - 14:17
Picture of user. AdminHunter Joined: 03/01/2013

I agree with Mike Neale about the standard errors.

Below is a version of the same model that is just easier (in my opinion) to read and write.


require(OpenMx)
nv <- 14
mv <- paste0('y', 1:nv)

# Make fake data for demonstration
phenfile0 <- as.data.frame(matrix(rnorm((nv+2)*100), nrow=100, ncol=nv+2))
names(phenfile0) <- c(letters[1:2], mv)

phenofacmodel <- mxModel(
name="PhenFactor",
manifestVars=list(exo=mv),
latentVars=list(exo=paste0('F', 1:5)),
type='LISREL',

# Data
mxData(observed=phenfile0[,3:16],type="raw"),

# Means
mxPath('one', mv, labels=paste0('m', 1:nv), values=colMeans(phenfile0[,3:16],na.rm=T)),

# Loadings
mxPath('F1', mv[12:14], labels=paste0('l1_', 12:14), values=1),
mxPath('F2', mv[9:11], labels=paste0('l2_', 9:11), values=1),
mxPath('F3', mv[6:8], labels=paste0('l3_', 6:8), values=1),
mxPath('F4', mv[3:5], labels=paste0('l4_', 3:5), values=1),
mxPath('F5', mv[1:2], labels=paste0('l5_', 1:2), values=1),

# Factor Variances
mxPath(paste0('F', 1:5), arrows=2, free=FALSE, values=1),

# Unique variances
mxPath(mv, arrows=2, labels=paste0('u', 1:14), lbound=0.0001,
values=diag(cov(phenfile0[,3:16],use="pair")))
)

# Loadings
mxEval(LX, phenofacmodel)

# Resid
mxEval(TD, phenofacmodel)

# Intercepts/means
mxEval(TX, phenofacmodel)

# Factor variances
mxEval(PH, phenofacmodel)

phenofacrun <- mxRun(phenofacmodel)
summary(phenofacrun, refModels=mxRefModels(phenofacrun, run = TRUE))

Log in or register to post comments
  Login

News