Twin Scripts
Message/Author
 George E. Higgins posted on Monday, June 05, 2006 - 7:54 am
Using twin data, I'm interested in producing the following model:

y1 --> y2 --> y3

I also need to take into account the heritability (hsquared) and the environmental influences (csquared). Are there any Mplus scripts available for such a model?
 Bengt O. Muthen posted on Monday, June 05, 2006 - 9:08 am
We don't have a script of this kind, but the analysis is possible. I take it that you mean that for each twin you have 3 outcomes, with an auto-regressive model imposed. So that the twin model would have 6 dependent variables and an ACE covariance structure. One question is if the ACE model should be imposed for y1 only, or 3 separate ACE models imposed for the 3 y's. The analogy would be a growth model for 3 time points where perhaps only the random intercept has an ACE model, but where the random slope could have an ACE model as well.
 George E. Higgins posted on Monday, June 05, 2006 - 10:12 am
I do mean for the 3 y's to have an ACE model. In the example, that you presented would there be example in the Mplus users guide for growth models that I may modify?
 Bengt O. Muthen posted on Monday, June 05, 2006 - 10:37 am
Not exactly, but I can guide you. You can start from example 7.29 and view the 4 indicators instead as one continuous outcome observed 4 times (3 times in your case). Here, the ACE model has the factor for the 2 twins as the dependent variables. In your case, the factor is the random intercept which implies that you should have unit loadings. And hold the intercept equal across time - and of course across twins.
 Joonmo Son posted on Thursday, May 03, 2007 - 8:33 am
Hi, I am working on twins study using ACE model script (Example 5.18 in Mplus 4.0 manual). But the model fit is very bad like below:

TESTS OF MODEL FIT

Chi-Square Test of Model Fit

Value 290.502
Degrees of Freedom 6
P-Value 0.0000

Chi-Square Contributions From Each Group

MZ 102.736
DZ 187.766

Chi-Square Test of Model Fit for the Baseline Model

Value 45.554
Degrees of Freedom 2
P-Value 0.0000

CFI/TLI

CFI 0.000
TLI -1.177

And the output presented error message:

WARNING: THE LATENT VARIABLE COVARIANCE MATRIX (PSI) IN GROUP MZ
IS NOT POSITIVE DEFINITE. THIS COULD INDICATE A NEGATIVE VARIANCE/
RESIDUAL VARIANCE FOR A LATENT VARIABLE, A CORRELATION GREATER OR EQUAL
TO ONE BETWEEN TWO LATENT VARIABLES, OR A LINEAR DEPENDENCY AMONG MORE
PROBLEM INVOLVING VARIABLE A2.

Would you tell me how I can fix these problems? Thank you in advance!
 Linda K. Muthen posted on Thursday, May 03, 2007 - 10:13 am
The warning is fine because you have fixed the covariance between a1 and a2 to one. You should have a similar warning to c2. You can ask for modification indices to see where the model misfit is.
 Joonmo Son posted on Thursday, May 03, 2007 - 12:08 pm
You are amazing! Yes, there was another warning regarding c2. But you mean that I can ignore the warnings, right?

Minimum M.I. value for printing the modification index 10.000

M.I. E.P.C. Std E.P.C. StdYX E.P.C.
Group MZ

BY Statements

A1 BY DON1 28.701 27.729 27.729 0.215
A2 BY DON1 28.701 54.757 54.757 0.424
C1 BY DON1 28.701 465.660 465.660 3.609
C2 BY DON1 999.000 0.000 0.000 0.000
E1 BY DON1 28.698 24.433 24.433 0.189

Would you let me know what I should do with the information above?

Model fit became far better after I excluded outliers from the data. Do you have any other recommendations that will help improve the model fit?

Thank you so much!
 Joonmo Son posted on Thursday, May 03, 2007 - 2:39 pm
I have another question.
Harold Snieder's 2000 piece on ACE model reported model fitting results, Va, Vc, Ve, and Vp.

Vc: Shared environmental variance
Ve: unique environmental variance
Vp: Total phenotypic variance

And

Vp=Va+Vc+Ve

Can I get such results from Mplus output?
If so, would you let me know how I can get them? Thanks a lot.
 Linda K. Muthen posted on Friday, May 04, 2007 - 8:18 am
Please send your input, data, output, and license number to support@statmodel.com so I can see the full model.

You can create vp using MODEL CONSTRAINT similar to the way heritability is created in Example 5.22.
 Joonmo Son posted on Wednesday, June 20, 2007 - 9:10 am
With your great help, I was able to run the ACE model. Now I am trying to plug in some control variables such as parents' education, R's education, and R's income in the ACE model. The outcome variable is a factor of volunteering and donating behaviors.

Could you let me know how I should construct the program based on Example 5.18 or 5.21? Thank you in advance!
 Linda K. Muthen posted on Wednesday, June 20, 2007 - 9:29 am
