Login

News

CI estimation of standardized path coefficients in common path way model

Posted on
No user picture. Liz Joined: 04/20/2016
Forums
Behavioral Genetics Models
Hi there,

Now I have a problem about CI estimation in a common pathway model. I can not estimate the CI of the path coeffcients of latent factors. I added following codes to the script:

matIl <- mxMatrix( type="Iden", nrow=nf, ncol=nf, name="Il")
invSDl <- mxAlgebra( expression=solve(sqrt(Il*V)), name="iSDl")
sta <- mxAlgebra(iSDl %*% al, name="sta") # standarized path coefficients of A
stc <- mxAlgebra(iSDl %*% cl, name="stc") # standarized path coefficients of C
ste <- mxAlgebra(iSDl %*% el, name="ste") # standarized path coefficients of E

The nf denotes the number of latant factors. If I use the nv (the number of all components in the model), the results will only show the parameters of components, not the latant factor which I'm interesed with. But now, it reply with "Error in runHelper(model, frontendStart, intervals, silent, suppressWarnings, : Setting improper index (1) from vector of Rf_length 0". I don't know how can I settle it (If I remove these codes, the script could run smoothly without estimating the standardized path coefficients). How can I reivise the codes? Many thanks.

Best,

Lizhi

Replied on Mon, 11/07/2016 - 05:01
No user picture. Liz Joined: 04/20/2016

No. How to constrained the latent variance to 1? Not the following codes?

matIl <- mxMatrix( type="Iden", nrow=nf, ncol=nf, name="Il")
invSDl <- mxAlgebra( expression=solve(sqrt(Il*V)), name="iSDl")

And why should the latents variance be fixed on 1. Actually, if I only estimatelatent path parameters, not their CI. The script could run smoothly. Many thanks.

Log in or register to post comments
Replied on Mon, 11/07/2016 - 10:27
Picture of user. AdminRobK Joined: 01/24/2014

It's pretty much impossible to troubleshoot without seeing the full script you're using.

Actually, if I only estimatelatent path parameters, not their CI. The script could run smoothly.

It's likely the problem is connected to your mxCI() statement. Seeing that in the context of the full script would help with troubleshooting.

And why should the latents variance be fixed on 1.

To identify the scale of the latent variables, there needs to be a constraint on their variances (or alternately, on the paths running from latent to manifest variables). However, if you're adapting an existing script, the constraint on the latent variances is probably implicitly built-in.

Log in or register to post comments
Replied on Mon, 11/07/2016 - 11:33
No user picture. Liz Joined: 04/20/2016

Many thanks, Rob. Attached is my full script. It is a two factor common path way model. The problem is the part of CI setup, if I deleted this part, the script could run smoothly:

#Standeraized parameters and CI
sta <- mxAlgebra(iSD %*% al, name="sta")
stc <- mxAlgebra(iSD %*% cl, name="stc")
ste <- mxAlgebra(iSD %*% el, name="ste")
StandardizedA <- mxAlgebra(A/V, name="StandardizedA")
StandardizedC <- mxAlgebra(C/V, name="StandardizedC")
StandardizedE <- mxAlgebra(E/V, name="StandardizedE")
CI <- mxCI(c('sta', 'stc', 'ste','StandardizedA', 'StandardizedC', 'StandardizedE'))

The matrix dimension of iSD (3x3x3) is different from that of "al“ (2x2x2), ”cl” or “el”. I'm still not clear why the multiply of "iSD" by "al" or "a" in the one trait ACE model indicate the standardization.

Full script as below:

nv <- 3 # number of variables
ntv <- nv*2 # number of total variables
selVars <- paste(Vars,c(rep(1,nv),rep(2,nv)),sep="_")

nf <- 2 # number of factors
nsub1 <- 1 # number of subphenotypes in factor 1
nsub2 <- 2 # number of subphenotypes in factor 2

# Select Data for Analysis
mzData <- subset(nl, zyg==1, selVars)
dzData <- subset(nl, zyg==2, selVars)
describe(mzData, skew=F)
describe(dzData, skew=F)
dim(mzData)
dim(dzData)

