Comparing ACE and ADE models

Posted on
No user picture. inooradd Joined: 03/11/2020
Hello,

I have a couple of questions regarding comparing ACE and ADE models.

1. Is there a quick method of comparing the fit of both of these models? As of right now I am running the scripts for each model, then comparing the two models using mxcompare:

<%
mxCompare( fitADENOTW , fitACENOTW ) %>

I know typically when we compare nested, we run something to the effect of...

<%
#Run ACE model without Twin Effect
modelACENOTW <- mxModel( fitACE, name="oneACENOTWba" )
modelACENOTW <- omxSetParameters( modelACENOTW, labels="tw11", free=FALSE, values=0 )
fitACENOTW <- mxTryHard( modelACENOTW, intervals=F )
mxCompare( fitACE, fitACENOTW ) %>

But I am unsure how to incorporate the change from the ACE model:

<%
covDZ <- mxAlgebra( expression= 0.5%x%A+ C+Tw, name="cDZ" )
covSIB <- mxAlgebra( expression= 0.5%x%A+ C, name="cSIB" ) %>

to the ADE model:

<%
covDZ <- mxAlgebra( expression= 0.5%x%A+ 0.25%x%D+Tw, name="cDZ" )
covSIB <- mxAlgebra( expression= 0.5%x%A+ 0.25%x%D, name="cSIB" )%>

in a short, correct manner that would preclude me from running two separate scrips.

2. Even when I run these models and compare their fit, I am unable to test whether differences in their fit are statistically significant. Notice the P value indicates NA. I know I can reference the AIC descriptives, with the ADE fitting better than ACE in this instance, but I also want to include the p value, if possible.

<%
< mxCompare( fitADENOTW , fitACENOTW )
base comparison ep minus2LL df AIC diffLL diffdf p
1 oneADENOTWba 8 4879.663 1740 1399.663 NA NA NA
2 oneADENOTWba oneACENOTWba 8 4882.399 1740 1402.399 2.735366 0 NA%>

This may have something to do with the similar degrees of freedom, since chi-square fit statistics uses df. However, I have used chi-square tests of fitness before (SEM) and, if I am not mistaken, it yielded a p value.

Any help would be greatly appreciated.

Replied on Thu, 05/28/2020 - 21:07
Picture of user. tbates Joined: 07/31/2009

If you use umx, then for univariate ACE models, you can use `umxReduce(m1)`

This will run ADE/ACE, AE, CE, and give you a comparison table and the relative probability of each in this set.

If you have a multivariate model, you could `umxModify` the model, and get a comparison (shown below). But it's probably more straightforward to use the dzAr parameter in `umxACE` to build the ADE model and compare the two.

FYI:

data(twinData)
twinData[, c("ht1", "ht2")] = twinData[, c("ht1", "ht2")] * 10
mzData = twinData[twinData$zygosity %in% "MZFF", ]
dzData = twinData[twinData$zygosity %in% "DZFF", ]

m1 = umxACE(selDVs = "ht", selCovs = "age", sep = "", dzData = dzData, mzData = mzData)

# Make ADE model
m2 = umxModify(m1, "dzCr_r1c1", value = .25, name="ADE", comp=TRUE)

ACE -2 × log(Likelihood) = 5944.845
Standardized solution

| | a1| d1| e1|
|:--|-----:|--:|----:|
|ht | 0.933| 0| 0.36|

Means: Intercept and (raw) betas from model$top$intercept and model$top$meansBetas

| | ht1| ht2|
|:---------|------:|------:|
|intercept | 16.446| 16.446|
|age | -0.005| -0.005|

|Model | EP|Δ -2LL |Δ df |p | AIC| Δ AIC|Compare with Model |
|:-----|--:|:---------|:----|:--|---------:|---------:|:------------------|
|ACE | 5| | | | -1843.155| 0.000000| |
|ADE | 5|-0.014679 |0 | | -1843.169| -0.014679|ACE |

Replied on Mon, 06/01/2020 - 18:04
No user picture. inooradd Joined: 03/11/2020

In reply to by tbates

I tried installing and using this function, but unfortunately its not available for R version 4.0.0.

Warning in install.packages :
package ‘umxReduce’ is not available (for R version 4.0.0)

While I wait for it to become available I will continue running two separate models then comparing them. Thank you!!

Replied on Mon, 06/01/2020 - 18:56
No user picture. inooradd Joined: 03/11/2020

In reply to by tbates

I was able to get it to run, but it seems to omit my sibling data set.

<%
m1 = umxACE(selDVs = "fem", selCovs=c("age","sex"), sep = "", mzData = mzData, dzData = dzData, sibData = sibData)
%>

and I get this message

<%
> m1 = umxACE(selDVs = "fem", selCovs=c("age","sex"), sep = "", mzData = mzData, dzData = dzData, sibData = sibData)
Error in umxACE(selDVs = "fem", selCovs = c("age", "sex"), sep = "", mzData = mzData, :
unused argument (sibData = sibData) %>

and it runs the model fitting with only two out of the three datasets.

<%
|Model | EP|Δ -2LL |Δ df |p | AIC| Δ AIC|Compare with Model |
|:-----|--:|:---------|:----|:-----|--------:|---------:|:------------------|
|ACE | 6| | | | 1059.607| 0.000000| |
|ADE | 6|1.294747 |0 | | 1060.901| 1.294747|ACE |
|CE | 5|0.3068545 |1 |0.580 | 1057.913| -1.693145|ACE |
|AE | 5|1.294747 |1 |0.255 | 1058.901| -0.705253|ACE | %>

Any idea how to fix this?

Replied on Tue, 06/02/2020 - 12:01
Picture of user. tbates Joined: 07/31/2009

Hi,
I've not implemented 5-group or extended family support in umxACE. It's not high on the priority list currently. I consider it from time to time but it would be a major development project, as IP, CP, etc all would have to grow in capability to make it worth doing.
Replied on Tue, 06/02/2020 - 13:16
No user picture. inooradd Joined: 03/11/2020

I am just excited to have used it. Ill stick with the multiple scripts for now. Thank you so much for this helpful piece of code I will hopefully use in my future studies.
Replied on Fri, 06/05/2020 - 13:45
Picture of user. AdminNeale Joined: 03/01/2013

In reply to by inooradd

In some ways, the difference between these models (which can be judged by AIC but not likelihood-ratio test) is moot when we consider unsigned variance components. Since C and D are confounded, both may be operating, even if rMZ=2rDZ, which would yield an AE model. The ADE and ACE will fit usually differently in a path coefficient specification because it won't allow D<0 or C<0. However, if c^2 or d^2 are allowed to go negative, the two models fit equally well. I find that this helps with the interpretation of these coefficients. When rMZ>2rDZ, a positive estimate of d^2 would be returned, but a negative estimate of c^2. Here the interpretability of the estimate counts more than the quantitative fit index or likelihood (which is the same for ACE and ADE). Note also that sampling error might make true world values of 40:10:50 A:C:E be estimated as 55:-5:50. Inspection of the confidence interval, or a likelihood ratio test of C=0 might reveal little loss of fit from fixing C to zero or to 10. Regardless, we should acknowledge that C and D are confounded and with twin data alone they are effectively estimated as a joint parameter combining non-addiitivity with shared environment (or the genetic effects of assortative mating which is also confounded with C in a classical twin study design).