You are here

Home › Forums › OpenMx 2.0 Discussion

Nested models in Multivariate Cholesky script - Error message for AE

3 posts / 0 new
Log in or register to post comments
Last post
Sun, 06/14/2015 - 15:08
#1
Samaneh2013's picture
Samaneh2013
Offline
Joined: 06/19/2013 - 17:00
Nested models in Multivariate Cholesky script - Error message for AE

Hi,

I have been running longitudinal multivariate cholesky decomposition. This seems o run ok but the problem is comparing nested models to ACE model and the Saturated model.
The error message is:
Error: The following error occurred while evaluating the subexpression 'solve(sqrt(CholAE.I * CholAE.C))' during the evaluation of 'iSDc' in model 'CholAE' : Lapack routine dgesv: system is exactly singular: U[1,1] = 0

The code:

-----------------------------------------------------------------------------

Cholesky Decomposition AE Model

------------------------------------------------------------------------------

CholAeModel <- mxModel( CholAceFit, name="CholAE" )
CholAeModel <- omxSetParameters( CholAeModel, labels=cLabs, free=FALSE, values=0 )
CholAeFit <- mxRun(CholAeModel)
mxCompare(CholAceFit, CholAeFit)

Top
Mon, 06/15/2015 - 18:37
neale's picture
neale
Offline
Joined: 07/31/2009 - 15:14
Matrix is singular

The issue here is that the CholAE.I *CholAE.C matrix cannot be inverted, because it has a zero in the first element (U[1,1]). This error message is, I agree, not very transparent. It seems to stem from the attempted standardization of the A C and E matrices so that confidence intervals can be calculated on the proportions of variance and the genetic/C/E correlations. These are reasonable things to do, but lead to problems on fitting submodels because inverting a matrix with zeroes in the diagonal (as is the case with C when it has had its elements fixed to zero) is not going to work - it involves dividing by zero which is a Bad Thing.

To solve your problem you could wipe out the standardization parts of the model

# Algebra to comopute Correlations
invSDa<- mxAlgebra( expression=solve(sqrt(I*A)), name="iSDa")
invSDc<- mxAlgebra( expression=solve(sqrt(I*C)), name="iSDc")
invSDe<- mxAlgebra( expression=solve(sqrt(I*E)), name="iSDe")
Rph <- mxAlgebra( expression=iSD  %&% V , name="Phcor")
Rg  <- mxAlgebra( expression=iSDa %&% A , name="Acor")
Rc  <- mxAlgebra( expression=iSDc %&% C , name="Ccor")
Re  <- mxAlgebra( expression=iSDe %&% E , name="Ecor")

and subsequent references to these objects in the model:

# Print the rG, rC, rE matrices (above, in the summary output you will find them with 95% CI)
CholAceFit$Acor
CholAceFit$Ccor
CholAceFit$Ecor
CholAceFit$Phcor
CholAceFit$scalar

Should you want them to inspect, you could use mxEval() to obtain them later, e.g.,

mxEval( expression=iSDa %&% A , model=CholAceFit, compute=T)

If, however, you wanted to see the confidence intervals (they've been commented out in your script anyway so I somewhat doubt it) then you could add those algebras expected to work. So A,C and E for the ACE model, or just A and E for the AE model with something like

CholAeModel   <- mxModel( CholAceFit, name="CholAE" )
CholAeModel   <- omxSetParameters( CholAeModel, labels= cLabs, free = FALSE, values=0 )
CholAeModel   <- mxModel(CholAeModel, 
                              mxAlgebra( expression=solve(sqrt(I*A)), name="iSDa"),
                              mxAlgebra( expression=solve(sqrt(I*E)), name="iSDe"),
                              mxAlgebra( expression=iSDa %&% A , name="Acor"),
                              mxAlgebra( expression=iSDe %&% E , name="Ecor"),
                              mxCI (c("Acor", "Ecor")
                                        )
CholAeFit     <- mxRun(CholAeModel, intervals=T)
Top
Tue, 06/16/2015 - 12:43
(Reply to #2)
Samaneh2013's picture
Samaneh2013
Offline
Joined: 06/19/2013 - 17:00
Problem resolved

Hi Neale,

Thank you very much for the detailed answer, this works now and the problem is resolved.

Top

Copyright © 2007-2024 The OpenMx Project