# Generate Descriptive Statistics
colMeans(mzData,na.rm=TRUE)
colMeans(dzData,na.rm=TRUE)
cov(mzData,use="complete")
cov(dzData,use="complete")
cor(mzData,use="complete")
cor(dzData,use="complete")

# Create Labels for Lower Triangular Matrices (fancy shorthand)
# Labels for a, c & e pathways from A, C & E latent factors to latent INTECEPT & SLOPE
AlLabs <- paste("al", do.call(c, sapply(seq(1, nf), function(x){ paste(x:nf, x,sep="") })), sep="")
ClLabs <- paste("cl", do.call(c, sapply(seq(1, nf), function(x){ paste(x:nf, x,sep="") })), sep="")
ElLabs <- paste("el", do.call(c, sapply(seq(1, nf), function(x){ paste(x:nf, x,sep="") })), sep="")

# Labels for a, c & e pathways from A, C & E residual to observed variables
AsLabs <- paste("as",1:nv,1:nv,sep="_")
CsLabs <- paste("cs",1:nv,1:nv,sep="_")
EsLabs <- paste("es",1:nv,1:nv,sep="_")

meMZVals <- rep(0,ntv) # start value for MZ means
meDZVals <- rep(0,ntv) # start value for DZ means
meanMZLabs <- paste("mMZ",1:ntv,sep="_")
meanDZLabs <- paste("mDZ",1:ntv,sep="_")

comFLabs1 <- paste("f",1,1:nsub1,sep="_")
comFLabs2 <- paste("f",2,1:nsub2,sep="_")

pathAl <- mxMatrix(name = "al", type = "Lower", nrow = 2, ncol = 2, labels = AlLabs, free=T, value=0.5)
pathCl <- mxMatrix(name = "cl", type = "Lower", nrow = 2, ncol = 2, labels = ClLabs, free=T, value=0.5)
pathEl <- mxMatrix(name = "el", type = "Lower", nrow = 2, ncol = 2, labels = ElLabs, free=T, value=0.5)

# Matrix f for factor loadings on latent phenotype
pathF1 <- mxMatrix( type="Full", nrow=nsub1, ncol=1, free=F, values=1, #lbound=-1, ubound=1,
labels=comFLabs1, name="fl1" )
pathF2 <- mxMatrix( type="Full", nrow=1, ncol=1, free=T, values=0.5, #lbound=-1, ubound=1,
labels=comFLabs2, name="fl2" ) #nrow= nsub2-1
pathTmp1 <- mxMatrix( type="Full", nrow=nsub2, ncol=1, free=F, values=0, name="patht1")
pathTmp2 <- mxMatrix( type="Full", nrow=nsub1, ncol=1, free=F, values=0, name="patht2")

pathF <- mxAlgebra( expression=cbind(rbind(fl1,patht1),rbind(patht2,fl1,fl2)), name="fl" )

# Matrices as, cs, and es to store a, c, and e path coefficients for specific factors
pathAs <- mxMatrix( type="Diag", nrow=nv, ncol=nv, free=TRUE, values=.5, labels=AsLabs,name="as")
pathCs <- mxMatrix( type="Diag", nrow=nv, ncol=nv, free=TRUE, values=.5, labels=CsLabs,name="cs")
pathEs <- mxMatrix( type="Diag", nrow=nv, ncol=nv, free=TRUE, values=.5, labels=EsLabs,name="es" )

# Matrices A, C, and E compute variance components
covA <- mxAlgebra( expression=fl %&% (al %*% t(al)) + as %*% t(as), name="A" )
covC <- mxAlgebra( expression=fl %&% (cl %*% t(cl)) + cs %*% t(cs), name="C" )
covE <- mxAlgebra( expression=fl %&% (el %*% t(el)) + es %*% t(es), name="E" )

# Algebra to compute total variances and standard deviations (diagonal only)
covP <- mxAlgebra( expression=A+C+E, name="V" )
matI <- mxMatrix( type="Iden", nrow=nv, ncol=nv, name="I")
invSD <- mxAlgebra( expression=solve(sqrt(I*V)), name="iSD")
corPh <- mxAlgebra(name ="rPh", expression = iSD%*%V%*%iSD)

