Hi all,

After attending the Workshop in March, we now tried to work with OpenMx for the first time without help... Not completely succesful yet. We only mildly adjusted a script provided at the workshop but it gives an error message: 'The observed data associated with the FIML objective in model 'MZ' does not contain dimnames'. Please find the script below, it seems to us that we did define dimnames (dimnames=selVars). What do we do wrong?

Many thanks in advance for your efforts.

Cathy & Marcel

require(OpenMx)

require(psych)

source("http://www.vipbg.vcu.edu/~vipbg/Tc24/GenEpiHelperFunctions.R")

# Prepare Data

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

# data(mean_schousbou_aus_mzm)

# describe(mean_schousbou_aus_mzm)

Vars <- 'bmi'

nv <- 1

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

ntv <- nv*2

mzData <- mean_schousbou_aus_mzm

dzData <- mean_schousbou_aus_dzm

# Print Descriptive Statistics

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

describe(mean_schousbou_aus_mzm)

colMeans(mean_schousbou_aus_mzm,na.rm=TRUE)

cov(mzData,use="complete")

describe(dzData)

colMeans(dzData,na.rm=TRUE)

cov(dzData,use="complete")

# Fit ACE Model with RawData and Matrices Input

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

univACEModel <- mxModel("univACE",

mxModel("ACE",

# Matrices a, c, and e to store a, c, and e path coefficients

mxMatrix( type="Lower", nrow=nv, ncol=nv, free=TRUE, values=.6, label="a11", name="a" ), #X

mxMatrix( type="Lower", nrow=nv, ncol=nv, free=TRUE, values=.6, label="c11", name="c" ), #Y

mxMatrix( type="Lower", nrow=nv, ncol=nv, free=TRUE, values=.6, label="e11", name="e" ), #Z

# 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(IV)), name="iSD"),

# 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

## Note that the rest of the mxModel statements do not change for bivariate/multivariate case

# Matrix & Algebra for expected means vector

mxMatrix( type="Full", nrow=1, ncol=nv, free=TRUE, values= 20, label="mean", name="Mean"),

mxAlgebra( expression= cbind(Mean,Mean), name="expMean"),

# 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)

),

mxModel("DZ",

mxData( observed=dzData, type="raw" ),

mxFIMLObjective( covariance="ACE.expCovDZ", means="ACE.expMean", dimnames=selVars)

),

mxAlgebra( expression=MZ.objective + DZ.objective, name="minus2sumll" ),

mxAlgebraObjective("minus2sumll")

)

univACEFit <- mxRun(univACEModel)

univACESumm <- summary(univACEFit)

univACESumm

Hi Marcel,

The dimnames argument in the mxFIMLObjective() function specifies the row and column names of the expected covariance matrix and the expected means matrix. The error message is reporting that the observed data needs column names. Once the observed data has column names, then the observed data is matched to the covariance and means matrices.

OpenMx does not assume that the data columns are in the same order for the observed data and the expected matrices. This allows the observed data to exist in an arbitrary order, and it allows the observed data to contain a superset of the data used in the model fitting process.

the error might perhaps say

"The observed data "mzData" provided to mxData in model 'MZ' does not contain dimnames"

instead of

'The observed data associated with the FIML objective in model 'MZ' does not contain dimnames'.

??

the error might perhaps say

"The observed data "mzData" provided to mxData in model 'MZ' does not contain dimnames"

instead of

'The observed data associated with the FIML objective in model 'MZ' does not contain dimnames'.

??

Or perhaps, to disambiguate between the cases where dimnames are missing in the data, and those where dimnames were not provided in the mxFIMLObjective() function call:

The observed data do not contain a variable with the name "bmi1" (or whatever has been found not to have a match - best of all a list of those variables not matched, perhaps with "none of the variables in the list provided to the mxFIMLObjective function could be found in the observed data".

There is an opportunity to be really nice to the user here: if there are no dimnames for the data, tell them that. If there are no dimnames in the Objective function call, tell them that. If both exist, but some of the variables in the Objective function cannot be found in the data, tell them what is in BOTH lists, and which variable(s) cannot be found. Possibly, this last list might only have one element in it, if the error checking and everything else stops when the first missing variable name is detected, in which case "Variable bmi1, at least, could not be found in the data. I can only analyze variables which are in both the data and the model, so please try to get these things to agree, and don't forget that they are of course case-sensitive.

Many thanks for your help! It works! We are happy!

Cathy and Marcel

One more question though; with column names added to the original data the script runs. We subsequently fitted ACE and AE models to the data, which provided us with path estimates and variance components for a, c and e. Unfortunately I can't recall how to get OpenMx to provide us with confidence intervals for the variance components as well. Could you help us out?

Thanks again for your help.

Cathy and Marcel

Yes, see thread http://openmx.psyc.virginia.edu/thread/153#comment-1054 - should do what you want.

Thanks. We added these commands into our script but we now get the following error:

"In model 'twinACE' NPSOL returned a non-zero status code 6. model does not satisfy the first-order optimality conditions to the required accuracy, and no improved point for the merit function could be found during the final linesearch (Mx status RED)".

We understand this suggests the model somehow doesn't fit the data properly. However, the model would run to calculate estimates of A, C, E, just not to calculate confidence intervals. Any suggestions?

Thanks again for your help,

Cathy and Marcel

Actually, status code 6 does not mean that the model does not fit the data. This code can occur whether the model fits brilliantly well or extremely badly or anywhere in between. What code 6 means is that the optimizer is not comfortable that it is at a global minimum. The usual ways to check this out are: i) to try optimization from different starting values to see if you get the same answer; ii) try optimizing again from the solution thus far. This second part is quite easy to set up:

run1 <- mxRun(mymodel)

run2 <- mxRun(run1)

If the solution improves (-2lnL moves in a negative direction) then it may be an idea to try

run3 <- mxRun(run2)

You may find similar suggestions in other threads - try searching for IFAIL or status code 6.

With a little fiddling with R it should be fairly easy to supply random numbers for starting values of parameters, and to try this several times in a loop. If you do this, please be aware that this procedure can start optimization in some really bad places, so many of the solutions may be worse than the one found originally. It is only when a solution which fits better that optimization has been found lacking. In my experience, more than half code 6's are spurious - the optimizer had found the global solution, but it was not happy that the gradients etc. were flat enough.