You are here

Confidence intervals for shared environmental correlation are (-1.00, 1.00)

13 posts / 0 new
Last post
JuliaJoplin's picture
Offline
Joined: 10/12/2014 - 16:44
Confidence intervals for shared environmental correlation are (-1.00, 1.00)

In a ACE bivariate Cholesky, any ideas why confidence intervals for a shared environmental correlation would be returned as (-1.00, 1.00)? My syntax is included below. The model converged without any errors.

bivACEModel <- mxModel("bivACE",
    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=.8, labels=c("a11","a21","a22"), name="a" ),
        mxMatrix( type="Lower", nrow=nv, ncol=nv, free=TRUE, values=.8, labels=c("c11","c21","c22"),name="c" ),
        mxMatrix( type="Lower", nrow=nv, ncol=nv, free=TRUE, values=.8, labels=c("e11","e21","e22"),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="iSD"),
        mxAlgebra(iSD %*% a, name="sta"),
        mxAlgebra(iSD %*% c, name="stc"),
        mxAlgebra(iSD %*% e, name="ste"),
        mxAlgebra(A/V, name="StandardizedA"),
        mxAlgebra(C/V, name="StandardizedC"),
        mxAlgebra(E/V, name="StandardizedE"),
        mxAlgebra(solve(sqrt(I*A)) %&% A, name="CorA"),
        mxAlgebra(solve(sqrt(I*C)) %&% C, name="CorC"),
        mxAlgebra(solve(sqrt(I*E)) %&% E, name="CorE"),
        mxAlgebra(solve(sqrt(I*V)) %&% V, name="CorP"),
   ## 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, 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="minus2sumloglikelihood" ),
    mxAlgebraObjective("minus2sumloglikelihood"),
    mxCI(c('ACE.sta', 'ACE.stc', 'ACE.ste')),
    mxCI(c('ACE.StandardizedA', 'ACE.StandardizedC', 'ACE.StandardizedE')),
    mxCI(c('ACE.CorA', 'ACE.CorC', 'ACE.CorE', 'ACE.CorP'))
)
bivACEFit <- mxRun(bivACEModel)
bivACESumm <- summary(bivACEFit)
bivACESumm
 
jpritikin's picture
Offline
Joined: 05/24/2012 - 00:35
version?

Please report the output of mxVersion()

JuliaJoplin's picture
Offline
Joined: 10/12/2014 - 16:44
So, I discovered that it was

So, I discovered that it was an old version of MX -- 1.4-3060-- running on an old version of R -- 3.00. Here are the confidence intervals for the syntax I ran above:

confidence intervals:
                              lbound    estimate       ubound