# Matrix and Algebra for constraint on variance of latent phenotype
covarLP <- mxAlgebra( expression= al %*% t(al) + cl %*% t(cl) + el %*% t(el), name="CovarLP" )
varLP <- mxAlgebra( expression= diag2vec(CovarLP), name="VarLP" )
unit <- mxMatrix( type="Unit", nrow=2, ncol=1, name="Unit")

#standardise the f parameters. Note the above iSD is just the invers of D_y
matI3 <- mxMatrix( type="Iden", nrow=nf, ncol=nf, name="I3")
standF <- mxAlgebra(name ="Sfl", expression = iSD %*% fl %*% sqrt(I3*CovarLP) )

#Standeraized parameters and CI
sta <- mxAlgebra(iSD %*% al, name="sta")
stc <- mxAlgebra(iSD %*% cl, name="stc")
ste <- mxAlgebra(iSD %*% el, name="ste")
StandardizedA <- mxAlgebra(A/V, name="StandardizedA")
StandardizedC <- mxAlgebra(C/V, name="StandardizedC")
StandardizedE <- mxAlgebra(E/V, name="StandardizedE")

CI <- mxCI(c('sta', 'stc', 'ste','StandardizedA', 'StandardizedC', 'StandardizedE'))

# Algebra for expected Mean and Variance/Covariance Matrices in MZ & DZ twins
meanMZ <- mxMatrix( type="Full", nrow=1, ncol=ntv, free=T, values=0, labels=meanMZLabs, name="expMeanMZ" )
meanDZ <- mxMatrix( type="Full", nrow=1, ncol=ntv, free=T, values=0, labels=meanDZLabs, name="expMeanDZ" )
covMZ <- mxAlgebra( expression=
rbind( cbind(A+C+E , A+C),
cbind(A+C , A+C+E)), name="expCovMZ" )
covDZ <- mxAlgebra( expression=
rbind( cbind(A+C+E , 0.5%x%A+C),
cbind(0.5%x%A+C , A+C+E)), name="expCovDZ" )

# Algebra for expected variance/covariance matrix in MZ
expCovMZ <- mxAlgebra(name = "expCovMZ",
expression = rbind (cbind(A+C+E, A+C),
cbind(A+C, A+C+E)))

# Algebra for expected variance/covariance matrix in DZ
expCovDZ <- mxAlgebra(name = "expCovDZ",
expression = rbind (cbind(A+C+E, 0.5%x%A+C),
cbind(0.5%x%A+C, A+C+E)))

# Data objects for Multiple Groups
dataMZ <- mxData( observed=mzData, type="raw" )
dataDZ <- mxData( observed=dzData, type="raw" )

# Objectives for MZ and DZ groups
objMZ <- mxExpectationNormal( covariance="expCovMZ", means="expMeanMZ", dimnames=selVars )
objDZ <- mxExpectationNormal( covariance="expCovDZ", means="expMeanDZ", dimnames=selVars )
funML <- mxFitFunctionML()

# Combine Groups
pars <- list(pathAl, pathCl, pathEl, pathF1, pathF2, pathTmp1, pathTmp2, pathF, pathAs, pathCs, pathEs, matI3, standF, covA, covC, covE, covP, matI, invSD, corPh, covarLP, varLP, unit, meanMZ, meanDZ,sta, stc, ste,StandardizedA, StandardizedC, StandardizedE)
modelMZ <- mxModel( pars, covMZ, dataMZ, objMZ, funML, name="MZ" )
modelDZ <- mxModel( pars, covDZ, dataDZ, objDZ, funML,name="DZ" )
fitML <- mxFitFunctionMultigroup(c("MZ.fitfunction","DZ.fitfunction") )
ComAceModel <- mxModel( "ComACE", pars, modelMZ, modelDZ, fitML, CI)

# Run CholACE model
ComAceFit <- mxRun(ComAceModel, intervals=TRUE)

