Hello

I have been checking whether the equal variances assumption can be applied to my data for twins 1 and 2 and for MZ and DZ twins. It appears that I can't make this assumption.

1. If I am not able to constrain the variances for twins 1 and 2, how do I deal with this? Do I just leave the model so that the variances for each twin are included separately in creating the covariance (ie as in the model below)? (script below)

2. If I come up with above (ie a significant difference in the -2LL for constraining of variances for twins 1 and 2 compared with the last model), how is it best to look at whether or not MZ/DZ variances can be constrained? Do we still just compare with the constraining of variances for Twin 1/Twin 2? Or go back to compare with the constrainining of the means (ie model twinSatSub2 below)?

3. If I needed to have the variances for MZ/DZ twins different how do I do this? Is it again like with the model used to set this up (see below) of keeping the MZ and DZ variances specified separately?

Many, many thanks!

Script

DASS
mzData
dzData

stMean stVar stCovMZ stCovDZ

nv selVars

twinSatSphereModel mxModel("MZ",

mxMatrix( type="Full", nrow=nv, ncol=2, free=TRUE,

values=stMean, label=c("MZm1", "MZm2"), name="interceptMZ" ),

mxMatrix( type="Full", nrow=1, ncol=1, free=TRUE,

values=-0.058, label="bsex", name="b_sex" ),

mxMatrix( type="Full", nrow=1, ncol=2, free=FALSE,

label="data.Sex", name="obs_sex" ),

mxMatrix( type="Full", nrow=1, ncol=1, free=TRUE,

values=-.124, label="bAge", name="b_age" ),

mxMatrix( type="Full", nrow=1, ncol=2, free=FALSE,

label=c("data.AgeA", "data.AgeB"), name="obs_age" ),

mxMatrix( type="Full", nrow=1, ncol=1, free=TRUE,

values=-.052, label="bedu", name="b_edu" ),

mxMatrix( type="Full", nrow=1, ncol=2, free=FALSE,

label=c("data.EducationA", "data.EducationB"), name="obs_edu" ),

mxMatrix( type="Full", nrow=1, ncol=1, free=TRUE,

values=1.041, label="bmarital", name="b_marital" ),

mxMatrix( type="Full", nrow=1, ncol=2, free=FALSE,

label=c("data.MaritalA", "data.MaritalB"), name="obs_marital" ),

mxMatrix( type="Full", nrow=1, ncol=1, free=TRUE,

values=10.091, label="bpsych", name="b_psych" ),

mxMatrix( type="Full", nrow=1, ncol=2, free=FALSE,

label=c("data.PsychexcludeA", "data.PsychexcludeB"), name="obs_psych" ),

mxMatrix( type="Full", nrow=1, ncol=1, free=TRUE,

values=14.538, label="bsphere", name="b_sphere" ),

mxMatrix( type="Full", nrow=1, ncol=2, free=FALSE,

label=c("data.SphereOnlyA", "data.SphereOnlyB"), name="obs_sphere" ),

mxAlgebra( expression= interceptMZ + b_sex %x% obs_sex + b_age %x% obs_age + b_edu %x% obs_edu + b_marital %x% obs_marital + b_psych %x% obs_psych + b_sphere %x% obs_sphere, name="expMeanMZ" ),

mxMatrix( type="Full", nrow=nv, ncol=nv, free=TRUE, values=stVar,

label=c("VMZ1"),name="VarMZtw1" ),

mxMatrix(type="Full", nrow=nv, ncol=nv, free=TRUE, values=stVar,

label=c("VMZ2"),name="VarMZtw2"),

mxMatrix(type="Full", nrow=nv, ncol=nv, free=TRUE, values=stCovMZ, name="CovMZ"),

mxAlgebra( expression= rbind ( cbind(VarMZtw1, t(CovMZ)),

cbind(CovMZ, VarMZtw2) ), name="expCovMZ" ),

mxData(observed=mzData, type="raw"),

mxFIMLObjective(covariance="expCovMZ", means="expMeanMZ", dimnames=selVars) ),

