You are here

Interpret >1 and negative results from variance parametrisation of twin model

8 posts / 0 new
Last post
JuanJMV's picture
Offline
Joined: 07/20/2016 - 13:13
Interpret >1 and negative results from variance parametrisation of twin model

Hi All,

I am using the new parametrization for a bivariate model, as suggested in "Type I Error Rates and Parameter Bias in Multivariate Behavioral Genetic Models".

The results (see below) include an rA correlation that is higher than one.

Am I doing something wrong or I have to interpret the rA correlation as 1?

The phenotypic correlation is .2. Is it be correct to say that is due to 1, -.96 and .2 for A, D and E respectively?

There are also negative values: there is no way to interpret these results, right?

Thank you so much in advance

$V
mxAlgebra 'V' 
$formula:  VA + VD + VE 
$result:
           [,1]       [,2]
[1,] 0.26015570 0.08651586
[2,] 0.08651586 0.64434885
dimnames: NULL
 
$iSD
mxAlgebra 'iSD' 
$formula:  solve(sqrt(I * V)) 
$result:
          [,1]      [,2]
[1,] 1.9605744 0.0000000
[2,] 0.0000000 1.2457746
dimnames: NULL
 
$rA
mxAlgebra 'rA' 
$formula:  solve(sqrt(I * VA)) %&% VA 
$result:
          [,1]      [,2]
[1,] 1.0000000 3.1877784
[2,] 3.1877784 1.0000000
dimnames: NULL
 
$rD
mxAlgebra 'rD' 
$formula:  solve(sqrt(I * VD)) %&% VD 
$result:
            [,1]        [,2]
[1,]  1.00000000 -0.95732896
[2,] -0.95732896  1.00000000
dimnames: NULL
 
$rE
mxAlgebra 'rE' 
$formula:  solve(sqrt(I * VE)) %&% VE 
$result:
           [,1]       [,2]
[1,] 1.00000000 0.20271077
[2,] 0.20271077 1.00000000
dimnames: NULL
 
$US
mxAlgebra 'US' 
$formula:  cbind(VA, VD, VE, VA/V, VD/V, VE/V) 
$result:
            VA          VA           VD           VD          VE          VE         SA          SA         SD          SD
US 0.036015032 0.130102084  0.053615331 -0.096971671 0.170525336 0.053385447 0.13843645 1.503794611  0.2060894 -1.12085427
US 0.130102084 0.046249649 -0.096971671  0.191371989 0.053385447 0.406727210 1.50379461 0.071777344 -1.1208543  0.29700059
           SE         SE
US 0.65547415 0.61705966
US 0.61705966 0.63122206

Hope that you are all well during this difficult situation

AdminRobK's picture
Offline
Joined: 01/24/2014 - 12:15
The simple answer is that,

The simple answer is that, under the direct-symmetric parameterization, the A, C, and E matrices are not necessarily covariance matrices, and therefore when standardized, are not necessarily interpretable as correlation matrices.

JuanJMV's picture
Offline
Joined: 07/20/2016 - 13:13
Thank you so much for your

Thank you so much for your response. Super useful as always!

In that case, what would be your recommendation for interpreting the etiological correlations?

Thank you so much for your patience.

With all good wishes

AdminNeale's picture
Offline
Joined: 03/01/2013 - 14:09
Big genetic correlation > 1

Hi Juan

The large genetic correlation, rg greater than 1, is probably just some sampling error around an estimate of rg that is close to 1. Note that the genetic variances are small so rg may have confidence intervals that overlap 1.0 - perhaps by a lot. You could try a constraint to keep the correlation within the -1 to 1 bound and see how much difference it makes to the -2lnL.

To interpret it, I would inspect the correlations in the data. Likely what is happening is that the MZ cross-correlation between the traits is a bit bigger than expected, given how big the within-variable across MZ correlations are. So if looking at a correlation matrix that goes T1V1 T1V2 T2V1 T2V2 for columns and rows, it may be that rT1V1,T2V2 and/or rT1V2,T2V1 are a bit bigger than rT1V1,rT2V1 and rT1V2,rT2V2 (check their geometric mean).

JuanJMV's picture
Offline
Joined: 07/20/2016 - 13:13
Dear professor Neale,

Dear professor Neale,

Thank you so much for your response.

