You are here

Standardizing factor loadings

8 posts / 0 new
Last post
szellers's picture
Offline
Joined: 04/04/2018 - 14:52
Standardizing factor loadings

I am working with a common path model and I am having difficulty standardizing my factor loadings. In the past, I have identified the model by constraining the variance of the latent phenotype to 1, then I standardize factor loadings by using a matrix of standard deviations (SDs on diagonal, 0s on off-diagonal) and multiplying that by the matrix of the unstandardized factor loadings (I can attach code for this if necessary)

In the model I am currently working with, I have identified the model by fixing the first factor loading to 1 and I am finding that the method of standardizing factor loadings I've used before doesn't seem to be working properly (I get standardized loadings greater than 1). What method should I be using to standardize loadings when the first loading is fixed to 1?

Thanks!

Stephanie

AdminRobK's picture
Offline
Joined: 01/24/2014 - 12:15
standardized loadings

One way or another, you need to multiply each loading by the standard deviation of the common factor, and divide it by the standard deviation of the corresponding observable variable. It's analogous to how you'd standardize a linear regression coefficient.

Previously, you were probably premultiplying the column vector of loadings by a diagonal matrix the diagonal elements of which were the reciprocals of the observable variables' standard deviations. Try Kronecker-multiplying the column of loadings by the latent factor's standard deviation, and then premultiply the resulting rescaled column vector by the same diagonal matrix as before. You can use parentheses to control order-of-operations.

I don't think I can be more specific without seeing the script you're working from.

szellers's picture
Offline
Joined: 04/04/2018 - 14:52
Thank you! I think I was not

Thank you! I think I was not considering the standard deviation of the common factor, as that would have just been 1 in previous models when the variance of the factor was constrained to 1, but that is not the case in this model. Below is my code, I am trying to standardize flBDMN and flBDCO.

###################################################################
#base models
###################################################################
nl <- 3 # number of latent factors
# Matrices ac, cc, and ec to store a, c, and e path coefficients for latent phenotype(s)
XMN <- mxMatrix(type="Lower", nrow=nl, ncol=nl, free=TRUE,
                values=c(.6, 
                         .1, .1, 
                         .1, .1, 0),
                labels=c("x11MN",
                         "x21MN", "x22MN", 
                         "x31MN", "x32MN", "x33MN"),
                lbound=.00001, name="XMN" )
 
YMN <- mxMatrix(type="Lower", nrow=nl, ncol=nl, free=TRUE,
                values=c(.6,
                         .1, .1, 
                         .1, .1, .1),
                labels=c("y11MN",
                         "y21MN", "y22MN",
                         "y31MN", "y32MN", "y33MN"),
                lbound=.00001, name="YMN" )
 
ZMN <- mxMatrix(type="Lower", nrow=nl, ncol=nl, free=TRUE,
                values=c(.6,
                         .1, .1,
                         .1, .1, .1),
                labels=c("z11MN",
                         "z21MN", "z22MN",
                         "z31MN", "z32MN", "z33MN"),
                lbound=.00001, name="ZMN" )
 
XCO <- mxMatrix(type="Lower", nrow=nl, ncol=nl, free=TRUE,
                values=c(.6, 
                         .1, .1, 
                         .1, .1, 0),
                labels=c("x11CO",
                         "x21CO", "x22CO", 
                         "x31CO", "x32CO", "x33CO"),
                lbound=.00001, name="XCO" )
 
YCO <- mxMatrix(type="Lower", nrow=nl, ncol=nl, free=TRUE,
                values=c(.6,
                         .1, .1, 
                         .1, .1, .1),
                labels=c("y11CO",
                         "y21CO", "y22CO",
                         "y31CO", "y32CO", "y33CO"),
                lbound=.00001, name="YCO" )
 
ZCO <- mxMatrix(type="Lower", nrow=nl, ncol=nl, free=TRUE,
                values=c(.6,
                         .1, .1,
                         .1, .1, .1),
                labels=c("z11CO",
                         "z21CO", "z22CO",
                         "z31CO", "z32CO", "z33CO"),
                lbound=.00001, name="ZCO" )
 
 
AlMN <- mxAlgebra(XMN %*% t(XMN), name="AlMN")
ClMN <- mxAlgebra(YMN %*% t(YMN), name="ClMN")
ElMN <- mxAlgebra(ZMN %*% t(ZMN), name="ElMN")
 
AlCO <- mxAlgebra(XCO %*% t(XCO), name="AlCO")
ClCO <- mxAlgebra(YCO %*% t(YCO), name="ClCO")
ElCO <- mxAlgebra(ZCO %*% t(ZCO), name="ElCO")
 