ACE.sta[1,1]           -6.896004e-01 -0.57868872 -0.456512579
ACE.sta[1,2]            0.000000e+00  0.00000000  0.000000000
ACE.sta[2,1]           -4.364658e-01 -0.22016369 -0.001328027
ACE.sta[2,2]           -4.997451e-01  0.34031281  0.499745056
ACE.stc[1,1]           -7.081026e-01 -0.62012270 -0.500041821
ACE.stc[1,2]            0.000000e+00  0.00000000  0.000000000
ACE.stc[2,1]           -3.084933e-01 -0.12736917  0.066596127
ACE.stc[2,2]           -4.546382e-01 -0.25391401  0.454638223
ACE.ste[1,1]           -5.612665e-01 -0.52968595 -0.499881726
ACE.ste[1,2]            0.000000e+00  0.00000000  0.000000000
ACE.ste[2,1]           -1.537588e-01 -0.09065141 -0.027862200
ACE.ste[2,2]           -9.033611e-01 -0.86417720 -0.824994808
ACE.StandardizedA[1,1]  2.084154e-01  0.33488064  0.475600204
ACE.StandardizedA[1,2]  1.800016e-03  0.50079589  1.020476196
ACE.StandardizedA[2,1]  1.800016e-03  0.50079589  1.020476196
ACE.StandardizedA[2,2]  1.620320e-03  0.16428486  0.304279511
ACE.StandardizedC[1,1]  2.500418e-01  0.38455216  0.501438707
ACE.StandardizedC[1,2] -1.585491e-01  0.31046452  0.753865370
ACE.StandardizedC[2,1] -1.585491e-01  0.31046452  0.753865370
ACE.StandardizedC[2,2]  2.303712e-17  0.08069523  0.240496725
ACE.StandardizedE[1,1]  2.498817e-01  0.28056720  0.315021278
ACE.StandardizedE[1,2]  5.976999e-02  0.18873960  0.322220089
ACE.StandardizedE[2,1]  5.976999e-02  0.18873960  0.322220089
ACE.StandardizedE[2,2]  6.880836e-01  0.75501991  0.825169268
ACE.CorA[1,1]           1.000000e+00  1.00000000  1.000000000
ACE.CorA[1,2]           1.527377e-03  0.54318396  1.000000000
ACE.CorA[2,1]           1.528116e-03  0.54318396  1.000000000
ACE.CorA[2,2]           1.000000e+00  1.00000000  1.000000000
ACE.CorC[1,1]           1.000000e+00  1.00000000  1.000000000
ACE.CorC[1,2]          -1.000000e+00  0.44837395  1.000000000
ACE.CorC[2,1]          -1.000000e+00  0.44837395  1.000000000
ACE.CorC[2,2]           1.000000e+00  1.00000000  1.000000000
ACE.CorE[1,1]           1.000000e+00  1.00000000  1.000000000
ACE.CorE[1,2]           3.214664e-02  0.10432667  0.175575168
ACE.CorE[2,1]           3.214664e-02  0.10432667  0.175575167
ACE.CorE[2,2]           1.000000e+00  1.00000000  1.000000000
ACE.CorP[1,1]           1.000000e+00  1.00000000  1.000000000
ACE.CorP[1,2]           2.137637e-01  0.25440754  0.294307088
ACE.CorP[2,1]           2.137637e-01  0.25440754  0.294307089
ACE.CorP[2,2]           1.000000e+00  1.00000000  1.000000000

Note the confidence intervals I'm confused about:

ACE.CorC[1,2]          -1.000000e+00  0.44837395  1.000000000
ACE.CorC[2,1]          -1.000000e+00  0.44837395  1.000000000

I updated both Open Mx and R and got this:

confidence intervals:
                              lbound    estimate     ubound note
ACE.sta[1,1]           -6.895056e-01  0.57868882  0.6895018     
ACE.sta[2,1]            1.754875e-03  0.22016381  0.4358245     
ACE.sta[1,2]                      NA  0.00000000         NA  !!!
ACE.sta[2,2]           -4.997114e-01 -0.34036930  0.4997087     
ACE.stc[1,1]            5.002152e-01  0.62012279  0.7080538     
ACE.stc[2,1]           -6.621980e-02  0.12736973  0.3082047     
ACE.stc[1,2]                      NA  0.00000000         NA  !!!
ACE.stc[2,2]           -4.546012e-01 -0.25385084  0.4546005     
ACE.ste[1,1]            4.998826e-01  0.52968573  0.5612627     
ACE.ste[2,1]            2.789288e-02  0.09065103  0.1537259     
ACE.ste[1,2]                      NA  0.00000000         NA  !!!
ACE.ste[2,2]           -9.031761e-01 -0.86417344 -0.8250037     
ACE.StandardizedA[1,1]  2.085394e-01  0.33488075  0.4752662     
ACE.StandardizedA[2,1]  3.922934e-03  0.50079578  1.0180140     
ACE.StandardizedA[1,2]  3.922934e-03  0.50079578  1.0180140     
ACE.StandardizedA[2,2]  1.624890e-03  0.16432336  0.3042418     
ACE.StandardizedC[1,1]  2.503027e-01  0.38455228  0.5013411     
ACE.StandardizedC[2,1] -1.565228e-01  0.31046565  0.7522242     
ACE.StandardizedC[1,2] -1.565228e-01  0.31046565  0.7522242     
ACE.StandardizedC[2,2]  6.078903e-38  0.08066330  0.2404413     
ACE.StandardizedE[1,1]  2.498834e-01  0.28056697  0.3150154     
ACE.StandardizedE[2,1]  5.990391e-02  0.18873857  0.3220464     
ACE.StandardizedE[1,2]  5.990391e-02  0.18873857  0.3220464     
ACE.StandardizedE[2,2]  6.881153e-01  0.75501334  0.8249501     
ACE.CorA[1,1]                     NA  1.00000000         NA  !!!
ACE.CorA[2,1]           5.431206e-01  0.54312061  1.0000000  !!!
ACE.CorA[1,2]           5.431206e-01  0.54312061  1.0000000  !!!
ACE.CorA[2,2]                     NA  1.00000000         NA  !!!
ACE.CorC[1,1]                     NA  1.00000000         NA  !!!
ACE.CorC[2,1]          -1.741573e-02  0.44846468  1.0000000     
ACE.CorC[1,2]          -6.765430e-01  0.44846468  1.0000000     
ACE.CorC[2,2]                     NA  1.00000000         NA  !!!
ACE.CorE[1,1]                     NA  1.00000000         NA  !!!
ACE.CorE[2,1]           3.218793e-02  0.10432669  0.1755334     
ACE.CorE[1,2]           3.218793e-02  0.10432669  0.1755334     
ACE.CorE[2,2]                     NA  1.00000000         NA  !!!
ACE.CorP[1,1]                     NA  1.00000000         NA  !!!
ACE.CorP[2,1]           2.137777e-01  0.25440777  0.2942948     
ACE.CorP[1,2]           2.137777e-01  0.25440777  0.2942948     
ACE.CorP[2,2]                     NA  1.00000000         NA  !!!

