Hi everyone, I am new here and learning OpenMX. I have a problem of getting confidence interval for standardized path coefficients. Could anyone please provide me the OpenMx language about this? I also provide my scripts below and hope someone can let me know what the problem is (the script about 'ciACE <- mxCI("VC[1:5, 16:30]")' ) or how to read the output about standardized CI.

# Load Data

nl <- read.csv(file="C:/R/TwinPairACE.csv", header=TRUE, sep=",")

# Recode Data for Analysis - Rescale variables to have variances around 1.0

nl$PC1 <- nl$PAIRC1/9.6

nl$PC2 <- nl$PAIRC2/9.6

nl$BDS1 <- nl$BDS1/2.4

nl$BDS2 <- nl$BDS2/2.4

nl$CHEXI_RE1 <- nl$CHEXI_RE1/18.5

nl$CHEXI_RE2 <- nl$CHEXI_RE2/18.5

nl$CVK1 <- nl$CVK1/9.7

nl$CVK2 <- nl$CVK2/9.7

nl$EVK1 <- nl$EVK1/12

nl$EVK2 <- nl$EVK2/12

# Select Variables for Analysis

vars <- c('PC','BDS','CHEXI_RE','CVK','EVK')

nv <- 5 # number of variables

ntv <- nv*2 # number of total variables

selVars <- paste(vars,c(rep(1,nv),rep(2,nv)),sep="")

# Select Covariates for Analysis

nl[,'Ageyr_Ind'] <- nl[,'Ageyr_Ind']/100

# nl <- nl[-which(is.na(nl$Ageyr_Ind)),] # missing value

covVars <- 'Ageyr_Ind'

nc <- nv # number of covariates

# Select Data for Analysis

mzData <- subset(nl, ZygTwin==1, c(selVars, covVars))

dzData <- subset(nl, ZygTwin==2, c(selVars, covVars))

# Set Starting Values

svMe <- c(7,5,5,1,1) # start value for means

svPa <- .4 # start value for path coefficient

svPaD <- vech(diag(svPa,nv,nv)) # start values for diagonal of covariance matrix

svPe <- .8 # start value for path coefficient for e

svPeD <- vech(diag(svPe,nv,nv)) # start values for diagonal of covariance matrix

lbPa <- .0001 # start value for lower bounds

lbPaD <- diag(lbPa,nv,nv) # lower bounds for diagonal of covariance matrix

lbPaD[lower.tri(lbPaD)] <- -10 # lower bounds for below diagonal elements

lbPaD[upper.tri(lbPaD)] <- NA # lower bounds for above diagonal elements

# Create Labels

labMe <- paste("mean",vars,sep="_")

labBe <- labFull("beta",nc,1)

# ------------------------------------------------------------------------------

# PREPARE MODEL

# ACE Model

# Create Matrices for Covariates and linear Regression Coefficients

defAge <- mxMatrix( type="Full", nrow=1, ncol=1, free=FALSE, labels=c("data.Ageyr_Ind"), name="Age" )

pathB <- mxMatrix( type="Full", nrow=nc, ncol=1, free=TRUE, values=.01, label=labBe, name="b" )

# Create Algebra for expected Mean Matrices

meanG <- mxMatrix( type="Full", nrow=1, ncol=ntv, free=TRUE, values=svMe, labels=labMe, name="meanG" )

expMean <- mxAlgebra( expression= meanG + cbind(t(b%*%Age),t(b%*%Age)), name="expMean" )

# Create Matrices for Path Coefficients

pathA <- mxMatrix( type="Lower", nrow=nv, ncol=nv, free=TRUE, values=svPaD, label=labLower("a",nv), lbound=lbPaD, name="a" )

pathC <- mxMatrix( type="Lower", nrow=nv, ncol=nv, free=TRUE, values=svPaD, label=labLower("c",nv), lbound=lbPaD, name="c" )

pathE <- mxMatrix( type="Lower", nrow=nv, ncol=nv, free=TRUE, values=svPeD, label=labLower("e",nv), lbound=lbPaD, name="e" )

# Create Algebra for Variance Components

covA <- mxAlgebra( expression=a %*% t(a), name="A" )

covC <- mxAlgebra( expression=c %*% t(c), name="C" )

covE <- mxAlgebra( expression=e %*% t(e), name="E" )

# Create Algebra for expected Variance/Covariance Matrices in MZ & DZ twins

covP <- mxAlgebra( expression= A+C+E, name="V" )

covMZ <- mxAlgebra( expression= A+C, name="cMZ" )

covDZ <- mxAlgebra( expression= 0.5%x%A+ C, name="cDZ" )

expCovMZ <- mxAlgebra( expression= rbind( cbind(V, cMZ), cbind(t(cMZ), V)), name="expCovMZ" )

expCovDZ <- mxAlgebra( expression= rbind( cbind(V, cDZ), cbind(t(cDZ), V)), name="expCovDZ" )

# Create Algebra for Standardization

matI <- mxMatrix( type="Iden", nrow=nv, ncol=nv, name="I")

invSD <- mxAlgebra( expression=solve(sqrt(I*V)), name="iSD")

# Calculate genetic and environmental correlations

corA <- mxAlgebra( expression=solve(sqrt(I*A))%&%A, name ="rA" )

corC <- mxAlgebra( expression=solve(sqrt(I*C))%&%C, name ="rC" )

corE <- mxAlgebra( expression=solve(sqrt(I*E))%&%E, name ="rE" )

# Create Data Objects for Multiple Groups

dataMZ <- mxData( observed=mzData, type="raw" )