# Matrices as, cs, and es to store a, c, and e path coefficients for specific factors
pathAsMN <- mxMatrix(type="Diag", nrow=nvMN, ncol=nvMN, free=TRUE, values=.3,
                     labels=c("a11MN", "a22MN", "a33MN", "a44MN", "a55MN", "a66MN", "a77MN", "a88MN"),
                     lbound=.00001, name="asMN" )
pathCsMN <- mxMatrix(type="Diag", nrow=nvMN, ncol=nvMN, free=TRUE, values=.3,
                     labels=c("c11MN","c22MN", "c33MN", "c44MN", "c55MN", "c66MN", "c77MN", "c88MN"),
                     lbound=.00001, name="csMN" )
pathEsMN <- mxMatrix(type="Diag", nrow=nvMN, ncol=nvMN, free=TRUE, values=.3,
                     labels=c("e11MN", "e22MN", "e33MN", "e44MN", "e55MN", "e66MN", "e77MN", "e88MN"),
                     lbound=.00001, name="esMN" )
 
pathAsCO <- mxMatrix(type="Diag", nrow=nvCO, ncol=nvCO, free=TRUE, values=.3,
                     labels=c("a11CO", "a22CO", "a33CO", "a44CO", "a55CO", "a66CO"),
                     lbound=.00001, name="asCO" )
pathCsCO <- mxMatrix(type="Diag", nrow=nvCO, ncol=nvCO, free=TRUE, values=.3,
                     labels=c("c11CO","c22CO", "c33CO", "c44CO", "c55CO", "c66CO"),
                     lbound=.00001, name="csCO" )
pathEsCO <- mxMatrix(type="Diag", nrow=nvCO, ncol=nvCO, free=TRUE, values=.3,
                     labels=c("e11CO", "e22CO", "e33CO", "e44CO", "e55CO", "e66CO"),
                     lbound=.00001, name="esCO" )
 
 
 
 
pathFlGrowthMN <- mxMatrix(type="Full", nrow=5, ncol=nl, free=FALSE, 
                           labels=c("data.zero", "data.intercept", "data.age14twinpair",
                                    "data.zero", "data.intercept", "data.age17twinpair",
                                    "data.zero", "data.intercept", "data.age20twinpair",
                                    "data.zero", "data.intercept", "data.age24twinpair",
                                    "data.zero", "data.intercept", "data.age29twinpair"), name="flGrowthMN", byrow=T)
 
pathFlGrowthCO <- mxMatrix(type="Full", nrow=3, ncol=nl, free=FALSE, 
                           labels=c("data.zero", "data.intercept", "data.age1twinpairC",
                                    "data.zero", "data.intercept", "data.age2twinpairC",
                                    "data.zero", "data.intercept", "data.age3twinpairC"), name="flGrowthCO", byrow=T)
 
 
flFree   <- c(FALSE, FALSE, FALSE,
              TRUE, FALSE, FALSE,
              TRUE, FALSE, FALSE)
 
flValues <- c(1, 0, 0,
              .6, 0, 0,
              .6, 0, 0)
 
#set first loading to 1 and make it fixed for BD to manifest variables
pathFlBDMN <- mxMatrix(type="Full", nrow=3, ncol=nl, free=flFree, 
                       values=flValues, lbound=.000001, labels=c("fBD1MN", NA, NA,
                                                                 "fBD2MN", NA, NA,
                                                                 "fBD3MN", NA, NA), name="flBDMN", byrow=T)
 
pathFlBDCO <- mxMatrix(type="Full", nrow=3, ncol=nl, free=flFree, 
                       values=flValues, lbound=.000001, labels=c("fBD1CO", NA, NA,
                                                                 "fBD2CO", NA, NA,
                                                                 "fBD3CO", NA, NA), name="flBDCO", byrow=T)
 
#combine factor loading matrices
pathflMN <- mxAlgebra(rbind(flBDMN, flGrowthMN), name="flMN")
 
pathflCO <- mxAlgebra(rbind(flBDCO, flGrowthCO), name="flCO")
 
 
 
 
covAMN <- mxAlgebra( expression=flMN %*% AlMN %*% t(flMN) + asMN %*% t(asMN), name="AMN" )
covCMN <- mxAlgebra( expression=flMN %*% ClMN %*% t(flMN) + csMN %*% t(csMN), name="CMN" )
covEMN <- mxAlgebra( expression=flMN %*% ElMN %*% t(flMN) + esMN %*% t(esMN), name="EMN" )
 
covACO <- mxAlgebra( expression=flCO %*% AlCO %*% t(flCO) + asCO %*% t(asCO), name="ACO" )
covCCO <- mxAlgebra( expression=flCO %*% ClCO %*% t(flCO) + csCO %*% t(csCO), name="CCO" )
covECO <- mxAlgebra( expression=flCO %*% ElCO %*% t(flCO) + esCO %*% t(esCO), name="ECO" )
 