It seems that covariates are related to the ACE components. See the Harden or Prescott papers on our website under Genetic Modeling to see how this is done.
 MAH posted on Friday, July 18, 2008 - 1:07 pm
Hi,

I am trying to learn to do twin analyses, so thank you for posting the genetic scipts. I have used your 4-group script "rawvc2a" translated from MX. I have a question concerning your use of sex as a covariate. In the first model statement, sex is covaried out of bmi. Why don't you get an error message that 'sex1' and 'sex2' has no variance for the MZM, DZM, MZF, and DZF group?

Thank you.
 MAH posted on Friday, July 18, 2008 - 3:40 pm
One more thing. In the output, the A,C,and E parameter estimates from the examples are, I think, unstandardized. How do I obtain standardized components or components that are interpretable as the % variance in the outcome that can be attributed to A, C, and E.

Thank you!
 Bengt O. Muthen posted on Friday, July 18, 2008 - 4:53 pm
Re: post of 1:07

This is an artificial example. In a real data example you should get such a message.
 Bengt O. Muthen posted on Friday, July 18, 2008 - 4:55 pm
Re: post of 3:40

See how "h" is used in Model Constraint in UG Ex 5.21.
 Joonmo Son posted on Tuesday, December 16, 2008 - 12:22 am
Hi, I got a question on how to make a new variable of correlation between additive genetic factors of males and females in the model constraint shown below:

MODEL CONSTRAINT:
NEW(af cf ef vf af2 cf2 ef2 staf stcf stef
am cm em vm am2 cm2 em2 stam stcm stem);

varf = af**2 + cf**2 + ef**2;
covmzf = af**2 + cf**2;
covdzf = 0.5*af**2 + cf**2;
vf = varf;
af2 = af**2;
cf2 = cf**2;
ef2 = ef**2;
staf = af2/vf;
stcf = cf2/vf;
stef = ef2/vf;

varm = am**2 + cm**2 + em**2;
covmzm = am**2 + cm**2;
covdzm = 0.5*am**2 + cm**2;
vm = varm;
am2 = am**2;
cm2 = cm**2;
em2 = em**2;
stam = am2/vm;
stcm = cm2/vm;
stem = em2/vm;

Would you let me know how I can generate a new variable of correlation coefficient between "af" (genetic factor for males) and "am" (genetic factor for females)? Eventually, I will set the correlation to unity to identify whether the genetic impact for both genders is the same. Thank you in advance.
 Linda K. Muthen posted on Tuesday, December 16, 2008 - 11:18 am
I would test that the genetic impact for the genders is the same by doing a chi-square or loglikelihood difference test. I would test the model that you have specified versus a model where am is changed to af in the following part of MODEL CONSTRAINT while the rest remains as is:

varm = am**2 + cm**2 + em**2;
covmzm = am**2 + cm**2;
covdzm = 0.5*am**2 + cm**2;
vm = varm;
am2 = am**2;
cm2 = cm**2;
em2 = em**2;
stam = am2/vm;
stcm = cm2/vm;
stem = em2/vm;
 Joonmo Son posted on Tuesday, December 16, 2008 - 7:35 pm
Thanks. But I had already done that as below. And there was no significant difference from the unconstrained model according to chi2 test:

varm = am**2 + cm**2 + em**2;
covmzm = am**2 + cm**2;
covdzm = 0.5*am**2 + cm**2;
vm = varm;
am2 = am**2;
cm2 = cm**2;
em2 = em**2;
stam = am2/vm;
stcm = cm2/vm;
stem = em2/vm;

am=af;

What I am now trying to figure out is a way to get the correlation coefficient between af and am (e.g., Rfm=.80) even though I already noticed that they are pretty much similar to each other. And I may set it to be 1 (i.e., Rfm=1) to test the equality of them in a different manner even though I am sure that I will get the same outcome. But I found it difficult to build such command for correlation coefficient in the model constraint. Would you help me?
 Bengt O. Muthen posted on Wednesday, December 17, 2008 - 10:02 am
If you want correlations, I think you have to choose the alternative way of doing ACE modeling, namely define ACE factors. See the Mplus UG.
 Erika Wolf posted on Friday, April 03, 2009 - 7:55 am
I'm running a twin model with two latent psychopathology factors and separate latent A C E factors for the two latent psychopathology factors. I wanted to run a model in which both latent psychopathology factors cross-load on both sets of latent ACE factors and the plan was then to drop paths that I don't think are needed from a theoretical standpoint and use a chi square difference testing approach. However, when I ran the parent model described above, I get the following message:

THE STANDARD ERRORS OF THE MODEL PARAMETER ESTIMATES COULD NOT BE COMPUTED. THE MODEL MAY NOT BE IDENTIFIED. CHECK YOUR MODEL. PROBLEM INVOLVING PARAMETER 26.

