Confidence Intervals for Univariate Ordinal ACE Model with 2 Thresholds

4 replies [Last post]
kspoon's picture
Offline
Joined: 06/17/2011

I was wondering if there was a reason I can't seem to get confidence intervals for an ordinal ACE model with 2 thresholds.

I made sure to specify the mxCI for my standardized variance components call in my model and to put the "intervals=T" in my mxRun line. In my model fit $output$confidenceIntervals exists, but does not have any values.

Is this normal? If not, is there a way to create confidence intervals in the ordinal case? If it's just a silly code error, my code is below.

Any help would be appreciated!

univACEOrdModel <- mxModel("univACEOrd",
mxModel("ACE",
# Matrices a, c, and e to store a, c, and e path coefficients
mxMatrix( type="Full", nrow=nv, ncol=nv, free=TRUE, values=.6, label="a11", name="a" ),
mxMatrix( type="Full", nrow=nv, ncol=nv, free=TRUE, values=.6, label="c11", name="c" ),
mxMatrix( type="Full", nrow=nv, ncol=nv, free=TRUE, values=.6, label="e11", name="e" ),
# Matrices A, C, and E compute variance components
mxAlgebra( expression=a %*% t(a), name="A" ),
mxAlgebra( expression=c %*% t(c), name="C" ),
mxAlgebra( expression=e %*% t(e), name="E" ),
# Algebra to compute total variances and standard deviations (diagonal only)
mxAlgebra( expression=A+C+E, name="V" ),
mxMatrix( type="Iden", nrow=nv, ncol=nv, name="I"),
mxAlgebra( expression=solve(sqrt(I*V)), name="sd"),
mxAlgebra( expression=cbind(A/VP,C/VP,E/VP),name="stndVCs"),
# Calculate 95% CIs here
mxAlgebra(A+C+E,name="VP"),
##Yes, it's repetitive, but I was desperate and the above was exactly how it ran in my continuous univariate ACE model
mxCI(c("stndVCs")),
# Constraint on variance of ordinal variables
mxConstraint(V == I, name="Var1"),
# Matrix & Algebra for expected means vector
mxMatrix( type="Zero", nrow=1, ncol=nv, name="M" ),
mxAlgebra( expression= cbind(M,M), name="expMean" ),
mxMatrix( type="Full", nrow=2, ncol=nv, free=TRUE, values=c(0.8,1.2), label=c("thre1","thre2"), name="T" ),
mxAlgebra( expression= cbind(T,T), dimnames=list(c('th1','th2'),selVars), name="expThre" ),
# Algebra for expected variance/covariance matrix in MZ
mxAlgebra( expression= rbind ( cbind(A+C+E , A+C),
cbind(A+C , A+C+E)), name="expCovMZ" ),
# Algebra for expected variance/covariance matrix in DZ, note use of 0.5, converted to 1*1 matrix
mxAlgebra( expression= rbind ( cbind(A+C+E , 0.5%x%A+C),
cbind(0.5%x%A+C , A+C+E)), name="expCovDZ" )
),
mxModel("MZ",
mxData( observed=mzData, type="raw" ),
mxFIMLObjective( covariance="ACE.expCovMZ", means="ACE.expMean", dimnames=selVars, thresholds="ACE.expThre" )
),
mxModel("DZ",
mxData( observed=dzData, type="raw" ),
mxFIMLObjective( covariance="ACE.expCovDZ", means="ACE.expMean", dimnames=selVars, thresholds="ACE.expThre" )
),
mxAlgebra( expression=MZ.objective + DZ.objective, name="min2sumll" ),
mxAlgebraObjective("min2sumll")
)

univACEOrdFit <- mxRun(univACEOrdModel,intervals=T)

mdewey's picture
Offline
Joined: 01/21/2011
Best to avoid T (and F)

Kelly
It is better to spell out TRUE and FALSE as one day you will use a variable T (or F) and wonder why your previously working example stops working

mspiegel's picture
Offline
Joined: 07/31/2009
Along the same lines, you

Along the same lines, you should avoid naming you MxMatrix or MxAlgebra object "T" or "F" or "c". The model will execute, but should you try to use mxEval() to examine components of your model, then the model names will override the R function or variable names. For example,
mxEval(c(MZ.objective[1,1], DZ.objective[1,1]), univACEOrdFit)
will no longer work because mxEval() beleives that 'c' is a MxMatrix object, not the c() function in R.

mspiegel's picture
Offline
Joined: 07/31/2009
It was only recently in the

It was only recently in the OpenMx 1.1 beta release that we started support confidence intervals specified in the submodels. You can try downloading the beta and see if the script works. Another alternative would be to specify mxCI("ACE.stndVCs") in the container model.

kspoon's picture
Offline
Joined: 06/17/2011
Thanks so much for the quick

Thanks so much for the quick reply.

That's the problem with being lazy and adapting other people's scripts. I removed the unnecessary submodel and everything is working perfectly now.