You are here

OpenMx sensitivity to start values in a,c,e?

5 posts / 0 new
Last post
tbates's picture
Offline
Joined: 07/31/2009 - 14:25
OpenMx sensitivity to start values in a,c,e?

running a tri-bivariate ACE model, and the results can be all over the show (with no convergence warnings), depending on the start values: for instance C wandering from .998 to .054 for one variable.
Any thoughts?
Here are the results, and script fragment

Results starting at a, c, and e at .6

               A1     A2    A3    C1    C2    C3    E1    E2    E3
var1     0.284     NA    NA 0.282    NA    NA 0.434    NA    NA
var2    -0.292  0.027    NA 0.998 0.434    NA 0.294 0.539    NA
var3     0.224 -0.028 0.206 0.366 0.458 0.196 0.411 0.570 0.598
[1] 7550.711

Now, starting at a=.6, c=.3, e=.4

              A1    A2    A3     C1     C2    C3    E1    E2    E3
var1     0.560    NA    NA  0.006     NA    NA 0.434    NA    NA
var2     0.652 0.343    NA  0.054  0.118    NA 0.294 0.539    NA
var3     0.636 0.477 0.277 -0.047 -0.047 0.124 0.411 0.570 0.598
[1] 7550.711

relevant fragment

<

pre>
share = mxModel("share",
mxMatrix("Lower", nrow=nVar, ncol=nVar, free=TRUE, values=.6, name="a", byrow=TRUE), # Additive genetic path coefficient
mxMatrix("Lower", nrow=nVar, ncol=nVar, free=TRUE, values=.6, name="c", byrow=TRUE), # Common environmental path coefficient
mxMatrix("Lower", nrow=nVar, ncol=nVar, free=TRUE, values=.6, name="e", byrow=TRUE), # Unique environmental path coefficient
# Multiply by each path coefficient by its inverse to get variance component
mxAlgebra(a %% t(a), name="A"), # additive genetic variance
mxAlgebra(c %
% t(c), name="C"), # common environmental variance
mxAlgebra(e %*% t(e), name="E") # unique environmental variance
)

mzModel = mxModel(share, name="MZ",
expMZMeans,
# MZ ACE algebra
mxAlgebra(rbind (cbind(A+C+E , A+C),
cbind(A+C , A+C+E)), dimnames = list(selVars, selVars), name="expCov"),
mxData(mzData, type="raw"),
mxFIMLObjective("expCov", "expMean")
)

dzModel = mxModel(share, name="DZ",
expDZMeans,
mxAlgebra(rbind (cbind(A+C+E , .5%x%A+C),
cbind(.5%x%A+C, A+C+E) ), dimnames = list(selVars, selVars), name="expCov"),
mxData(dzData, type="raw"),
mxFIMLObjective("expCov", "expMean")
)

Build and Run ACE model

twinACEModel = mxModel("twinACE",mzModel, dzModel, mxAlgebra(MZ.objective + DZ.objective, name="twin"), mxAlgebraObjective("twin") )
twinACEFit = mxRun(twinACEModel)

tbates's picture
Offline
Joined: 07/31/2009 - 14:25
Interestingly, rewriting the

Interestingly, rewriting the script so that it uses three-models so that the MZ and DZ sub-models both do their algebra off a single set of ace and ACE matrices and algebras contained in the third shared model eliminates this instability, with the same solution reliably emerging from diverse starts.

Steve's picture
Offline
Joined: 07/30/2009 - 14:03
Hmmm.... That is interesting.

Hmmm.... That is interesting. We're trying to do the consolidation of the algebra in a smart way. It must be eliminating either some multiplication or a place where your first model had two parameters plus a constraint. I'd need to work through the model step by step and see where the two sets of convariance algebra differ.

neale's picture
Offline
Joined: 07/31/2009 - 15:14
The output fragment you

The output fragment you provide indicates strongly that the model is not identified. Different parameter estimates emerge, but the -2lnL is the same. I expect that somehow the information from one group is not being used in the model-fitting.

tbates's picture
Offline
Joined: 07/31/2009 - 14:25
my error was to assume that

my error was to assume that there was a global space of matrix names across the supermodel and all contained groups so matrix names acted like labels equating all matrices with the same name to be the same.

The fix here was to use the name of the mz group in the algebra of the dz group, linking the two algebras to the same underlying matrices.

mzModel = mxModel(name="mz", 
    mxMatrix("Lower", nVar, nVar, TRUE, .6, name="a", byrow=TRUE), # Additive genetic path coefficient
    mxMatrix("Lower", nVar, nVar, TRUE, .6, name="c", byrow=TRUE), # Common environmental path coefficient
    mxMatrix("Lower", nVar, nVar, TRUE, .6, name="e", byrow=TRUE), # Unique environmental path coefficient
    # Multiply by each path coefficient by its inverse to get variance component
    mxAlgebra(a %*% t(a), name="A"), # additive genetic variance
    mxAlgebra(c %*% t(c), name="C"), # common environmental variance
    mxAlgebra(e %*% t(e), name="E"), # unique environmental variance

    expMZMeans,
    # MZ ACE algebra
    mxAlgebra(rbind (cbind(A+C+E , A+C),
    cbind(A+C   , A+C+E)), dimnames = list(selVars, selVars), name="expCov"),
    mxData(mzData, type="raw"),
    mxFIMLObjective("expCov", "expMeanMZ")
)


dzModel = mxModel(name="dz", 
  expDZMeans,
  mxAlgebra(rbind (cbind(mz.A+mz.C+mz.E   , .5%x%mz.A+mz.C),
                   cbind(.5%x%mz.A+mz.C,    mz.A+mz.C+mz.E)  ), dimnames = list(selVars, selVars), name="expCov"),
    mxData(dzData, type="raw"),
    mxFIMLObjective("expCov", "expMeanDZ")
)

#Build and Run ACE model
model = mxModel("twinACE",mzModel, dzModel, mxAlgebra(mz.objective + dz.objective, name="twin"), mxAlgebraObjective("twin") )
fit = mxRun(model)

It will important to clearly explain the scope of labels and names to users, and that information from different groups has to be imported by using the object's full name, including the originating group.

Easy enough once you know :-)