However, the confidence intervals are still strange for the shared environmental correlation.

ACE.CorC[2,1]          -1.741573e-02  0.44846468  1.0000000     
ACE.CorC[1,2]          -6.765430e-01  0.44846468  1.0000000   

Thanks for your help.

jpritikin's picture
Offline
Joined: 05/24/2012 - 00:35
CIs

You still haven't reported mxVersion() so we know which version you are using.

Version 2.2.4 has new diagnostics for investigation of profile CIs. I believe you can access them with model$compute$steps$CI$output$detail

JuliaJoplin's picture
Offline
Joined: 10/12/2014 - 16:44
Oops -- the first set of

Oops -- the first set of results I list above come from "1.4-3060". But I downloaded the latest OpenMX version for the second set of results that I include above, so I will try running this model$compute$steps$CI$output$detail when I get home tonight. Thank you for your help!

jpritikin's picture
Offline
Joined: 05/24/2012 - 00:35
syntax

Hm, on second thought, that syntax won't work for OpenMx 2.2.4. You'll need something like model$compute$steps[[2]]$output$detail

See if you can find the output slot of the mxComputeConfidenceInterval object.

AdminRobK's picture
Offline
Joined: 01/24/2014 - 12:15
BTW, what you'll see from

BTW, what you'll see from model$compute$steps[[2]]$output$detail is a data table. Each row corresponds to an attempt to find a confidence limit. The first column is the name of the matrix element for which the confidence limit is being calculated (e.g., ACE.StandardizedC[2,2]). The second column tells you the value of the limit that was found. The third column tells you whether it's a lower limit (1) or an upper limit (0). The fourth tells you the value of the fitfunction at the limit. The rest of the columns tell you the values of the other parameters at the limit.

JuliaJoplin's picture
Offline
Joined: 10/12/2014 - 16:44
Still having some problems with confidence intervals...

Hi all. I really appreciated the helpful comments on this thread and was satisfied that the very small C for the second variable explained this issue. However, I am continuing to have some issues with confidence intervals now that I have transitioned to OpenMX 2.2.4. I'm doing a similar analysis as listed above, only it's a scalar model, so raw phenotypic variance is allowed to vary across males and females.