mxModel("DZ",

mxMatrix( type="Full", nrow=nv, ncol=2, free=TRUE,

values=stMean, label=c("DZm1", "DZm2"), name="interceptDZ" ),

mxMatrix( type="Full", nrow=1, ncol=1, free=TRUE,

values=-0.058, label="MZ.bsex", name="b_sex" ),

mxMatrix( type="Full", nrow=1, ncol=2, free=FALSE,

label="data.Sex", name="obs_sex" ),

mxMatrix( type="Full", nrow=1, ncol=1, free=TRUE,

values=-.124, label="MZ.bAge", name="b_age" ),

mxMatrix( type="Full", nrow=1, ncol=2, free=FALSE,

label=c("data.AgeA", "data.AgeB"), name="obs_age" ),

mxMatrix( type="Full", nrow=1, ncol=1, free=TRUE,

values=-.052, label="MZ.bedu", name="b_edu" ),

mxMatrix( type="Full", nrow=1, ncol=2, free=FALSE,

label=c("data.EducationA", "data.EducationB"), name="obs_edu" ),

mxMatrix( type="Full", nrow=1, ncol=1, free=TRUE,

values=1.041, label="MZ.bmarital", name="b_marital" ),

mxMatrix( type="Full", nrow=1, ncol=2, free=FALSE,

label=c("data.MaritalA", "data.MaritalB"), name="obs_marital" ),

mxMatrix( type="Full", nrow=1, ncol=1, free=TRUE,

values=10.091, label="MZ.bpsych", name="b_psych" ),

mxMatrix( type="Full", nrow=1, ncol=2, free=FALSE,

label=c("data.PsychexcludeA", "data.PsychexcludeB"), name="obs_psych" ),

mxMatrix( type="Full", nrow=1, ncol=1, free=TRUE,

values=14.538, label="MZ.bsphere", name="b_sphere" ),

mxMatrix( type="Full", nrow=1, ncol=2, free=FALSE,

label=c("data.SphereOnlyA", "data.SphereOnlyB"), name="obs_sphere" ),

mxAlgebra( expression= interceptDZ + b_sex %x% obs_sex + b_age %x% obs_age + b_edu %x% obs_edu + b_marital %x% obs_marital + b_psych %x% obs_psych + b_sphere %x% obs_sphere, name="expMeanDZ" ),

mxMatrix(type="Full", nrow=nv, ncol=nv, free=TRUE, values=stVar,

label=c("VDZ1"), name="VarDZtw1"),

mxMatrix(type="Full", nrow=nv, ncol=nv, free=TRUE, values=stVar,

label=c("VDZ2"), name="VarDZtw2"),

mxMatrix(type="Full", nrow=nv, ncol=nv, free=TRUE, values=stCovDZ, name="CovDZ"),

mxAlgebra( expression= rbind ( cbind(VarDZtw1, t(CovDZ)),

cbind(CovDZ, VarDZtw2) ), name="expCovDZ" ),

mxData(observed=dzData, type="raw"),

mxFIMLObjective(covariance="expCovDZ", means="expMeanDZ", dimnames=selVars)),

mxAlgebra(MZ.objective + DZ.objective, name="min2LL"),

mxAlgebraObjective("min2LL"))

twinSatSphereFit (twinSatFitSphereSum

# Specify and Run Sub-Model 1 equating means across twin order within zygosity group

twinSatModelSub1
twinSatModelSub1$MZ$interceptMZ
twinSatModelSub1$DZ$interceptDZ
twinSatFitSub1
twinSatModelSub1Summ
twinSatModelSub1Summ

tableFitStatistics(twinSatSphereFit, twinSatFitSub1)

# Specify and Run Sub-Model 2 equating means across twin order AND zygosity group

twinSatModelSub2
twinSatModelSub2$MZ$interceptMZ
twinSatModelSub2$DZ$interceptDZ
twinSatFitSub2
twinSatModelSub2Summ
twinSatModelSub2Summ

tableFitStatistics(twinSatFitSub1, twinSatFitSub2)

tableFitStatistics(twinSatSphereFit, twinSatFitSub2)

# Specify and Run Sub-Model 3 equating Variances across twin order in MZ and DZ group

twinSatModelSub3
twinSatModelSub3$MZ$VarMZtw1
twinSatModelSub3$MZ$VarMZtw2
twinSatModelSub3$DZ$VarDZtw1
twinSatModelSub3$DZ$VarDZtw2
twinSatFitSub3
twinSatModelSub3Summ
twinSatModelSub3Summ

tableFitStatistics(twinSatFitSub2,twinSatFitSub3)

# Specify and Run Sub-Model 4 equating Variances across twin order AND zygosity group

twinSatModelSub4
twinSatModelSub4$MZ$VarMZtw1
twinSatModelSub4$MZ$VarMZtw2
twinSatModelSub4$DZ$VarDZtw1
twinSatModelSub4$DZ$VarDZtw2
twinSatFitSub4
twinSatModelSub4Summ
twinSatModelSub4Summ

tableFitStatistics(twinSatFitSub3,twinSatFitSub4)

For next time, please include your script as an attachment. It will be easier to read the replies to the post.

It's pretty odd that you can't constrain variances across twins, as which twin is 1 and which is 2 is usually random. If the assignment is not random (i.e., opposite sex twins always put the male as 1 and the female as 2), then you've discovered a potentially interesting effect. If the assignment truly is random, then this is just an "unhappy randomization", and you should be able to re-assign twin numbers. One of the more experienced BG folk on this forum should probably weigh in on common practice, however.

Thanks for that thought. Our assignment to twin 1 and 2 was not random as it was based on their birth order...... any thoughts?

Then it appears you've found an effect of birth order. Is that relevant to your research?

In my experience, there are few effects of birth order in twins. What usually causes failure of the equal variances of T1 & T2 assumption is non-normality of the data. I recommend you (@koak) test this.

Thank you so, so much. This has fixed the problem totally!

Much indebted!