dataDZ <- mxData( observed=dzData, type="raw" )

# Create Expectation Objects for Multiple Groups

expMZ <- mxExpectationNormal( covariance="expCovMZ", means="expMean", dimnames=selVars )

expDZ <- mxExpectationNormal( covariance="expCovDZ", means="expMean", dimnames=selVars )

funML <- mxFitFunctionML()

# Create Model Objects for Multiple Groups

pars <- list(pathB, meanG, matI, invSD,

pathA, pathC, pathE, covA, covC, covE, covP, corA, corC, corE)

defs <- list(defAge)

modelMZ <- mxModel( name="MZ", pars, defs, expMean, covMZ, expCovMZ, dataMZ, expMZ, funML )

modelDZ <- mxModel( name="DZ", pars, defs, expMean, covDZ, expCovDZ, dataDZ, expDZ, funML )

multi <- mxFitFunctionMultigroup( c("MZ","DZ") )

# Create Algebra for Variance Components

rowVC <- rep('VC',nv)

colVC <- rep(c('A','C','E','SA','SC','SE'),each=nv)

estVC <- mxAlgebra( expression=cbind(A,C,E,A/V,C/V,E/V), name="VC", dimnames=list(rowVC,colVC))

# Create Confidence Interval Objects

ciACE <- mxCI( "VC[1:5,16:30]" )

# Build Model with Confidence Intervals

modelACE <- mxModel( "mulACEEF", pars, modelMZ, modelDZ, multi, estVC, ciACE )

# Run ACE Model

fitACE <- mxRun( modelACE, intervals=T )

summary( fitACE, verbose=T )

# ACE Estimated Path Coefficients

matACEepaths <- c("a","c","e","iSD")

labACEepaths <- c("PathA","PathC","PathE","iSD")

formatOutputMatrices(fitACE, matACEepaths, labACEepaths, vars,4)

# ACE Standardized Path Coefficients (pre-multiplied by inverse of standard deviations)

matACEpaths <- c("iSD %*% a","iSD %*% c","iSD %*% e")

labACEpaths <- c("stPathA","stPathC","stPathE")

formatOutputMatrices(fitACE, matACEpaths, labACEpaths, vars,4)

Many thanks!

## summary

OpenMx's

`summary()`

function returns the SE (standard error) on free parameters. And a common basis for a 95% CI is to report 1.96 * SE + and - around parameter estimate.Two limitations of this are that parameter you want might be an algebra (for instance, in a twin model you might want a CI around

`A`

, but`A = a %*% t(a)`

. Second, the SE is always symmetric around the parameter, and is not the same as computing the parameter value at which fit declines by p= .05 in a confidence interval. But it is close enough for most purposes.## mxCI

OpenMx can compute confidence intervals around parameters. To do this, you add

`mxCI`

objects to the model. The help on mxCI has an example.trickyou not only have to add these to the model but, because they can take some time, request them run by adding`mxRun..., intervals = TRUE)`

to your mxRun call.These

`mxCI`

results appear in the`summary`

output in a CI section.## mxSE

You can also request SEs for derived parameters (including algebras) using

`mxSE`

These are very handy!## umx

Finally, the

`umx`

companion package (also on CRAN) supports CIs. most models, like`umxACE`

have CI objects added, andm1 = umxConfint(m1, run=TRUE) is designed to be a flexible way to report or request CIs from existing CI objects, or to add them and compute them.

You add new CIs via the

`parm`

parameter (like the generic R function`confint`

). Some models support "smart" CIs, where umx knows enough about the model type to monitor only the interesting free parameters (like`umxCP`

).Any suggestions for improvements? submit them to https://github.com/tbates/umx/issues

Many thanks for your quick reply! Sorry that I am a new hand of the OpenMx scripts.

If I want to build the multivariate Chelosky ACE model for 5 variables, is the script that 'mxCI ("VC[1:5, 16:30]")' correct? The estimated values of modelName.VC[1, 16] to VC[5, 30] should be the same as the standardized path coefficients? Actually they are not. Also it seems the values of modelName.VC[1,1] to VC[5, 15] are not the likelihood path coefficients.

If I used mxSE (labMe, model=modelACE), there is an error message:"Model does not have output and 'cov' argument is missing. I'm a doctor, not a bricklayer! as this model run with mxRun?" Anything wrong?

Anyway, I have been thankful for your reply. If you know my above problems, could you please let me know?

You are correct that

`VC`

does not contain path coefficients. As its name suggests, it contains variance-covariance components. If you want standardized path coefficients, you need to make an MxAlgebra that calculates them, and request CIs for that algebra. For instance:The value of argument

`model`

to`mxSE()`

has to be an MxModel objectthat has been run(as with`mxRun()`

). Try using`fitACE`

instead.Many thanks for your reply, and I have solved the problem!

I noticed the SE in the summary report. But we should not simply use these SE for calculating the CI for standardized paths, right? The SE is based on raw data and can be used for calculating the CI for the likelihood path coefficients; but the standardized paths used standardized scores (z scores), so the SE doesn't work for standardized CI. Is this right?

I'm not sure what you mean by "likelihood path coefficients", but you seem to have the right idea: it is not valid to use the SEs from

`summary()`

output to construct CIs for standardized path coefficients, because`summary()`

output only contains SEs for free parameters, and in your script, the standardized path coefficients are not free parameters, but are "algebras" i.e. functions of free parameters.You can use

`mxSE()`

to get SEs for algebras, but the profile-likelihood CIs are recommended due to their superior theoretical properties.