# Intra Class Corelations

# Select Variables for Analysis

vars <- "AG_Ln" # list of variables names

nv <- 1 # number of variables

ntv <- nv*2 # number of total variables

selVars <- paste(vars,c(rep(1,nv),rep(2,nv)),sep="")
# Select Covariates for Analysis

#twinData[,'age'] <- twinData[,'age']/100

#twinData <- twinData[-which(is.na(twinData$age)),]

covVars <- c("age","Sex1","Sex2")

# Select Data for Analysis

mzData <- subset(twinData, Zyg==1, c(selVars, covVars))

dzData <- subset(twinData, Zyg==2, c(selVars, covVars))

# Set Starting Values

svMe <- 4 # start value for means

svPa <- .4 # start value for path coefficient

svPe <- .8 # start value for path coefficient for e

lbPa <- .01 # start value for lower bounds

# ------------------------------------------------------------------------------

# PREPARE MODEL

# ACE Model

# Create Matrices for Covariates and linear Regression Coefficients

defAge <- mxMatrix( type="Full", nrow=1, ncol=1, free=FALSE, labels=c("data.age"), name="Age" )

pathB1 <- mxMatrix( type="Full", nrow=1, ncol=1, free=TRUE, values=.01, label="b11", name="b1" )

defSex <- mxMatrix( type="Full", nrow=1, ncol=2, free=FALSE, labels=c("data.Sex1", "data.Sex2"), name="Sex" )

pathB2 <- mxMatrix( type="Full", nrow=1, ncol=1, free=TRUE, values=.01, label="b12", name="b2" )

# Create Algebra for expected Mean Matrices

meanG <- mxMatrix( type="Full", nrow=1, ncol=ntv, free=TRUE, values=svMe, labels="xbmi", name="meanG" )

expMean <- mxAlgebra( expression= meanG + cbind(b1%*%Age,b1%*%Age)+b2%*%Sex, name="expMean" )

# Create Matrices for Path Coefficients

pathA <- mxMatrix( type="Lower", nrow=nv, ncol=nv, free=TRUE, values=svPa, label="a11", lbound=lbPa, name="a" )

pathC <- mxMatrix( type="Lower", nrow=nv, ncol=nv, free=TRUE, values=svPa, label="c11", lbound=lbPa, name="c" )

pathE <- mxMatrix( type="Lower", nrow=nv, ncol=nv, free=TRUE, values=svPe, label="e11", lbound=lbPa, name="e" )

# Create Algebra for Variance Components

covA <- mxAlgebra( expression=a %*% t(a), name="A" )

covC <- mxAlgebra( expression=c %*% t(c), name="C" )

covE <- mxAlgebra( expression=e %*% t(e), name="E" )

# Create Algebra for expected Variance/Covariance Matrices in MZ & DZ twins

covP <- mxAlgebra( expression= A+C+E, name="V" )

covMZ <- mxAlgebra( expression= A+C, name="cMZ" )

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

expCovMZ <- mxAlgebra( expression= rbind( cbind(V, cMZ), cbind(t(cMZ), V)), name="expCovMZ" )

expCovDZ <- mxAlgebra( expression= rbind( cbind(V, cDZ), cbind(t(cDZ), V)), name="expCovDZ" )

# Create Data Objects for Multiple Groups

dataMZ <- mxData( observed=mzData, type="raw" )

dataDZ <- mxData( observed=dzData, type="raw" )

# Create Expectation Objects for Multiple Groups

expMZ <- mxExpectationNormal( covariance="expCovMZ", means="expMean", dimnames=selVars )

expDZ <- mxExpectationNormal( covariance="expCovDZ", means="expMean", dimnames=selVars )

funML <- mxFitFunctionML()

# Create Model Objects for Multiple Groups

pars <- list(pathB1,pathB2, meanG, pathA, pathC, pathE, covA, covC, covE, covP)

defs <- list(defAge, defSex)

modelMZ <- mxModel( name="MZ", pars, defs, expMean, covMZ, expCovMZ, dataMZ, expMZ, funML )

modelDZ <- mxModel( name="DZ", pars, defs, expMean, covDZ, expCovDZ, dataDZ, expDZ, funML )

multi <- mxFitFunctionMultigroup( c("MZ","DZ") )

# Create Algebra for Variance Components

rowVC <- rep('VC',nv)

colVC <- rep(c('A','C','E','SA','SC','SE'),each=nv)

estVC <- mxAlgebra( expression=cbind(A,C,E,A/V,C/V,E/V), name="VC", dimnames=list(rowVC,colVC))

# Create Confidence Interval Objects

ciACE <- mxCI( "VC" )

# Build Model with Confidence Intervals

modelACE <- mxModel( "oneACEca", pars, modelMZ, modelDZ, multi, estVC, ciACE )

# ------------------------------------------------------------------------------

# RUN MODEL

# Run ACE Model

fitACE <- mxTryHard( modelACE, intervals=T, extraTries = 10, finetuneGradient = T, exhaustive = T )

sumACE <- summary( fitACE )

summary(fitACE)

# Compare with Saturated Model

mxCompare( fit, fitACE )

#lrtSAT(fitACE,-2ll,df)

