Download
1 / 22
220 likes | 403 Views
OpenMx. Frühling Rijsdijk. OpenMx Models. Steven Boker1 Michael Neale2 Hermine Maes2 Paras Mehta3 Michael Wilde4 Timothy Brick1 Jerey Spies1 Michael Spiegel1 Ryne Estabrook1 Sarah Kenny4 John Fox5 Timothy Bates6 1University of Virginia; 2Virginia Commonwealth University;
E N D
OpenMx Frühling Rijsdijk
OpenMx Models Steven Boker1 Michael Neale2 Hermine Maes2 Paras Mehta3 Michael Wilde4 Timothy Brick1 Jerey Spies1 Michael Spiegel1 Ryne Estabrook1 Sarah Kenny4 John Fox5 Timothy Bates6 • 1University of Virginia; • 2Virginia Commonwealth University; • 3University of Houston; • 4University of Chicago, Argonne National Labs; • 5McMasters University; 6University of Edinburgh
Open Mx is OpenMx is free open source software for fitting Structural Equation Models (SEM) to observed data. OpenMx offers the features you would expect in an SEM software package, but OpenMx works in ways that will make your modeling jobs easier and will allow you to do things that other SEM packages don't. Integration with R OpenMx works as an integral part of the R statistical software system. You have available the full power of the R statistical software system for data manipulation, graphics, simulation, and report generation http://openmx.psyc.virginia.edu
Open Mx is 1. A free, full-featured, open source SEM package. 2. Runs on Windows, Mac OS-X, and Linux. 3. Runs inside the R statistical programming environment. OpenMx features: 1. A new approach to model specification. 2. Allows both path-style and matrix-style scripting. 3. Web-based forums, tutorials, and a wiki. http://openmx.psyc.virginia.edu
Install OpenMx What computers run OpenMx? • You can run OpenMx on computers using Windows XP, Windows Vista, Mac (Intel or PPC) OS-X 10.5 or later, and varieties of Linux (32 and 64 bit). What do I need to do first? • In order to install OpenMx, you will need R version 2.9.x or 2.10.1. How do I install OpenMx? • Open up an R session and copy the following line into the R command line: • source('http://openmx.psyc.virginia.edu/getOpenMx.R') • execute line • A few lines of R output will scroll by and your OpenMx (and snow and Matrix) library will be installed
Biometrical Genetic Theory Observed Data model building Twin Model Summary Statistics Path Diagrams System of Linear Equations Covariance Algebra Path Tracing Rules Predicted Var/Cov of the Model Observed Var/Cov of the Data SEM
Model for MZ or DZ Pairs Reared Together 1 1 /.5 C A C E E A 1 1 1 1 1 1 e a e c a c PTwin 1 PTwin 2
Predicted Var-Cov Matrices Tw1 Tw2 Tw1 Tw2 Tw1 Tw2 Tw1 Tw2
Saturated Twin Model To get Twin correlations (MZ and DZ) Testing e.g. equality in means and variances across Twin 1 and 2 and MZ and DZ groups
Mean1 Mean2 Variance / Covariance Stand 2x2 r L1 L2 S1 S2 V1tw1 V1tw2 Diag 2x2
r L1 L2 S1 S2 V1tw1 V1tw2 Mean1 Mean2 Variance / Covariance * * V1T1 V1T2 S12 S1 * r * S2 V1T1 S22 V1T2 S1 * r * S2
Means and Variances Multinormal Probability Density Function: -|2πΣ |-n/2 e -.5((xi - μ) Σ-1 (xi - μ)’) make use of all available data get unbiased estimates if missing data are missing at random Use Maximum Likelihood to estimate free Parameters: 2 means, 2 variances, 1 covariance
ACEuniv.R require(OpenMx) source("GenEpiHelperFunctions.R") # ----------------------------------------------------------------------- # Prepare Data # ----------------------------------------------------------------------- NCdata <- read.table ('N-CortisolNA.csv', header=T, sep=',') names (NCdata) str(NCdata) nv <- 1 ntv <- nv*2 Vars <-('ncomp') selVars <- c('ncomp1','ncomp2') summary(NCdata) mzData <- subset(NCdata, zyg==1, selVars) dzData <- subset(NCdata, zyg==2, selVars) summary(mzData) summary(dzData) colMeans(mzData,na.rm=TRUE) cov(mzData,use="complete") colMeans(dzData,na.rm=TRUE) cov(dzData,use="complete")
Mean1 Mean2 Specify and Run Saturated Model with RawData and Matrix-style Input twinSatModel <- mxModel("twinSat", mxModel("MZ", mxMatrix(type="Full", nrow=1, ncol=ntv, free=T, values=c(30,30), name="expMeanMZ"), mxMatrix(type="Diag", nrow=ntv, ncol=ntv, free=TRUE, values=8, lbound=.001, name="expSDMZ" ), mxMatrix(type="Stand", nrow=ntv, ncol=ntv, free=TRUE, values=.8, lbound=-.99, ubound=.99, name="expCorMZ" ), mxAlgebra( expression= expSDMZ %*% expCorMZ %*% expSDMZ, name="expCovMZ" ), mxData(mzData, type="raw"), mxFIMLObjective("expCovMZ", "expMeanMZ", selVars)), mxModel("DZ", mxMatrix("Full", 1, 2, T, 30, name="expMeanDZ"), mxMatrix(type="Diag", nrow=ntv, ncol=ntv, free=TRUE, values=8, lbound=.001, name="expSDDZ" ), mxMatrix(type="Stand", nrow=ntv, ncol=ntv, free=TRUE, values=.6, lbound=-.99, ubound=.99, name="expCorDZ" ), mxAlgebra( expression= expSDDZ %*% expCorDZ %*% expSDDZ, name="expCovDZ" ), mxData(dzData, type="raw"), mxFIMLObjective("expCovDZ", "expMeanDZ", selVars)), mxAlgebra(MZ.objective + DZ.objective, name="-2sumLL"), mxAlgebraObjective("-2sumLL")) twinSatFit <- mxRun(twinSatModel) r L1 L2 S1 S2 V1tw1 V1tw2
TESTS ASSUMPTION OF EQUAL VARIANCES # ------------------------------------------------------------------------------------------------------------------ # Specify and Run Saturated SubModel 1 equating Variances ACROSS TWINS # ---------------------------------------------------------------------------------------------------------------- twinSatModelSub1 <- twinSatModel twinSatModelSub1$MZ$expSDMZ <- mxMatrix("Diag", 2, 2, T, 8, "sdMZ", name="expSDMZ") twinSatModelSub1$DZ$expSDDZ <- mxMatrix("Diag", 2, 2, T, 8, "sdDZ", name="expSDDZ") twinSatFitSub1 <- mxRun(twinSatModelSub1) #Check!! ExpSDMZ <- mxEval(MZ.expSDMZ, twinSatFitSub1) ExpSDMZ ExpSDDZ <- mxEval(DZ.expSDDZ, twinSatFitSub1) ExpSDDZ
TESTS ASSUMPTION OF EQUAL VARIANCES # ------------------------------------------------------------------------------------------------------------------ # Specify and Run Saturated SubModel 1 equating Variances ACROSS TWINS and ZYG # ---------------------------------------------------------------------------------------------------------------- twinSatModelSub1 <- twinSatModel twinSatModelSub1$MZ$expSDMZ <- mxMatrix("Diag", 2, 2, T, 8, "sd", name="expSDMZ") twinSatModelSub1$DZ$expSDDZ <- mxMatrix("Diag", 2, 2, T, 8, "sd", name="expSDDZ") twinSatFitSub1 <- mxRun(twinSatModelSub1) #Check!! ExpSDMZ <- mxEval(MZ.expSDMZ, twinSatFitSub1) ExpSDMZ ExpSDDZ <- mxEval(DZ.expSDDZ, twinSatFitSub1) ExpSDDZ
C A E 1 1 a e c PTwin 1 SPECIFY AND RUN ACE MODEL WITH RAWDATA AND MATRIX STYLE INPUT univACEModel <- mxModel("univACE", mxModel("ACE", # Matrix for expected means vector for MZ and DZ twins mxMatrix("Full", 1, 2, T, 20, "mean", name="expMean"), # Matrices a11, c11, and e11 to store the a, c, and e path coefficients mxMatrix("Full", nrow=1, ncol=1, free=TRUE, values=3, label="a11", name="a"), mxMatrix("Full", nrow=1, ncol=1, free=TRUE, values=3, label="c11", name="c"), mxMatrix("Full", nrow=1, ncol=1, free=TRUE, values=3, label="e11", name="e"), # Matrices A, C, and E to compute variance components mxAlgebra(a * t(a), name="A"), mxAlgebra(c * t(c), name="C"), mxAlgebra(e * t(e), name="E"), # Algebra to compute total variances and SD mxAlgebra( expression=A+C+E, name="V" ), mxMatrix( type="Iden", nrow=nv, ncol=nv, name="I"), mxAlgebra( expression=solve(sqrt(I*V)), name="sd"), # Algebra to model expected variance/covariance matrix in MZ mxAlgebra(rbind (cbind(A+C+E , A+C), cbind(A+C , A+C+E)), name="expCovMZ"), # Algebra to model expected variance/covariance matrix in DZ mxAlgebra(rbind (cbind(A+C+E , .5%x%A+C), cbind(.5%x%A+C , A+C+E)), name="expCovDZ") ),
mxModel("MZ", mxData(mzData, type="raw"), mxFIMLObjective("ACE.expCovMZ", "ACE.expMean",selVars) ), mxModel("DZ", mxData(dzData, type="raw"), mxFIMLObjective("ACE.expCovDZ", "ACE.expMean",selVars) ), mxAlgebra(MZ.objective + DZ.objective, name="-2sumLL"), mxAlgebraObjective("-2sumLL") ) # end of model univACE univACEFit <- mxRun(univACEModel) # ----------------------------------------------------------------------- # Generate ACE Model Output # ----------------------------------------------------------------------- LL_ACE <- mxEval(objective, univACEFit) LL_ACE
univTwinSatFit univTwinSatFit@algebras univTwinSatFit@algebras$’-2sumll’ univTwinSatFit@submodels univTwinSatFit@submodels$MZ =univTwinSatFit$MZ univTwinSatFit$MZ@matrices univTwinSatFit$MZ@matrices$CholMZ =univTwinSatFit$MZ$CholMZ univTwinSatFit$MZ@algebras univTwinSatFit$MZ@algebras$expCovMZ =univTwinSatFit$MZ$expCovMZ univTwinSatFit$MZ$objective univTwinSatFit$MZ$data all information in MxModel fitted object list of algebras of container model specific algebra list of all items in all child models list of all items in specific child model list of all matrices in specific child model specific matrix in specific child model list of all algebras in specific child model specific algebra in specific child model objective of specific child model data of specific child model Extract Information
R Matrix Operators http://openmx.psyc.virginia.edu/wiki/matrix-operators-and-functions For example: Mx OpenMx Inverse ~ solve() Transpose ‘ t() Row x col * %*% Dot Prod . * Kronecker Prod @ %x% Vert adhesion _ rbind() Hor adhesion | cbind()
GenEpiHelperFunctions.R • Function "parameterSpecifations()" prints labels of a MxMatrix with square brackets surrounding free parameters; returns a matrix of strings • Function "expectedMeansCovariances()" prints expected means and expected covariance matrices for all submodels • Function "formatOutputMatrices()" prints matrix with specified labels and number of decimals • Function "formatMatrix()" returns a matrix with specified dimnames and of decimal places • Function "tableFitStatistics()" prints fit statistics with labels for Full Model and list of Nested Models
