Login

News

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

Posted on
No user picture. JuliaJoplin Joined: 10/12/2014
Forums
Behavioral Genetics Models
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

Replied on Wed, 06/03/2015 - 20:32
No user picture. JuliaJoplin Joined: 10/12/2014

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.

Log in or register to post comments
Replied on Wed, 06/03/2015 - 23:32
Picture of user. jpritikin Joined: 05/23/2012

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

Log in or register to post comments
Replied on Thu, 06/04/2015 - 10:13
No user picture. JuliaJoplin Joined: 10/12/2014

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!
Log in or register to post comments
Replied on Thu, 06/04/2015 - 10:52
Picture of user. jpritikin Joined: 05/23/2012

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.

Log in or register to post comments
Replied on Thu, 06/04/2015 - 10:59
Picture of user. AdminRobK Joined: 01/24/2014

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.
Log in or register to post comments
Replied on Tue, 07/14/2015 - 23:16
No user picture. JuliaJoplin Joined: 10/12/2014

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.

Log in or register to post comments
Replied on Wed, 07/15/2015 - 10:41
Picture of user. AdminRobK Joined: 01/24/2014

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.

Log in or register to post comments
Replied on Thu, 07/16/2015 - 19:55
No user picture. JuliaJoplin Joined: 10/12/2014

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

Log in or register to post comments
Replied on Fri, 07/17/2015 - 11:23
Picture of user. AdminRobK Joined: 01/24/2014

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.

Log in or register to post comments
Replied on Thu, 06/04/2015 - 07:45
Picture of user. AdminNeale Joined: 03/01/2013

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.

Log in or register to post comments
  Login

News