Login

News

Decomposed ACE variances...now trying to run regressions

Posted on
No user picture. sppspp10 Joined: 02/05/2021
Forums
Behavioral Genetics Models
I am 100% new to twin data and OpenMx, and so far the tutorials and documentations have been **superbly** helpful!
I used this [genetic epi path specification annotated guide] (https://vipbg.vcu.edu/vipbg/OpenMx2/docs//OpenMx/latest/GeneticEpi_Path.html#ace-model-a-twin-analysis) to decompose the variances of ACE on a given phenotype (my DV).

Now, though, I need to run a regression and I can't figure out how to feed the ACE variances from above into the regressions. I assume I have to constrain the parameters somehow?

I have a categorical IV (an unshared environmental factor) ~ continuous Phenotype DV + covariates

And my research question is whether the E (unshared environmental factor) is significantly associated with the DV, after controlling for A and C confounders.

Is there any documentation I missed on this script?

Thank you so much for any help.

Replied on Sun, 03/07/2021 - 16:55
Picture of user. AdminRobK Joined: 01/24/2014

You should post your script as it currently stands.

Now, though, I need to run a regression and I can't figure out how to feed the ACE variances from above into the regressions. I assume I have to constrain the parameters somehow?

I have a categorical IV (an unshared environmental factor) ~ continuous Phenotype DV + covariates

And my research question is whether the E (unshared environmental factor) is significantly associated with the DV, after controlling for A and C confounders.

I don't follow this part. In particular, I think you mean "continuous Phenotype DV ~ categorical IV (an unshared environmental factor) + covariates"? Also, what makes you so sure the IV is nonshared-environmental? How many categories does it have, and are they ordered?

Log in or register to post comments
Replied on Wed, 03/10/2021 - 13:31
No user picture. sppspp10 Joined: 02/05/2021

THanks, AdminRobK.
Script is below, but so far it just includes the ACE variance breakdown. I do not understand how to write the regression script (which is where I am stuck). To clarify your
I don't follow this part. In particular, I think you mean "continuous Phenotype DV ~ categorical IV (an unshared environmental factor) + covariates"? Also, what makes you so sure the IV is nonshared-environmental? How many categories does it have, and are they ordered?

:
1. You are right, I meant that the continuous DV ~ categorical IV + covariates
2. The IV being non-shared is the crux of my research question. My research question asks if going outdoors into nature affects mental health. So, the behavior of going outside into nature (which may have some genetic influences, but is also the observed difference in unshared environmental, i.e. the behavior) is hypothesized to be associated with improved mental health (the phenotype DV).
3. The IV has 5 categories, ordered.

Any help feeding my ACE model summary into the regression model is needed.

# Select Variables for Analysis
selVars <- c('cohenpss_score_A','cohenpss_score_B') #Selecting these columns for analysis
aceVars <- c("A1","C1","E1","A2","C2","E2") #creating latent variables

# Select Data for Analysis: Creating two subsets for MZ and DZ to generate descriptive stats for each group type
mzData <- subset(df, zyg_A==1, selVars)
dzData <- subset(df, zyg_A==0, selVars)

# Generate Descriptive Statistics
colMeans(mzData,na.rm=TRUE) #mean ~11
colMeans(dzData,na.rm=TRUE)
cov(mzData,use="complete") #cov ~45
cov(dzData,use="complete")

# Path objects for Multiple Groups
manifestVars=selVars #Manifest = the Stress phenotype
latentVars=aceVars

# variances of latent variables
latVariances <- mxPath( from=aceVars, arrows=2,
free=FALSE, values=1 )

# means of latent variables... AND SPECIFYING SINGLE=HEADED ARROWS FROM TRIANGLE (WHICH IS FIXED AT 1) TO EACH OF THE LATENT VARIABLES, FIXED AT 0

latMeans <- mxPath( from="one", to=aceVars, arrows=1,
free=FALSE, values=0 )

# means of observed variables

obsMeans <- mxPath( from="one", to=selVars, arrows=1,
free=TRUE, values=10, labels="mean" ) #The mean cohen_pss_score is around 11

# path coefficients for twin 1
pathAceT1 <- mxPath( from=c("A1","C1","E1"), to="cohenpss_score_A", arrows=1,
free=TRUE, values=.5, label=c("a","c","e") )
# path coefficients for twin 2
pathAceT2 <- mxPath( from=c("A2","C2","E2"), to="cohenpss_score_B", arrows=1,
free=TRUE, values=.5, label=c("a","c","e") )

#The common environmental factors are the same for both twins (by definition), so setting the fixed correlation of 1 for C1 and C2

# covariance between C1 & C2
covC1C2 <- mxPath( from="C1", to="C2", arrows=2,
free=FALSE, values=1 )

# covariance between A1 & A2 in MZ twins
covA1A2_MZ <- mxPath( from="A1", to="A2", arrows=2,
free=FALSE, values=1 )
# covariance between A1 & A2 in DZ twins
covA1A2_DZ <- mxPath( from="A1", to="A2", arrows=2,
free=FALSE, values=.5 )

# Data objects for Multiple Groups
dataMZ <- mxData( observed=mzData, type="raw" )
dataDZ <- mxData( observed=dzData, type="raw" )

# Combine Groups
paths <- list( latVariances, latMeans, obsMeans,
pathAceT1, pathAceT2, covC1C2 )
modelMZ <- mxModel(model="MZ", type="RAM", manifestVars=selVars,
latentVars=aceVars, paths, covA1A2_MZ, dataMZ )
modelDZ <- mxModel(model="DZ", type="RAM", manifestVars=selVars,
latentVars=aceVars, paths, covA1A2_DZ, dataDZ )

obj <- mxFitFunctionMultigroup(c("MZ", "DZ"))
modelACE <- mxModel(model="ACE", modelMZ, modelDZ, obj )

# Run Model
fitACE <- mxRun(modelACE)
sumACE <- summary(fitACE)
```

Log in or register to post comments
Replied on Wed, 03/10/2021 - 14:45
Picture of user. AdminRobK Joined: 01/24/2014

Let me be sure I understand. Your present script biometrically decomposes the variance in your continuous DV. You want to also regress the DV onto your ordinal IV and some nuisance covariates. You further want to biometrically decompose the ordinal IV. Is all that correct?
Log in or register to post comments
Replied on Sun, 03/14/2021 - 17:18
No user picture. Leo Joined: 01/09/2020

How about a simple bivariate genetic model? You can set this up as a Cholesky with covariates. There are many scripts here, but if that’s all you want to do, you might want to consider umx. The function umxACE will take care of this. If you go down the openmx road, be sure to to check out joint continuous ordinal examples with definition variables.

Subsequently you can observe the e cross path and do some testing with it. I.e. does model fit decrease significantly if you fix the path to zero?

Log in or register to post comments
Replied on Thu, 03/18/2021 - 15:59
Picture of user. AdminRobK Joined: 01/24/2014

OK, I'll assume that the two "nuisance" covariates are age and sex. In that case, since the biometrical variance decompositions of both the DV and the IV are of interest, the IV will need to be adjusted for age and sex as well.

Give the attached script a try. Obviously, I can't actually test it, since I don't have your dataset.

Note that there are alternate ways in which the location and scale of the latent liability underlying the ordinal IV could be identified.

Edit: my mistake--you'll need to use `latentVars=c(aceVars,dummyVars)` when you construct `modelMZ` and `modelDZ`.

File attachments
210318.R
Log in or register to post comments
Replied on Tue, 03/30/2021 - 23:36
No user picture. sppspp10 Joined: 02/05/2021

Thank you very much for your helpful script. Any suggestions on the best way to fix a covariance matrix that is not a definite positive?

fitACE <- mxRun(modelACE)
Running ACE with 16 parameters
constraint-adjusted standard errors could not be calculated because the coefficient matrix of the quadratic form was uninvertibleIn model 'ACE' Optimizer returned a non-zero status code 10. Starting values are not feasible. Consider mxTryHard()Error: The job for model 'ACE' exited abnormally with the error message: fit is not finite (0: The continuous part of the model implied covariance (loc3) is not positive definite in data 'DZ.data' row 89. Detail:
covariance = matrix(c( # 6x6
0.7575, 0.37875, 0, 0, 0, 0
, 0.37875, 0.7575, 0, 0, 0, 0
, 0, 0, 0, 0, 0, 0
, 0, 0, 0, 0, 0, 0
, 0, 0, 0, 0, 0, 0
, 0, 0, 0, 0, 0, 0), byrow=TRUE, nrow=6, ncol=6)

1: The continuous part of the model implied covariance (loc3) is not positive definite in data 'DZ.data' row 82. Detail:
covariance = matrix(c( # 6x6
0.7575, 0.37875, 0, 0, 0, 0
, 0.37875, 0.7575, 0, 0, 0, 0
, 0, 0, 0, 0, 0, 0
, 0, 0, 0, 0, 0, 0
, 0, 0, 0, 0, 0, 0
, 0, 0, 0, 0, 0, 0), byrow=TRUE, nrow=6, ncol=6)

Log in or register to post comments
Replied on Wed, 03/31/2021 - 12:36
Picture of user. AdminRobK Joined: 01/24/2014

How is the model-expected covariance matrix 6x6? I thought there were only two endogenous variables per twin, so it should be 4x4.

At any rate, it appears some paths that should be free are instead fixed to zero.

Log in or register to post comments
  Login

News