bivACEModel <- mxModel("bivACE",
    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=.8, labels=c("am11","am21","am22"), name="am" ),
        mxMatrix( type="Lower", nrow=nv, ncol=nv, free=TRUE, values=.8, labels=c("cm11","cm21","cm22"),name="cm" ),
        mxMatrix( type="Lower", nrow=nv, ncol=nv, free=TRUE, values=.8, labels=c("em11","em21","em22"),name="em" ),
        mxMatrix( type="Lower", nrow=nv, ncol=nv, free=TRUE, values=.8, labels=c("am11","am21","am22"), name="af" ),
        mxMatrix( type="Lower", nrow=nv, ncol=nv, free=TRUE, values=.8, labels=c("cm11","cm21","cm22"),name="cf" ),
        mxMatrix( type="Lower", nrow=nv, ncol=nv, free=TRUE, values=.8, labels=c("em11","em21","em22"),name="ef" ),
        mxMatrix(type="Diag", 4, 4, labels=rep("k",4), free=TRUE, values=1, name="scalar"),
    # Matrices A, C, and E compute variance components
        mxAlgebra( expression=am %*% t(am), name="Am" ),
        mxAlgebra( expression=cm %*% t(cm), name="Cm" ),
        mxAlgebra( expression=em %*% t(em), name="Em" ),
        mxAlgebra( expression=af %*% t(af), name="Af" ),
        mxAlgebra( expression=cf %*% t(cf), name="Cf" ),
        mxAlgebra( expression=ef %*% t(ef), name="Ef" ),
    # Algebra to compute total variances and standard deviations (diagonal only)
        mxAlgebra( expression=Am+Cm+Em, name="Vm" ),
        mxAlgebra( expression=Af+Cf+Ef, name="Vf" ),
        mxMatrix( type="Iden", nrow=nv, ncol=nv, name="I"),
        mxAlgebra( expression=solve(sqrt(I*Vm)), name="iSDm"),
        mxAlgebra( expression=solve(sqrt(I*Vf)), name="iSDf"),
        mxAlgebra(iSDm %*% am, name="stam"),
        mxAlgebra(iSDm %*% cm, name="stcm"),
        mxAlgebra(iSDm %*% em, name="stem"),
        mxAlgebra(iSDf %*% af, name="staf"),
        mxAlgebra(iSDf %*% cf, name="stcf"),
        mxAlgebra(iSDf %*% ef, name="stef"),
        mxAlgebra(Am/Vm, name="StandardizedAm"),
        mxAlgebra(Cm/Vm, name="StandardizedCm"),
        mxAlgebra(Em/Vm, name="StandardizedEm"),
        mxAlgebra(Af/Vf, name="StandardizedAf"),
        mxAlgebra(Cf/Vf, name="StandardizedCf"),
        mxAlgebra(Ef/Vf, name="StandardizedEf"),
        mxAlgebra(solve(sqrt(I*Am)) %&% Am, name="CorAm"),
        mxAlgebra(solve(sqrt(I*Cm)) %&% Cm, name="CorCm"),
        mxAlgebra(solve(sqrt(I*Em)) %&% Em, name="CorEm"),
        mxAlgebra(solve(sqrt(I*Am)) %&% Af, name="CorAf"),
        mxAlgebra(solve(sqrt(I*Cf)) %&% Cf, name="CorCf"),
        mxAlgebra(solve(sqrt(I*Ef)) %&% Ef, name="CorEf"),
        mxAlgebra(solve(sqrt(I*Vm)) %&% Vm, name="CorPm"),
        mxAlgebra(solve(sqrt(I*Vf)) %&% Vf, name="CorPf"),
       ## Note that the rest of the mxModel statements do not change for bivariate/multivariate case
    # Matrix & Algebra for expected means vector
        mxMatrix("Full", nrow=1, ncol=nv, free=TRUE, values= .6, label="meanM", name="Mm"), 
        mxMatrix("Full", nrow=1, ncol=nv, free=TRUE, values= .6, label="meanF", name="Mf"),
        mxAlgebra(cbind(Mm,Mm), name="expMeanM"),
        mxAlgebra(cbind(Mf,Mf), name="expMeanF"),
# Expected covariance matrix in MZ, for males and females
mxAlgebra( scalar %*% (rbind( cbind(Am+Cm+Em , Am+Cm),
                                cbind(Am+Cm    , Am+Cm+Em))), name="expCovMZm" ),
 
mxAlgebra(
rbind(
cbind(Af+Cf+Ef, Af+Cf),
cbind(Af+Cf, Af+Cf+Ef)
),
name="expCovMZf"
),
 
# Expected covariance matrix in DZ, for males and females
mxAlgebra(  scalar %*% (rbind( cbind(Am+Cm+Em    , 0.5%x%Am+Cm),
                                cbind(0.5%x%Am+Cm , Am+Cm+Em))), name="expCovDZm" ),
 
mxAlgebra(
rbind(
cbind(Af+Cf+Ef,   (0.5 %x% Af)+Cf),
cbind((0.5 %x% Af)+Cf, Af+Cf+Ef)
),
name="expCovDZf"
)
),       
 