Log in or register to post comments
Replied on Mon, 11/14/2016 - 10:50
Picture of user. AdminRobK Joined: 01/24/2014

We've had a tough time troubleshooting the problem you're reporting, since we don't have the data to run your script. Can you share your data, or maybe generate data that's similar to it?
Log in or register to post comments
Replied on Tue, 11/15/2016 - 06:42
No user picture. Liz Joined: 04/20/2016

Many thanks for your kind help, Rob! I have figured it out yestarday. It's about the matrix setting. For the "V" contains the varience of all the phenotype and latent factor, hence I make a new matrix which only contains the variance of latent factors, and then multiply this matrix with I (2x2) and al. The script runs successfully and the results seems correct.

Many thanks and best wishes,

Lizhi

Log in or register to post comments
Replied on Tue, 11/15/2016 - 11:18
Picture of user. AdminRobK Joined: 01/24/2014

There's something I forgot to mention. I think you might still be missing something crucial to model identification here. You don't actually have the MxConstraint needed for model identification. I guess it would be something like

identifyingConstraint <- mxConstraint(VarLP == Unit, name="ic")

which you would need to then put into your MxModel.
Log in or register to post comments
Replied on Tue, 11/15/2016 - 15:44
No user picture. Liz Joined: 04/20/2016

You mean the constraint of latent valence? Are the following codes enough, or I should add the code you suggested to it?

# Matrix and Algebra for constraint on variance of latent phenotype
covarLP <- mxAlgebra( expression= al %*% t(al) + cl %*% t(cl) + el %*% t(el), name="CovarLP" )
varLP <- mxAlgebra( expression= diag2vec(CovarLP), name="VarLP" )
unit <- mxMatrix( type="Unit", nrow=2, ncol=1, name="Unit")

Log in or register to post comments
Replied on Tue, 11/15/2016 - 16:12
Picture of user. AdminRobK Joined: 01/24/2014

Something like this:

# Matrix and Algebra for constraint on variance of latent phenotype
covarLP <- mxAlgebra( expression= al %*% t(al) + cl %*% t(cl) + el %*% t(el), name="CovarLP" )
varLP <- mxAlgebra( expression= diag2vec(CovarLP), name="VarLP" )
unit <- mxMatrix( type="Unit", nrow=2, ncol=1, name="Unit")
identifyingConstraint <- mxConstraint(VarLP == Unit, name="ic")

And you would need to put indentifyingConstraint somewhere in your MxModel.
Log in or register to post comments
Replied on Tue, 11/08/2016 - 17:56
Picture of user. neale Joined: 07/31/2009

I believe that:

pathAl <- mxMatrix(name = "al", type = "Lower", nrow = 2, ncol = 2, labels = AlLabs, free=T, value=0.5)
pathCl <- mxMatrix(name = "cl", type = "Lower", nrow = 2, ncol = 2, labels = ClLabs, free=T, value=0.5)
pathEl <- mxMatrix(name = "el", type = "Lower", nrow = 2, ncol = 2, labels = ElLabs, free=T, value=0.5)

should really be

pathAl <- mxMatrix(name = "al", type = "Lower", nrow = nf, ncol = nf, labels = AlLabs, free=T, value=0.5)
pathCl <- mxMatrix(name = "cl", type = "Lower", nrow = nf, ncol = nf, labels = ClLabs, free=T, value=0.5)
pathEl <- mxMatrix(name = "el", type = "Lower", nrow = nf, ncol = nf, labels = ElLabs, free=T, value=0.5)

so that the matrices end up the same dimension as iSDI etc.

I note also that the standardization approach isn't ideal - and could cause problems with submodels if, e.g., one of the elements of V, which is used as a denominator in a division, is zero. This would happen if A C and E covariances (but not variances) were set to zero.

StandardizedA <- mxAlgebra(A/V, name="StandardizedA")
StandardizedC <- mxAlgebra(C/V, name="StandardizedC")
StandardizedE <- mxAlgebra(E/V, name="StandardizedE")

Log in or register to post comments
  Login

News