matIMN <- mxMatrix( type="Iden", nrow=3, ncol=3, name="IMN")
covPMN <- mxAlgebra( expression=AMN[1:3, 1:3]+CMN[1:3, 1:3]+EMN[1:3, 1:3], name="VMN" )
invSDMN <- mxAlgebra( expression=solve(sqrt(IMN*VMN)), name="iSDMN")
 
matICO <- mxMatrix( type="Iden", nrow=3, ncol=3, name="ICO")
covPCO <- mxAlgebra( expression=ACO[1:3, 1:3]+CCO[1:3, 1:3]+ECO[1:3, 1:3], name="VCO" )
invSDCO <- mxAlgebra( expression=solve(sqrt(ICO*VCO)), name="iSDCO")
 
meanGMN <- mxMatrix( type="Full", nrow=1, ncol=nvMN, free=TRUE, values=.3, 
                     labels=c('m1MN', 'm2MN', 'm3MN', 'm4MN', 'm5MN', 'm6MN', 'm7MN', 'm8MN'), name="expMeanMN" )
regCoefMN <- mxMatrix(type="Full", nrow=1, ncol=nvMN, free=TRUE, values=0, 
                      labels=c("beta1MN", "beta2MN", "beta3MN", "beta4MN", "beta5MN", "beta6MN", "beta7MN", "beta8MN"), name="betaMN")
sexCovMN <- mxMatrix(type="Full", nrow=1, ncol=nvMN, free=FALSE, labels="data.sex", name="sexMN")
MuMN <- mxAlgebra(expression= expMeanMN+(betaMN*sexMN), name="MuMN")
meanGMN2 <- mxAlgebra(cbind(MuMN, MuMN), name="expMeanMN2")
 
meanGCO <- mxMatrix( type="Full", nrow=1, ncol=nvCO, free=TRUE, values=.3, 
                     labels=c('m1CO', 'm2CO', 'm3CO', 'm4CO', 'm5CO', 'm6CO'), name="expMeanCO" )
regCoefCO <- mxMatrix(type="Full", nrow=1, ncol=nvCO, free=TRUE, values=0, 
                      labels=c("beta1CO", "beta2CO", "beta3CO", "beta4CO", "beta5CO", "beta6CO"), name="betaCO")
sexCovCO <- mxMatrix(type="Full", nrow=1, ncol=nvCO, free=FALSE, labels="data.sex", name="sexCO")
MuCO <- mxAlgebra(expression= expMeanCO+(betaCO*sexCO), name="MuCO")
meanGCO2 <- mxAlgebra(cbind(MuCO, MuCO), name="expMeanCO2")
 
 
covMZMN <- mxAlgebra( expression= rbind( cbind(AMN+CMN+EMN,       AMN+CMN),
                                         cbind(   AMN+CMN,    AMN+CMN+EMN)), name="expCovMZMN" )
covDZMN <- mxAlgebra( expression= rbind( cbind(AMN+CMN+EMN, 0.5%x%AMN+CMN),
                                         cbind(0.5%x%AMN+CMN, AMN+CMN+EMN)), name="expCovDZMN" )
 
covMZCO <- mxAlgebra( expression= rbind( cbind(ACO+CCO+ECO,       ACO+CCO),
                                         cbind(   ACO+CCO,    ACO+CCO+ECO)), name="expCovMZCO" )
covDZCO <- mxAlgebra( expression= rbind( cbind(ACO+CCO+ECO, 0.5%x%ACO+CCO),
                                         cbind(0.5%x%ACO+CCO, ACO+CCO+ECO)), name="expCovDZCO" )
 
 
### Combine Groups
objMZMN <- mxExpectationNormal( covariance="expCovMZMN", means="expMeanMN2", dimnames=selVarsMN )
objDZMN <- mxExpectationNormal( covariance="expCovDZMN", means="expMeanMN2", dimnames=selVarsMN )
 
objMZCO <- mxExpectationNormal( covariance="expCovMZCO", means="expMeanCO2", dimnames=selVarsCO )
objDZCO <- mxExpectationNormal( covariance="expCovDZCO", means="expMeanCO2", dimnames=selVarsCO )
 
parsMN <- list( XMN, YMN, ZMN, AlMN, ClMN, ElMN, pathAsMN, pathCsMN, pathEsMN, pathflMN, covAMN, covCMN, 
                covEMN, meanGMN, meanGMN2, regCoefMN, sexCovMN, MuMN, pathFlGrowthMN, pathFlBDMN, covPMN, matIMN, invSDMN)
 
parsCO <- list( XCO, YCO, ZCO, AlCO, ClCO, ElCO, pathAsCO, pathCsCO, pathEsCO, pathflCO, covACO, covCCO,
                covECO, meanGCO, meanGCO2, regCoefCO, sexCovCO, MuCO, pathFlGrowthCO, pathFlBDCO, covPCO, matICO, invSDCO)
 