THE CONDITION NUMBER IS -0.469D-16.

The parameter it is referring to is one that I know from prior work with this dataset is approximately 0, so I'm wondering if this could be the explanation for the error message (b/c the estimate is so close to the boundary). If I set that path to 0, I get the same type of message about another path that, again, from prior work, I know approximates 0. If I then set that second parameter to 0 as well, my model runs fine. Can you confirm my interpretation of this? Thanks.
 Bengt O. Muthen posted on Friday, April 03, 2009 - 11:03 am
Is it the variance of c?
 Erika Wolf posted on Friday, April 03, 2009 - 11:21 am
In the first instance of this, it is actually the cross-path from the second genetic factor (A) to the first latent psychopathology factor (I know from prior analyses that the cross-path actually works the other way, the first genetic factor, A, exerts effects on both psychopathology factors, which is also theoretically appropriate).

In the second instance of this type of error message, it is the cross path from the second non-shared environment factor (E) to the first latent psychopathology factor (and I'd imagine there is not much variance left over in the first psychopathology factor, given what is already explained by the first A and E factors associated with it).

The variances of all the ACE factors are fixed at 1.

Any advice is greatly appreciated! Thanks!
 Bengt O. Muthen posted on Saturday, April 04, 2009 - 4:15 pm
Please send the input, output and data for this with your license number to support@statmodel.com so we can take a closer look at it.
 Joonmo Son posted on Friday, May 08, 2009 - 7:28 am
Hi, I need your advice on estimator of twin modeling and chi-square test.

I used ACE modeling to see whether there is genetic effect in a prosocial behavior (i.e., volunteering). The default ML estimator reported that there is gender difference in genetic component because genetic effect on volunteering was significantly greater for females than for males based on the chi-square model fit comparison test. But the endogenous variable (volunteering hours per month) is nonnormal because it is highly skewed to the left. So I tested the same model with GLS and MLM (to use Satorra-Bentler scaled chi-square). GLS estimator reported very similar outcome I had with ML. But Satorra-Bentler scaled chi-square test informed me that there is no gender difference between male and female twins in terms of genetic impact. I would like to ask you which estimator you would recommend for chi-square test in this kind of case. By the way, the sample size is 501 with four subgroups by zygosity. Thank you in advance.
 Linda K. Muthen posted on Saturday, May 09, 2009 - 10:53 am
The differences are likely due to the fact that your variables are skewed. I would use MLM to protect against non-normality.
 Joonmo Son posted on Sunday, May 10, 2009 - 4:36 am
 Erika Wolf posted on Thursday, September 17, 2009 - 9:37 am
Hi,
I ran an ACE model on two common factors (2 sets of distinct ACE components for each factor) and now need to calculate residual ACE components for each of the indicators of my 2 factors. Can you give me a sense of how I would go about doing this? Would I create new variables and regress the indicators on the new variables, constraining the paths, given the MZ or DZ status of the twin pair? Thanks!
 Bengt O. Muthen posted on Thursday, September 17, 2009 - 5:19 pm
Do you mean doing ACE modeling of the residuals for the indicators? If so, you can replace the residuals by factors which in turn will be decomposed into ACE.
 Erika Wolf posted on Friday, September 18, 2009 - 8:18 am
Yes, that is what I mean--thanks! I'm working on the syntax now!
 Erika Wolf posted on Friday, September 18, 2009 - 9:47 am
This is a follow-up to my post from yesterday:
I'm running into a problem with one non-shared environmental component of an indicator's residual, so I'd like to make sure I did this generally correctly.

In order to decompose the residuals of the indicators into the ACE components, I added in one set of ACE factors for each indicator using BY statements (e.g., A5 BY var1, C5 BY var1, E5 by var1) and then set the correlations between the ACE factors as you normally would (e.g., r = 1 for the As for MZ, r = .5 for the A factors for the DZ, etc) and set cross correlations between ACE factors for different indicators to 0. I also set the ACE factor variances to 1 and the ACE factor means to 0, as I had previously done for the ACE factors for my latent variables in the model. Is this what you meant by replacing the residuals with factors? Thanks for your help!
 Bengt O. Muthen posted on Friday, September 18, 2009 - 10:10 am
First redo your original model using factors for the residuals. Check that you get the same number of parameters and logL.

Then take the next step.
 Erika Wolf posted on Tuesday, September 22, 2009 - 6:35 am
This is an updated question to my last post (above) of 9/18: When I replaced the residuals with factors, I ran into many problems, including inability to compute standard errors and chi square. Perhaps this was because some of the residuals are quite close to the boundary of 0? At any rate, I then re-wrote my syntax using the model constraint approach (rather than A1 BY factor 1, etc) and I was able to run ACE modeling of the residuals of the indicators without a problem. My question is this: for subsequent analyses, it will be easier if I can return to the alternate (non model constraint) approach for modeling the main ACE factors. I'm wondering if it is acceptable to run that portion of the model using the model comand and model the ACE components of the residuals of the indicators using the model constraint approach in the same script? I tested this out (the hybrid script approach vs just the model constraint approach) and got equivalent results (loadings), but of course the # of freely estimated parameters differs somewhat. Thanks for your help!
 Erika Wolf posted on Tuesday, September 22, 2009 - 6:42 am
One additional question, I noticed when I added in ACE modeling of the residuals of the indicators, one of my main genetic paths to the latent phenotype decreased somewhat compared to the model without the ACE modeling of the residuals of the indicators; in addition two of the indicator loadings on the latent phenotype changed (one increased, one decreased). I was surprised by this as I thought that the ACE modeling of the residual factors should just divy up the residual variance of the indicators, after accounting for the latent factor, and should not affect the indicator loadings or other genetic paths in the model. Any thoughts on this?
 Bengt O. Muthen posted on Tuesday, September 22, 2009 - 10:54 am
Did you do the first step that I mentioned? It should not have any problems but replicate your LL (and #par.'s). If you have problems with it, then you are not doing it right.
 Erika Wolf posted on Tuesday, September 22, 2009 - 11:12 am
I do get the same # of free parameters and chi square value (I'm using WLSMV estimator as I have three-level categorical data), but I get the following message (sorry I quoted the wrong message above) for both my MZ and DZ groups for the same variable (I don't get this message in the prior step without the factors for the residuals of the indicators):

WARNING: THE RESIDUAL COVARIANCE MATRIX (THETA) IN GROUP MZ IS NOT
POSITIVE DEFINITE. THIS COULD INDICATE A NEGATIVE VARIANCE/RESIDUAL
VARIANCE FOR AN OBSERVED VARIABLE, A CORRELATION GREATER OR EQUAL TO ONE
BETWEEN TWO OBSERVED VARIABLES, OR A LINEAR DEPENDENCY AMONG MORE THAN TWO
PROBLEM INVOLVING VARIABLE ALCH.

For what it's worth, I set the factor means to 0 and the factor variances to 1 for the factors replacing the residuals of the indicators. Thanks again for all your help. This website is invaluable!
 Erika Wolf posted on Tuesday, September 22, 2009 - 12:55 pm
Also of note, although the loadings on the factors (that are replacing the residual variances) should be 1 as they are all single indicators and I set them to be equal across twin pairs, I noticed that in the standardized output, I get values other than 1 for the DZ group only (often taking out of range values greater than 1 on many variables).
 Bengt O. Muthen posted on Tuesday, September 22, 2009 - 1:43 pm
If you have categorical outcomes, you don't have the possibility of estimating parameters related to the residuals. If you look at a regular factor analysis of binary outcomes for example, you see that there is no free parameter related to the residual variance. The residual variance is a "remainder", making the total y* variance add to one. How this affects ACE modeling is explained further in the Prescott (2004) article on our web site under Genetics.
 Amber Singh posted on Wednesday, September 30, 2009 - 12:40 pm
I'm trying run a sex limitation model (rawVC3a) and I receive the following error message:

*** ERROR
One or more variables have a variance of zero.
Check your data and format statement.

From an earlier post, I understand that I should get this message (the variance of sex for MZM, for example, is zero). But the model then terminates...how can I get around this?
 Linda K. Muthen posted on Thursday, October 01, 2009 - 9:51 am
If you are using sex as a grouping variable, you cannot also use it as a covariate in the model which is what it sounds like you are doing. You can add VARIANCES=NOCHECK; to the DATA command and see if that helps.
 Amber Singh posted on Thursday, October 01, 2009 - 1:59 pm
Thanks. Sex is not the grouping variable, zygosity is the grouping variable. The VARIANCES=NOCHECK seemed to do the trick.
 Lulu Dong posted on Friday, April 01, 2011 - 3:57 pm
I am following the common pathway twin model script on the mplus website. In my output, Mplus did not give the Standardized Coefficients (under STDYX etc.) for the common and As, Cs, Es, as well as Xs, Ys, and Zs. I wonder why this is and how I can get standardized output for these parameters.

Thanks!
 Linda K. Muthen posted on Saturday, April 02, 2011 - 8:43 am
 Erika Wolf posted on Thursday, September 29, 2011 - 12:38 pm
If I run a bivariate Cholesky model using the latent variable approach (as in Carol Prescott's example 5), shouldn't I get relatively equivalent estimates of heritability (and the ACE paths) as I would if I run the same Cholesky model using the model constraint approach, as in the mplus example rawvc4a?
 Linda K. Muthen posted on Friday, September 30, 2011 - 8:32 am