# Print Goodness-of-fit Statistics & Parameter Estimates

fitGofs(fitACE)

fitEsts(fitACE)

# ------------------------------------------------------------------------------

# RUN SUBMODELS

# Test Significance of Covariate

modelCov <- mxModel( fitACE, name="testCov" )

modelCov <- omxSetParameters( modelCov, label=c("b11","b12"), free=FALSE, values=0 )

fitCov <- mxRun( modelCov )

mxCompare( fitACE, fitCov )

# Run AE model

modelAE <- mxModel( fitACE, name="oneAEca" )

modelAE <- omxSetParameters( modelAE, labels="c11", free=FALSE, values=0 )

fitAE <- mxRun( modelAE, intervals=T )

mxCompare( fitACE, fitAE )

fitGofs(fitAE)

fitEsts(fitAE)

summary(fitAE)

# Run CE model

modelCE <- mxModel( fitACE, name="oneCEca" )

modelCE <- omxSetParameters( modelCE, labels="a11", free=FALSE, values=0 )

fitCE <- mxRun( modelCE, intervals=T )

mxCompare( fitACE, fitCE )

fitGofs(fitCE)

fitEsts(fitCE)

# Run E model

modelE <- mxModel( fitAE, name="oneEca" )

modelE <- omxSetParameters( modelE, labels="a11", free=FALSE, values=0 )

fitE <- mxRun( modelE, intervals=T )

mxCompare( fitAE, fitE )

fitGofs(fitE)

fitEsts(fitE)

`# Print Comparative Fit Statistics`

mxCompare( fitACE, nested <- list(fitCov, fitAE, fitCE, fitE) )

round(rbind(fitACE$VC$result,fitAE$VC$result,fitCE$VC$result,fitE$VC$result),4)

Thank you so much.

## cov2cor()?

`cov2cor()`

, and what results does it give you?Log in or register to post comments

In reply to cov2cor()? by AdminRobK

## Hi,

I have tried with this:

mxEval(cov2cor(V), fitACE, T) the result is: [,1] [1,] 1

and with: mxEval(cov2cor(expCovDZ), fitACE, T) the result is:

Error: The following error occurred while evaluating the expression 'cov2cor(expCovDZ)' in model 'oneACEca' : 'V' is not a square numeric matrix

Thank you so much!

Log in or register to post comments

## Ahhh

shouldproduce a 1x1 matrix with a 1 in it.The second one is a namespace issue. `mxEval` is confused about the R object called `expCovDZ` which is an `MxAlgebra` and the OpenMx object given the same name which can be interpreted as a numeric matrix. Name the R object something else and you should get a 2x2 correlation matrix.

# change this line

expCovDZ <- mxAlgebra( expression= rbind( cbind(V, cDZ), cbind(t(cDZ), V)), name="expCovDZ" )

# to something like

rexpCovDZ <- mxAlgebra( expression= rbind( cbind(V, cDZ), cbind(t(cDZ), V)), name="expCovDZ" )

`# then put 'rexpCovDZ' into the model.`

Log in or register to post comments

In reply to Ahhh by mhunter

## Thank you so much for your

Log in or register to post comments

In reply to Thank you so much for your by JuanJMV

## Are you having a new problem?

Log in or register to post comments

In reply to Thank you so much for your by JuanJMV

## my guess

I'm going to take a wild guess that all three of the following lines of code will do what you want for the DZ twins (and with the obvious changes, what you want for the MZ twins):

mxEval(cov2cor(expCovDZ), fitACE$DZ, T)

cov2cor(mxGetExpected(fitACE$DZ, "covariance"))

mxEval(cov2cor(DZ.expCovDZ), fitACE, T)

Am I correct?

Log in or register to post comments

In reply to my guess by AdminRobK

## Yeah!!! you are right. Thank

Thanks for your help I really appreciate it.

Log in or register to post comments

In reply to my guess by AdminRobK

## Just one thing more I did the

[,1] [,2] [,3] [,4] [,5] [,6]

[1,] 1.00000000 0.66182596 0.37644256 0.62216985 0.54028297 0.42317649

[2,] 0.66182596 1.00000000 0.32231799 0.54028297 0.70316889 0.32361090

[3,] 0.37644256 0.32231799 1.00000000 0.42317649 0.32361090 0.78676051

[4,] 0.62216985 0.54028297 0.42317649 1.00000000 0.66182596 0.37644256

[5,] 0.54028297 0.70316889 0.32361090 0.66182596 1.00000000 0.32231799

[6,] 0.42317649 0.32361090 0.78676051 0.37644256 0.32231799 1.00000000

So I should have 3 intra-class correlations but I don't know what value is for each variable.

Also, I would like to know if the intra-class correlations should be the same in the univariate analyses and trivariate analysis?

Thanks in advance

Log in or register to post comments

In reply to Just one thing more I did the by JuanJMV

## between-twin, within-trait?

Log in or register to post comments

In reply to between-twin, within-trait? by AdminRobK

## Thank you so much!!! and the

Log in or register to post comments

In reply to Thank you so much!!! and the by JuanJMV

## Right.

Log in or register to post comments

In reply to Right. by AdminRobK

## Thank you so much!!!

Log in or register to post comments