funML <- mxFitFunctionML()
 
mtfsModelMZ <- mxModel( covMZMN, mtfsDataMZ, objMZMN, funML, name="MZMN", parsMN)
mtfsModelDZ <- mxModel( covDZMN, mtfsDataDZ, objDZMN, funML, name="DZMN", parsMN)
 
caddModelMZ <- mxModel( covMZCO, caddDataMZ, objMZCO, funML, name="MZCO", parsCO)
caddModelDZ <- mxModel( covDZCO, caddDataDZ, objDZCO, funML, name="DZCO", parsCO)
 
fitMLMN <- mxFitFunctionMultigroup(c("MZMN.fitfunction","DZMN.fitfunction"))
LGMN <- mxModel( "LinearGrowthACE", mtfsModelMZ, mtfsModelDZ, fitMLMN)
mxOption(NULL,"Default optimizer","SLSQP")
LGFitMN <- mxRun(LGMN) ### Run model
RefModLGMN <- mxRefModels(LGFitMN, run=TRUE)
 
summary(LGFitMN, refModels=RefModLGMN)
LGFitMN$output$status$code
round(LGFitMN$output$algebras$MZMN.AlMN, 5)
round(LGFitMN$output$algebras$MZMN.ClMN, 5)
round(LGFitMN$output$algebras$MZMN.ElMN, 5)
 
fitMLCO <- mxFitFunctionMultigroup(c("MZCO.fitfunction","DZCO.fitfunction"))
LGCO <- mxModel("LinearGrowthACE", caddModelMZ, caddModelDZ, fitMLCO)
mxOption(NULL,"Default optimizer","SLSQP")
LGFitCO <- mxRun(LGCO) ### Run model
RefModLGCO <- mxRefModels(LGFitCO, run=TRUE)
 
summary(LGFitCO, refModels=RefModLGCO)
LGFitCO$output$status$code
round(LGFitCO$output$algebras$MZCO.AlCO, 5)
round(LGFitCO$output$algebras$MZCO.ClCO, 5)
round(LGFitCO$output$algebras$MZCO.ElCO, 5)
 
round(LGFitMN$output$algebras$MZMN.iSDMN %*% LGFitMN$output$matrices$MZMN.flBDMN[,1], 2)
round(LGFitCO$output$algebras$MZCO.iSDCO %*% LGFitCO$output$matrices$MZCO.flBDCO[,1], 2)
AdminRobK's picture
Offline
Joined: 01/24/2014 - 12:15
Explain script?

Could you explain what sort of model the script is meant to fit? I'm actually pretty confused. The script in your post doesn't even define certain variables, like nvMN or selVarsCO.

Also, could you provide the MxAlgebra you used previously to standardize "flBDMN" and "flBDCO"?

szellers's picture
Offline
Joined: 04/04/2018 - 14:52
Whoops!

Ah, I had left out the data prep code for this script. I am using two longitudinal twin samples and fitting a common path model with three indicators and a linear growth model with three indicators in CO sample and five indicators in MN sample. I am also allowing the common path latent factor to correlate with the slope and intercept of the linear growth model.

Previously to standardize flBDMN and flBDCO I had specified, inside the model, the matrices (separately in each sample, below is CO example)

matICO <- mxMatrix( type="Iden", nrow=3, ncol=3, name="ICO")
covPCO <- mxAlgebra( expression=ACO[1:3, 1:3]+CCO[1:3, 1:3]+ECO[1:3, 1:3], name="VCO" )
invSDCO <- mxAlgebra( expression=solve(sqrt(ICO*VCO)), name="iSDCO")

Then after running the model, I used the output to run the algebra. I took the unstandardized loadings and the iSDCO matrix to calculate standardized values using this command

LGFitCO$output$algebras$MZCO.iSDCO %*% LGFitCO$output$matrices$MZCO.flBDCO
 

Full script below

########################
### Biometric Models ###
########################
###
### Select observed variables
VarsCO <- c("CDres", "ADres", "ADHDres", "FREQ1", "FREQ2", "FREQ3")
nvCO <- length(VarsCO)
ntvCO <- nvCO*2
selVarsCO <- paste(VarsCO, c(rep(0,nvCO), rep(1,nvCO)), sep="")
 
VarsMN <- c("CDres", "ADres", "ADHDres", "mj_12m_freq_14", "mj_12m_freq_17", "mj_12m_freq_20", "mj_12m_freq_24", "mj_12m_freq_29")
nvMN <- length(VarsMN)
ntvMN <- nvMN*2
selVarsMN <- paste(VarsMN, c(rep(0,nvMN), rep(1,nvMN)), sep="")
 