mxModel("MZm",
mxData(observed=mzmData , type="raw"),
mxFIMLObjective(covariance="ACE.expCovMZm", means="ACE.expMeanM", dimnames=selVars)
),
 
mxModel("DZm",
mxData(observed=dzmData , type="raw"),
mxFIMLObjective(covariance="ACE.expCovDZm", means="ACE.expMeanM", dimnames=selVars)
),
 
mxModel("MZf",
mxData(observed=mzfData , type="raw"),
mxFIMLObjective(covariance="ACE.expCovMZf", means="ACE.expMeanF", dimnames=selVars)
),
 
mxModel("DZf",
mxData(observed=dzfData , type="raw"),
mxFIMLObjective(covariance="ACE.expCovDZf", means="ACE.expMeanF", dimnames=selVars)
),
 
mxAlgebra(MZm.objective + DZm.objective + MZf.objective + DZf.objective, name="minus2sumll"),
mxAlgebraObjective("minus2sumll"),
 
    mxCI(c('ACE.stam', 'ACE.stcm', 'ACE.stem','ACE.staf', 'ACE.stcf', 'ACE.stef')),
    mxCI(c('ACE.StandardizedAm', 'ACE.StandardizedCm', 'ACE.StandardizedEm', 'ACE.StandardizedAf', 'ACE.StandardizedCf', 'ACE.StandardizedEf')),
    mxCI(c('ACE.CorAm', 'ACE.CorCm', 'ACE.CorEm', 'ACE.CorAf', 'ACE.CorCf', 'ACE.CorEf', 'ACE.CorPf', 'ACE.CorPm'))
)
bivACEFit <- mxRun(bivACEModel, intervals=T)
bivACESumm <- summary(bivACEFit)
bivACESumm

When I run this model (it converges, I get the OpenMX 2.2.4 message about changing mxFIMLObjective), this is what I get:

free parameters:
    name     matrix row col    Estimate  Std.Error A
1   am11     ACE.am   1   1 0.528644361 0.05390055  
2   am21     ACE.am   2   1 0.326612961 0.08777112  
3   am22     ACE.am   2   2 0.474737615 0.08227078  
4   cm11     ACE.cm   1   1 0.565182277 0.05006764  
5   cm21     ACE.cm   2   1 0.187113198 0.07721006  
6   cm22     ACE.cm   2   2 0.164046483 0.19358661  
7   em11     ACE.em   1   1 0.484450835 0.01340149  
8   em21     ACE.em   2   1 0.078741665 0.02363708  
9   em22     ACE.em   2   2 0.612407656 0.01728468  
10     k ACE.scalar   1   1 1.246218021 0.05237848  
11 meanM     ACE.Mm   1   1 0.001183142 0.03170260  
12 meanF     ACE.Mf   1   1 0.004545240 0.02652372  
 
confidence intervals:
                               lbound   estimate    ubound note