I have tried to constraint the rG <1 but I cannot make the model to work.

Here the lines that I have used and the error message.

MZconstraint <- mxConstraint(rA[2,1] < 1)
MZconstraint2 <- mxConstraint(rA[1,2] < 1)
modelADECONS <- mxModel( fitADE, name="ADECONS", MZconstraint,MZconstraint2)
fitADECONS <- mxTryHard( modelADECONS, intervals=F )

Warning message:
In mxTryHard(modelADECONS, intervals = F) :
Error : The following error occurred while evaluating the subexpression 'solve(sqrt(ADECONS.I * ADECONS.VA))' during the evaluation of 'rA' in model 'ADECONS' : Lapack routine dgesv: system is exactly singular: U[1,1] = 0

Regarding the correlations in the data, I think rT1V1,T2V2 are bigger than rT1V1,rT2V1 so probably that is the problem. I have calculated the geometric mean and there are differences.

Here the outputs in case useful.

$expCovMZ
mxAlgebra 'expCovMZ'
$formula: rbind(cbind(V, cMZ), cbind(t(cMZ), V))
$result:
[,1] [,2] [,3] [,4]
[1,] 0.260155920 0.086515864 0.089630573 0.033130345
[2,] 0.086515864 0.644348796 0.033130345 0.237621696
[3,] 0.089630573 0.033130345 0.260155920 0.086515864
[4,] 0.033130345 0.237621696 0.086515864 0.644348796
dimnames: NULL

> mxEval(cov2cor(expCovMZ), fitADE$MZ, T)
[,1] [,2] [,3] [,4]
[1,] 1.000000000 0.211309187 0.344526362 0.080918643
[2,] 0.211309187 1.000000000 0.080918643 0.368778055
[3,] 0.344526362 0.080918643 1.000000000 0.211309187
[4,] 0.080918643 0.368778055 0.211309187 1.000000000

Thank you so much in advance

AdminNeale's picture
Offline
Joined: 03/01/2013 - 14:09
Only one constraint needed

Hi Juan

It is not necessary to add both constraints - they do the same thing because rA is symmetric. Just one is better (the two are perfectly correlated).

Second, VA for variable 1 is set to, starting at, or fixed to zero, so when it tries to calculate rA, (remembering a correlation is a covariance divided b.y the square root of the product of the variances). In the present case with U[1,1]=0, it is trying to divide by zero, which isn't going to work.

Third, I would make the constraint a bit more effective (prevent <-1 as well as >1) with

MZconstraint <- mxConstraint(abs(rA[2,1]) < 1)

But in the end, I think there's something wrong with your VA parameters.

JuanJMV's picture
Offline
Joined: 07/20/2016 - 13:13
Dear Michael,

Dear Michael,

Thank you so much for your response. I have tried using the lines that you suggested but I got the same error message.

Error : The following error occurred while evaluating the subexpression 'solve(sqrt(ADECONS.I * ADECONS.VA))' during the evaluation of 'rA' in model 'ADECONS' : Lapack routine dgesv: system is exactly singular: U[1,1] = 0

I have also tried changing the starting values for VA but It does not work either.

Here you can see the values for VA from the ADE model without constraints:

$VA
SymmMatrix 'VA'

$labels
[,1] [,2]
[1,] "VA11" "VA21"
[2,] "VA21" "VA22"

$values
[,1] [,2]
[1,] 0.036015727 0.130102805
[2,] 0.130102805 0.046250591

Could it be because VA[2,1] is higher than VA[1,1] and VA[2,2]?

Thank you so much in advance for your help, I am totally stuck with this.

With all good wishes.

JuanJMV's picture
Offline
Joined: 07/20/2016 - 13:13
Hi again,

Hi again,

I have been trying to fix it but I cannot make it work.

I have checked the starting values for VA but I am not sure why VA is set to, starting at, or fixed to zero.

Here the lines for VA:

svPa <- .2 # start value for path coefficient
covA <- mxMatrix( type="Symm", nrow=nv, ncol=nv, free=TRUE, values=valDiag(svPa,nv), label=labLower("VA",nv), name="VA" )

So, if I am not wrong, VA should start at .2. Am I doing something wrong?

I have also tried to find U[1,1] = 0 in the model but I cannot figure out where is that allocated.

Thank you so much for your help with this.
With all good wishes.