############
#definition variables
DefVarsCO <- c("age1twinpairC", "age2twinpairC", "age3twinpairC", "intercept", "sex", "zero")
DefVarsMN <- c("age14twinpair", "age17twinpair", "age20twinpair", "age24twinpair", "age29twinpair", "intercept", "sex", "zero")
 
 
 
### Select Data for Analysis
caddMzData    <- subset(cadd.w, zygosity==1, select=c(selVarsCO, DefVarsCO))
caddDzData    <- subset(cadd.w, zygosity==2, select=c(selVarsCO, DefVarsCO))
 
caddDataMZ   <- mxData( observed=caddMzData, type="raw" )
caddDataDZ   <- mxData( observed=caddDzData, type="raw" )
 
mtfsMzData    <- subset(mtfs.w, zygosity==1, select=c(selVarsMN, DefVarsMN))
mtfsDzData    <- subset(mtfs.w, zygosity==2, select=c(selVarsMN, DefVarsMN))
 
mtfsDataMZ   <- mxData( observed=mtfsMzData, type="raw" )
mtfsDataDZ   <- mxData( observed=mtfsDzData, type="raw" )
 
 
 
####################################################################
nl <- 3 # number of latent factors
# Matrices ac, cc, and ec to store a, c, and e path coefficients for latent phenotype(s)
XMN <- mxMatrix(type="Lower", nrow=nl, ncol=nl, free=TRUE,
                values=c(.6, 
                         .1, .1, 
                         .1, .1, 0),
                labels=c("x11MN",
                         "x21MN", "x22MN", 
                         "x31MN", "x32MN", "x33MN"),
                lbound=.00001, name="XMN" )
 
YMN <- mxMatrix(type="Lower", nrow=nl, ncol=nl, free=TRUE,
                values=c(.6,
                         .1, .1, 
                         .1, .1, .1),
                labels=c("y11MN",
                         "y21MN", "y22MN",
                         "y31MN", "y32MN", "y33MN"),
                lbound=.00001, name="YMN" )
 
ZMN <- mxMatrix(type="Lower", nrow=nl, ncol=nl, free=TRUE,
                values=c(.6,
                         .1, .1,
                         .1, .1, .1),
                labels=c("z11MN",
                         "z21MN", "z22MN",
                         "z31MN", "z32MN", "z33MN"),
                lbound=.00001, name="ZMN" )
 
XCO <- mxMatrix(type="Lower", nrow=nl, ncol=nl, free=TRUE,
                values=c(.6, 
                         .1, .1, 
                         .1, .1, 0),
                labels=c("x11CO",
                         "x21CO", "x22CO", 
                         "x31CO", "x32CO", "x33CO"),
                lbound=.00001, name="XCO" )
 
YCO <- mxMatrix(type="Lower", nrow=nl, ncol=nl, free=TRUE,
                values=c(.6,
                         .1, .1, 
                         .1, .1, .1),
                labels=c("y11CO",
                         "y21CO", "y22CO",
                         "y31CO", "y32CO", "y33CO"),
                lbound=.00001, name="YCO" )
 
ZCO <- mxMatrix(type="Lower", nrow=nl, ncol=nl, free=TRUE,
                values=c(.6,
                         .1, .1,
                         .1, .1, .1),
                labels=c("z11CO",
                         "z21CO", "z22CO",
                         "z31CO", "z32CO", "z33CO"),
                lbound=.00001, name="ZCO" )
 
 
AlMN <- mxAlgebra(XMN %*% t(XMN), name="AlMN")
ClMN <- mxAlgebra(YMN %*% t(YMN), name="ClMN")
ElMN <- mxAlgebra(ZMN %*% t(ZMN), name="ElMN")
 
AlCO <- mxAlgebra(XCO %*% t(XCO), name="AlCO")
ClCO <- mxAlgebra(YCO %*% t(YCO), name="ClCO")
ElCO <- mxAlgebra(ZCO %*% t(ZCO), name="ElCO")
 
 
# Matrix and Algebra for constraint on variance of latent phenotype
covarLPMN <- mxAlgebra( expression= AlMN[1,1] + ClMN[1,1] + ElMN[1,1], name="CovarLPMN" )
varLPMN <- mxAlgebra( expression= diag2vec(CovarLPMN), name="VarLPMN" )
unitMN <- mxMatrix( type="Unit", nrow=1, ncol=1, name="UnitMN")
varLP1MN <- mxConstraint( expression=VarLPMN == UnitMN, name="varLP1MN")
 