ACE.stam[1,1]           -0.6896907881 0.57901225 0.6896905     
ACE.stam[2,1]           -0.5679600624 0.37095527 0.5681496     
ACE.stam[1,2]                      NA 0.00000000        NA  !!!
ACE.stam[2,2]           -0.6266624793 0.53918993 0.6266620     
ACE.stcm[1,1]            0.4990879855 0.61903141 0.7070824     
ACE.stcm[2,1]            0.0260882338 0.21251645 0.3781078     
ACE.stcm[1,2]                      NA 0.00000000        NA  !!!
ACE.stcm[2,2]           -0.4358926357 0.18631810 0.4358972     
ACE.stem[1,1]            0.5007497566 0.53060808 0.5622438     
ACE.stem[2,1]            0.0370386511 0.08943195 0.1423264     
ACE.stem[1,2]                      NA 0.00000000        NA  !!!
ACE.stem[2,2]            0.6586055309 0.69555062 0.7339534     
ACE.staf[1,1]           -0.6896907881 0.57901225 0.6896905     
ACE.staf[2,1]           -0.5679600624 0.37095527 0.5681496     
ACE.staf[1,2]                      NA 0.00000000        NA  !!!
ACE.staf[2,2]           -0.6266624793 0.53918993 0.6266620     
ACE.stcf[1,1]            0.4990879855 0.61903141 0.7070824     
ACE.stcf[2,1]            0.0260882338 0.21251645 0.3781078     
ACE.stcf[1,2]                      NA 0.00000000        NA  !!!
ACE.stcf[2,2]           -0.4358926357 0.18631810 0.4358972     
ACE.stef[1,1]            0.5007497566 0.53060808 0.5622438     
ACE.stef[2,1]            0.0370386511 0.08943195 0.1423264     
ACE.stef[1,2]                      NA 0.00000000        NA  !!!
ACE.stef[2,2]            0.6586055309 0.69555062 0.7339534     
ACE.StandardizedAm[1,1]  0.2088288309 0.33525519 0.4756665     
ACE.StandardizedAm[2,1]  0.2476788423 0.54542966 0.8643968     
ACE.StandardizedAm[1,2]  0.2476788423 0.54542966 0.8643968     
ACE.StandardizedAm[2,2]  0.2365927972 0.42833359 0.5379178     
ACE.StandardizedCm[1,1]  0.2492653686 0.38319988 0.4999650     
ACE.StandardizedCm[2,1]  0.0381772694 0.33406785 0.6048082     
ACE.StandardizedCm[1,2]  0.0381772694 0.33406785 0.6048082     
ACE.StandardizedCm[2,2]  0.0006711864 0.07987768 0.2498342     
ACE.StandardizedEm[1,1]  0.2507500486 0.28154493 0.3161177     
ACE.StandardizedEm[2,1]  0.0504972567 0.12050249 0.1937924     
ACE.StandardizedEm[1,2]  0.0504972567 0.12050249 0.1937924     
ACE.StandardizedEm[2,2]  0.4408910293 0.49178873 0.5476654     
ACE.StandardizedAf[1,1]  0.2088288309 0.33525519 0.4756665     
ACE.StandardizedAf[2,1]  0.2476788423 0.54542966 0.8643968     
ACE.StandardizedAf[1,2]  0.2476788423 0.54542966 0.8643968     
ACE.StandardizedAf[2,2]  0.2365927972 0.42833359 0.5379178     
ACE.StandardizedCf[1,1]  0.2492653686 0.38319988 0.4999650     
ACE.StandardizedCf[2,1]  0.0381772694 0.33406785 0.6048082     
ACE.StandardizedCf[1,2]  0.0381772694 0.33406785 0.6048082     
ACE.StandardizedCf[2,2]  0.0006711864 0.07987768 0.2498342     
ACE.StandardizedEf[1,1]  0.2507500486 0.28154493 0.3161177     
ACE.StandardizedEf[2,1]  0.0504972567 0.12050249 0.1937924     
ACE.StandardizedEf[1,2]  0.0504972567 0.12050249 0.1937924     
ACE.StandardizedEf[2,2]  0.4408910293 0.49178873 0.5476654     
ACE.CorAm[1,1]                     NA 1.00000000        NA  !!!
ACE.CorAm[2,1]           0.5668008284 0.56680083 0.8572289  !!!
ACE.CorAm[1,2]           0.5668008284 0.56680083 0.8572289  !!!
ACE.CorAm[2,2]                     NA 1.00000000        NA  !!!
ACE.CorCm[1,1]                     NA 1.00000000        NA  !!!
ACE.CorCm[2,1]           0.7519342020 0.75193420 1.0000000  !!!
ACE.CorCm[1,2]           0.5593583429 0.75193420 1.0000000     
ACE.CorCm[2,2]                     NA 1.00000000        NA  !!!
ACE.CorEm[1,1]                     NA 1.00000000        NA  !!!
ACE.CorEm[2,1]           0.0530388967 0.12752738 0.2007366     
ACE.CorEm[1,2]           0.0530388967 0.12752738 0.2007366     
ACE.CorEm[2,2]                     NA 1.00000000        NA  !!!
ACE.CorAf[1,1]                     NA 1.00000000        NA  !!!
ACE.CorAf[2,1]           0.5668008284 0.56680083 0.8572289  !!!
ACE.CorAf[1,2]           0.5668008284 0.56680083 0.8572289  !!!
ACE.CorAf[2,2]                     NA 1.00000000        NA  !!!
ACE.CorCf[1,1]                     NA 1.00000000        NA  !!!
ACE.CorCf[2,1]           0.7519342020 0.75193420 1.0000000  !!!
ACE.CorCf[1,2]           0.5593583429 0.75193420 1.0000000     
ACE.CorCf[2,2]                     NA 1.00000000        NA  !!!
ACE.CorEf[1,1]                     NA 1.00000000        NA  !!!
ACE.CorEf[2,1]           0.0530388967 0.12752738 0.2007366     
ACE.CorEf[1,2]           0.0530388967 0.12752738 0.2007366     
ACE.CorEf[2,2]                     NA 1.00000000        NA  !!!
ACE.CorPf[1,1]                     NA 1.00000000        NA  !!!
ACE.CorPf[2,1]           0.3543287069 0.39379532 0.4319831     
ACE.CorPf[1,2]           0.3543287068 0.39379532 0.4319831     
ACE.CorPf[2,2]                     NA 1.00000000        NA  !!!
ACE.CorPm[1,1]                     NA 1.00000000        NA  !!!
ACE.CorPm[2,1]           0.3543287069 0.39379532 0.4319831     
ACE.CorPm[1,2]           0.3543287068 0.39379532 0.4319831     
ACE.CorPm[2,2]                     NA 1.00000000        NA  !!!

