You are here

ACE: negative variance component from H.Maes’ published script and slides

2 posts / 0 new
Last post
lemire's picture
Offline
Joined: 12/13/2019 - 12:35
ACE: negative variance component from H.Maes’ published script and slides

Hi,
could you help provide some insights as to why running the following script, which uses data(twinData), gives negative estimates for Vc in an ACE model?

https://ibg.colorado.edu/cdrom2020/maes/twinModeling/oneACEvc.R

These are the most relevant lines related to the model:

covA      <- mxMatrix( type="Symm", nrow=nv, ncol=nv, free=TRUE, values=svPa, label="VA11", name="VA" ) 
covC      <- mxMatrix( type="Symm", nrow=nv, ncol=nv, free=TRUE, values=svPa, label="VC11", name="VC" )
covE      <- mxMatrix( type="Symm", nrow=nv, ncol=nv, free=TRUE, values=svPa, label="VE11", name="VE" )
 
# Create Algebra for expected Variance/Covariance Matrices in MZ & DZ twins
covP      <- mxAlgebra( expression= VA+VC+VE, name="V" )
covMZ     <- mxAlgebra( expression= VA+VC, name="cMZ" )
covDZ     <- mxAlgebra( expression= 0.5%x%VA+ VC, 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" )

Running the script provides the following solution:

fitACE$output$confidenceIntervals
                      lbound    estimate        ubound
oneACEvc.US[1,1]  0.60785108  0.75506776  0.9254683654
oneACEvc.US[1,2] -0.29756250 -0.14471244 -0.0061102627
oneACEvc.US[1,3]  0.15055312  0.16935009  0.1913908419

These estimates are further reported on page 71 of the accompanying slides found at

https://ibg.colorado.edu/cdrom2020/maes/twinModeling/TwinModeling2020c.pdf

I don’t really understand why these results would make sense. Of course, just looking at the slides I don’t have the full narrative, which is why I am asking here.

If VA/VC/VE are known to be variances, wouldn't including lbound=0 in the mxMatrix functions be more appropriate? Was that just an oversight? The slides further discuss negative C in the next slide, so it doesn't appear it was.

Furthermore, if setting lbound=0 is appropriate, I noticed that the confidence intervals failed to be calculated:

                     lbound      estimate ubound
oneACEvc.US[1,1] 0.54832905 6.1730227e-01     NA
oneACEvc.US[1,2]         NA 1.3552200e-12     NA
oneACEvc.US[1,3]         NA 1.7304629e-01     NA

can this be fixed?

Thanks for your help,

Mathieu

AdminRobK's picture
Offline
Joined: 01/24/2014 - 12:15
boundary conditions

See this other post concerning the negative variance components. The matter of negative variance components has been discussed a number of times on these forums, so the OpenMx website's search feature may be useful for finding previous such discussion.

The Verhulst et al. paper linked in that other post addresses why some of the confidence limits aren't shown when you add the lbounds: boundary conditions. If there is one or more active bounds near the confidence-limit solution, the limit may be wrong in that it makes the CI too narrow. Therefore, OpenMx does not print it in the default summary() output. You can see a table of "CI details" if you pass argument verbose=TRUE to summary(), and I expect that table will show the confidence limits missing from the default output to be flagged as either "active box constraints" or "alpha level not reached".