covarLPCO <- mxAlgebra( expression= AlCO[1,1] + ClCO[1,1] + ElCO[1,1], name="CovarLPCO" )
varLPCO <- mxAlgebra( expression= diag2vec(CovarLPCO), name="VarLPCO" )
unitCO <- mxMatrix( type="Unit", nrow=1, ncol=1, name="UnitCO")
varLP1CO <- mxConstraint( expression=VarLPCO == UnitCO, name="varLP1CO")
 
 
 
 
# Matrices as, cs, and es to store a, c, and e path coefficients for specific factors
pathAsMN <- mxMatrix(type="Diag", nrow=nvMN, ncol=nvMN, free=TRUE, values=.3,
                     labels=c("a11MN", "a22MN", "a33MN", "a44MN", "a55MN", "a66MN", "a77MN", "a88MN"),
                     lbound=.00001, name="asMN" )
pathCsMN <- mxMatrix(type="Diag", nrow=nvMN, ncol=nvMN, free=TRUE, values=.3,
                     labels=c("c11MN","c22MN", "c33MN", "c44MN", "c55MN", "c66MN", "c77MN", "c88MN"),
                     lbound=.00001, name="csMN" )
pathEsMN <- mxMatrix(type="Diag", nrow=nvMN, ncol=nvMN, free=TRUE, values=.3,
                     labels=c("e11MN", "e22MN", "e33MN", "e44MN", "e55MN", "e66MN", "e77MN", "e88MN"),
                     lbound=.00001, name="esMN" )
 
pathAsCO <- mxMatrix(type="Diag", nrow=nvCO, ncol=nvCO, free=TRUE, values=.3,
                     labels=c("a11CO", "a22CO", "a33CO", "a44CO", "a55CO", "a66CO"),
                     lbound=.00001, name="asCO" )
pathCsCO <- mxMatrix(type="Diag", nrow=nvCO, ncol=nvCO, free=TRUE, values=.3,
                     labels=c("c11CO","c22CO", "c33CO", "c44CO", "c55CO", "c66CO"),
                     lbound=.00001, name="csCO" )
pathEsCO <- mxMatrix(type="Diag", nrow=nvCO, ncol=nvCO, free=TRUE, values=.3,
                     labels=c("e11CO", "e22CO", "e33CO", "e44CO", "e55CO", "e66CO"),
                     lbound=.00001, name="esCO" )
 
 
 
 
pathFlGrowthMN <- mxMatrix(type="Full", nrow=5, ncol=nl, free=FALSE, 
                           labels=c("data.zero", "data.intercept", "data.age14twinpair",
                                    "data.zero", "data.intercept", "data.age17twinpair",
                                    "data.zero", "data.intercept", "data.age20twinpair",
                                    "data.zero", "data.intercept", "data.age24twinpair",
                                    "data.zero", "data.intercept", "data.age29twinpair"), name="flGrowthMN", byrow=T)
 
pathFlGrowthCO <- mxMatrix(type="Full", nrow=3, ncol=nl, free=FALSE, 
                           labels=c("data.zero", "data.intercept", "data.age1twinpairC",
                                    "data.zero", "data.intercept", "data.age2twinpairC",
                                    "data.zero", "data.intercept", "data.age3twinpairC"), name="flGrowthCO", byrow=T)
 
 
flFree   <- c(FALSE, FALSE, FALSE,
              TRUE, FALSE, FALSE,
              TRUE, FALSE, FALSE)
 
flValues <- c(1, 0, 0,
              .6, 0, 0,
              .6, 0, 0)
 
#set first loading to 1 and make it fixed for BD to manifest variables
pathFlBDMN <- mxMatrix(type="Full", nrow=3, ncol=nl, free=flFree, 
                       values=flValues, lbound=.000001, labels=c("fBD1MN", NA, NA,
                                                                 "fBD2MN", NA, NA,
                                                                 "fBD3MN", NA, NA), name="flBDMN", byrow=T)
 
pathFlBDCO <- mxMatrix(type="Full", nrow=3, ncol=nl, free=flFree, 
                       values=flValues, lbound=.000001, labels=c("fBD1CO", NA, NA,
                                                                 "fBD2CO", NA, NA,
                                                                 "fBD3CO", NA, NA), name="flBDCO", byrow=T)
 
#combine factor loading matrices
pathflMN <- mxAlgebra(rbind(flBDMN, flGrowthMN), name="flMN")
 
pathflCO <- mxAlgebra(rbind(flBDCO, flGrowthCO), name="flCO")
 
 
 
 
covAMN <- mxAlgebra( expression=flMN %*% AlMN %*% t(flMN) + asMN %*% t(asMN), name="AMN" )
covCMN <- mxAlgebra( expression=flMN %*% ClMN %*% t(flMN) + csMN %*% t(csMN), name="CMN" )
covEMN <- mxAlgebra( expression=flMN %*% ElMN %*% t(flMN) + esMN %*% t(esMN), name="EMN" )
 