There seem to be a number of irregularities here. For example, look at the CIs around these selected path estimates.

ACE.stam[1,1]           -0.6896907881 0.57901225 0.6896905     
ACE.stam[2,1]           -0.5679600624 0.37095527 0.5681496     
ACE.stam[2,2]           -0.6266624793 0.53918993 0.6266620     
 
ACE.stcm[2,2]           -0.4358926357 0.18631810 0.4358972   

And then, for the genetic/shared environment correlations, also some strange results:

ACE.CorAm[2,1]           0.5668008284 0.56680083 0.8572289  !!!
 
ACE.CorCm[2,1]           0.7519342020 0.75193420 1.0000000  !!!
ACE.CorCm[1,2]           0.5593583429 0.75193420 1.0000000     

Unfortunately, when I use the code

model$compute$steps[[2]]$output$detail

I get "NULL." Any insights would be quite valuable for me -- I want to be able to report if the standardized paths are significant and have confidence intervals for the genetic and environmental correlations.

AdminRobK's picture
Offline
Joined: 01/24/2014 - 12:15
Unfortunately, when I use the
Unfortunately, when I use the code

model$compute$steps[[2]]$output$detail

I get "NULL."

Really? What do you see from str(bivACEFit$compute$steps) ?

There seem to be a number of irregularities here. For example, look at the CIs around these selected path estimates.

ACE.stam[1,1]           -0.6896907881 0.57901225 0.6896905     
ACE.stam[2,1]           -0.5679600624 0.37095527 0.5681496     
ACE.stam[2,2]           -0.6266624793 0.53918993 0.6266620     
 
ACE.stcm[2,2]           -0.4358926357 0.18631810 0.4358972 

There's nothing wrong here that jumps out at me...

And then, for the genetic/shared environment correlations, also some strange results:

ACE.CorAm[2,1]           0.5668008284 0.56680083 0.8572289  !!!
 
ACE.CorCm[2,1]           0.7519342020 0.75193420 1.0000000  !!!
ACE.CorCm[1,2]           0.5593583429 0.75193420 1.0000000     

Since elements [1,2] and [2,1] of 'CorCm' are supposed to be equal due to symmetry (right?), I think you can trust the CI provided for element [1,2]. The CI for [2,1] of 'CorAm' is a problem, though. There are a couple things you could try. For one, if you installed version 2.2.4 from virginia.edu (i.e., not from CRAN), you could run the line

mxOption(NULL, "Default optimizer", "NPSOL")

and then re-run your script. Alternately, if you only want to know if element [2,1] of 'CorAm' differs significantly from zero, you could make a new MxModel just like your old one, except with an MxConstraint that fixes the element to zero. Then, use mxCompare with your new and old model. Finally, if you really need to find the lower confidence limit, you could try running new MxModels that fix the element to various plausible values for the lower limit (again via an MxConstraint). This is rather labor-intensive, but you will have found an approximate lower limit if you find a candidate value that worsens the model's fit by about 3.84.

JuliaJoplin's picture
Offline
Joined: 10/12/2014 - 16:44
Here it is

Thank you so very much for all of these excellent ideas. I believe I did install from CRAN, but I will try installing from virginia.edu and using the mxOption, and then I'll report back