covACO <- mxAlgebra( expression=flCO %*% AlCO %*% t(flCO) + asCO %*% t(asCO), name="ACO" )
covCCO <- mxAlgebra( expression=flCO %*% ClCO %*% t(flCO) + csCO %*% t(csCO), name="CCO" )
covECO <- mxAlgebra( expression=flCO %*% ElCO %*% t(flCO) + esCO %*% t(esCO), name="ECO" )
 
matIMN <- mxMatrix( type="Iden", nrow=3, ncol=3, name="IMN")
covPMN <- mxAlgebra( expression=AMN[1:3, 1:3]+CMN[1:3, 1:3]+EMN[1:3, 1:3], name="VMN" )
invSDMN <- mxAlgebra( expression=solve(sqrt(IMN*VMN)), name="iSDMN")
 
matICO <- mxMatrix( type="Iden", nrow=3, ncol=3, name="ICO")
covPCO <- mxAlgebra( expression=ACO[1:3, 1:3]+CCO[1:3, 1:3]+ECO[1:3, 1:3], name="VCO" )
invSDCO <- mxAlgebra( expression=solve(sqrt(ICO*VCO)), name="iSDCO")
 
 
 
meanGMN <- mxMatrix( type="Full", nrow=1, ncol=nvMN, free=TRUE, values=.3, 
                     labels=c('m1MN', 'm2MN', 'm3MN', 'm4MN', 'm5MN', 'm6MN', 'm7MN', 'm8MN'), name="expMeanMN" )
regCoefMN <- mxMatrix(type="Full", nrow=1, ncol=nvMN, free=TRUE, values=0, 
                      labels=c("beta1MN", "beta2MN", "beta3MN", "beta4MN", "beta5MN", "beta6MN", "beta7MN", "beta8MN"), name="betaMN")
sexCovMN <- mxMatrix(type="Full", nrow=1, ncol=nvMN, free=FALSE, labels="data.sex", name="sexMN")
MuMN <- mxAlgebra(expression= expMeanMN+(betaMN*sexMN), name="MuMN")
meanGMN2 <- mxAlgebra(cbind(MuMN, MuMN), name="expMeanMN2")
 
meanGCO <- mxMatrix( type="Full", nrow=1, ncol=nvCO, free=TRUE, values=.3, 
                     labels=c('m1CO', 'm2CO', 'm3CO', 'm4CO', 'm5CO', 'm6CO'), name="expMeanCO" )
regCoefCO <- mxMatrix(type="Full", nrow=1, ncol=nvCO, free=TRUE, values=0, 
                      labels=c("beta1CO", "beta2CO", "beta3CO", "beta4CO", "beta5CO", "beta6CO"), name="betaCO")
sexCovCO <- mxMatrix(type="Full", nrow=1, ncol=nvCO, free=FALSE, labels="data.sex", name="sexCO")
MuCO <- mxAlgebra(expression= expMeanCO+(betaCO*sexCO), name="MuCO")
meanGCO2 <- mxAlgebra(cbind(MuCO, MuCO), name="expMeanCO2")
 
 
 
 
covMZMN <- mxAlgebra( expression= rbind( cbind(AMN+CMN+EMN,       AMN+CMN),
                                         cbind(   AMN+CMN,    AMN+CMN+EMN)), name="expCovMZMN" )
covDZMN <- mxAlgebra( expression= rbind( cbind(AMN+CMN+EMN, 0.5%x%AMN+CMN),
                                         cbind(0.5%x%AMN+CMN, AMN+CMN+EMN)), name="expCovDZMN" )
 
covMZCO <- mxAlgebra( expression= rbind( cbind(ACO+CCO+ECO,       ACO+CCO),
                                         cbind(   ACO+CCO,    ACO+CCO+ECO)), name="expCovMZCO" )
covDZCO <- mxAlgebra( expression= rbind( cbind(ACO+CCO+ECO, 0.5%x%ACO+CCO),
                                         cbind(0.5%x%ACO+CCO, ACO+CCO+ECO)), name="expCovDZCO" )
 
 
### Combine Groups
objMZMN <- mxExpectationNormal( covariance="expCovMZMN", means="expMeanMN2", dimnames=selVarsMN )
objDZMN <- mxExpectationNormal( covariance="expCovDZMN", means="expMeanMN2", dimnames=selVarsMN )
 
objMZCO <- mxExpectationNormal( covariance="expCovMZCO", means="expMeanCO2", dimnames=selVarsCO )
objDZCO <- mxExpectationNormal( covariance="expCovDZCO", means="expMeanCO2", dimnames=selVarsCO )
 
parsMN <- list( XMN, YMN, ZMN, AlMN, ClMN, ElMN, pathAsMN, pathCsMN, pathEsMN, pathflMN, covAMN, covCMN, 
                covEMN, meanGMN, meanGMN2, regCoefMN, sexCovMN, MuMN, pathFlGrowthMN, pathFlBDMN, covPMN, matIMN, invSDMN)
 
parsCO <- list( XCO, YCO, ZCO, AlCO, ClCO, ElCO, pathAsCO, pathCsCO, pathEsCO, pathflCO, covACO, covCCO,
                covECO, meanGCO, meanGCO2, regCoefCO, sexCovCO, MuCO, pathFlGrowthCO, pathFlBDCO, covPCO, matICO, invSDCO)
 
funML <- mxFitFunctionML()
 
mtfsModelMZ <- mxModel( covMZMN, mtfsDataMZ, objMZMN, funML, name="MZMN", parsMN)
mtfsModelDZ <- mxModel( covDZMN, mtfsDataDZ, objDZMN, funML, name="DZMN", parsMN)
 
caddModelMZ <- mxModel( covMZCO, caddDataMZ, objMZCO, funML, name="MZCO", parsCO)
caddModelDZ <- mxModel( covDZCO, caddDataDZ, objDZCO, funML, name="DZCO", parsCO)
 
fitMLMN <- mxFitFunctionMultigroup(c("MZMN.fitfunction","DZMN.fitfunction"))
LGMN <- mxModel( "LinearGrowthACE", mtfsModelMZ, mtfsModelDZ, fitMLMN)
mxOption(NULL,"Default optimizer","SLSQP")
LGFitMN <- mxRun(LGMN) ### Run model
RefModLGMN <- mxRefModels(LGFitMN, run=TRUE)
 
fitMLCO <- mxFitFunctionMultigroup(c("MZCO.fitfunction","DZCO.fitfunction"))
LGCO <- mxModel("LinearGrowthACE", caddModelMZ, caddModelDZ, fitMLCO)
mxOption(NULL,"Default optimizer","SLSQP")
LGFitCO <- mxRun(LGCO) ### Run model
RefModLGCO <- mxRefModels(LGFitCO, run=TRUE)
 
 
summary(LGFitMN, refModels=RefModLGMN)
summary(LGFitCO, refModels=RefModLGCO)
 
 
### Get standardized loadings
LGFitMN$output$algebras$MZMN.iSDMN %*% LGFitMN$output$matrices$MZMN.flBDMN
LGFitCO$output$algebras$MZCO.iSDCO %*% LGFitCO$output$matrices$MZCO.flBDCO
AdminRobK's picture
Offline
Joined: 01/24/2014 - 12:15
I think I get it

So, for (say) the Minnesota cohort, "A1MN" is the name of the additive-genetic covariance matrix of the 3 common factors, and "asMN" is the name of the unique additive-genetic covariance matrix of the observable phenotypes, right? If so, then I guess

AlMN[1,1] + ClMN[1,1] + ElMN[1,1]

would be the variance of the first common factor. So, you could create additional algebras,

sd_L1 <- mxAlgebra( sqrt(AlMN[1,1] + ClMN[1,1] + ElMN[1,1]), name="sdL1")
flBDMN_std <- mxAlgebra( iSDMN %*% (sdL1 %x% flBDMN), name="flBDMNstd")

and the algebra named "flBDMNstd" would be the standardized loadings you want for the Minnesota cohort. Does that sound right?

szellers's picture
Offline
Joined: 04/04/2018 - 14:52
I think I got it

Hello,

I tried to go through the steps you'd originally suggested above, working with the output from my model, I ran:

### Get standardized loadings
LGFitMN$output$algebras$MZMN.VMN
VarMN <- LGFitMN$output$algebras$MZMN.AlMN[1,1] + LGFitMN$output$algebras$MZMN.ClMN[1,1] + LGFitMN$output$algebras$MZMN.ElMN[1,1]
SDMN <- sqrt(VarMN)
MNStep1 <- SDMN %x% LGFitMN$output$matrices$MZMN.flBDMN 
MNStep2 <- LGFitMN$output$algebras$MZMN.iSDMN %*% MNStep1
MNStep2
 
LGFitCO$output$algebras$MZCO.VCO
VarCO <- LGFitCO$output$algebras$MZCO.AlCO[1,1] + LGFitCO$output$algebras$MZCO.ClCO[1,1] + LGFitCO$output$algebras$MZCO.ElCO[1,1]
SDCO <- sqrt(VarCO)
COStep1 <- SDCO %x% LGFitCO$output$matrices$MZCO.flBDCO 
COStep2 <- LGFitCO$output$algebras$MZCO.iSDCO %*% COStep1
COStep2

These loadings seem to be in agreement with what I'd expect them to be, given previous results.
The full script specifying all the matrices is the same that I'd posted in my above reply.

AdminRobK's picture
Offline
Joined: 01/24/2014 - 12:15
OK, looks good. We agree

OK, looks good. We agree about how to go about it.

Log in or register to post comments