With the first set of CIs that I pasted above, even though there were no errors, I was surprised to see that for these paths, the upper and lower limits have almost exactly the same absolute values. This just seemed very curious to me!

ACE.stam[1,1]           -0.6896907881 0.57901225 0.6896905     
ACE.stam[2,1]           -0.5679600624 0.37095527 0.5681496     
ACE.stam[2,2]           -0.6266624793 0.53918993 0.6266620     
 
ACE.stcm[2,2]           -0.4358926357 0.18631810 0.4358972

Also, I made a silly mistake when I was using model$compute$steps[[2]]$output$detail, and now I've got the output that you had previously described in an earlier post. What I don't understand is why the upper and lower limits are different than what appears in my model summary. For example, looking at the correlations, this is what appears from model$compute$steps[[2]]$output$detail (leaving out the values of other parameters at the limit for the sake o legibility):

99           ACE.CorAm[2,1]  0.8572289122     0 11230.44  
100          ACE.CorAm[2,1]  0.8179721008     1 11230.44 
101          ACE.CorAm[1,2]  0.8572289120     0 11230.44 
102          ACE.CorAm[1,2]  0.8204799792     1 11230.44
 
107          ACE.CorCm[2,1]  1.0000000000     0 11230.44  
108          ACE.CorCm[2,1]  0.8368893663     1 11230.44  
109          ACE.CorCm[1,2]  1.0000000000     0 11230.44  
110          ACE.CorCm[1,2]  0.5593583429     1 11230.44  
AdminRobK's picture
Offline
Joined: 01/24/2014 - 12:15
I believe I did install from
I believe I did install from CRAN, but I will try installing from virginia.edu and using the mxOption, and then I'll report back.

Instructions for installing from the virginia.edu repository are here.

With the first set of CIs that I pasted above, even though there were no errors, I was surprised to see that for these paths, the upper and lower limits have *almost* exactly the same absolute values.

Your MxAlgebra 'stam' is just the Cholesky factor 'am' with its diagonal elements divided by standard deviations. The signs of some of the parameters in the Cholesky factors are indeterminate. For instance, consider the [2,2] element of 'am', which you've labeled 'am22'. It only appears in the 2x2 variance component 'Am' as squared; specifically, the [2,2] element of 'Am' equals am21^2 + am22^2. The CIs where the upper and lower limit have approximately the same absolute magnitude should be correct, given the parameterization of this model. If they still make you uneasy, you could experiment with setting lower bounds of zero on some parameters, or just use a different parameterization.

Also, I made a silly mistake when I was using model$compute$steps[[2]]$output$detail, and now I've got the output that you had previously described in an earlier post. What I don't understand is why the upper and lower limits are different than what appears in my model summary.

Notice that the lower limits in the 'detail' output you highlighted are actually greater than the point estimate. I'm not completely sure, but I suspect that summary() doesn't report confidence limits that are so obviously wrong, and instead just uses the point estimate as the value for the problematic confidence limit, and flags it with "!!!" to alert the user that something is amiss.

AdminNeale's picture
Offline
Joined: 03/01/2013 - 14:09
Because the estimate of C 2 2 is very small

From your output, the estimate of C for the second variable is very small, and not significantly different from zero. Therefore there is very little information to estimate the C correlation with other variables, so I expect the result that you got is the right answer: rC could be fixed to any value between -1 and 1 without significantly worsening the fit.

ACE.StandardizedC[2,2]  2.303712e-17  0.08069523  0.240496725

This brings up the point that estimates of correlations between latent factors are not distributed the same as correlations between variables. Their precision depends on estimates of other parameters in the model. It is why two-step analyses (estimate say a genetic correlation matrix and then plug that into factor analysis to explore its structure) can only be exploratory. Standard tests of, e.g., the number of factors etc. would be misleading. The single-step approach, which you are using, is the right way to go to test such hypotheses.

JuliaJoplin's picture
Offline
Joined: 10/12/2014 - 16:44
Ah -- this makes sense. Thank

Ah -- this makes sense. Thank you for this helpful answer. In terms of reporting these results in a coherent way, then, I could include the (-1, 1) confidence interval for the shared environmental correlation and then demonstrate that the model fits as well by eliminating the c21 path?