Multivariate Means Var Model - help

2 replies [Last post]
pgseye's picture
Offline
Joined: 10/13/2009

Hi,

I'm unfamiliar with Mx and am still learning SEM concepts - apologies if I'm missing something fundamental with the following.

I've tried adapting a script (one of Hermines) used in the recent twins course for comparing means and variances for a multivariate trait across right and left eyes.

I thought commenting out the 4 MxAlgebra functions across twin order and zygosity would allow me to compare grand means and variances across R and L, but I receive an error message saying that an argument is missing.

Can someone tell me where I'm going wrong?

Thanks heaps,

Paul Sanfilippo

> require(OpenMx)
> require(psych)
> source("http://www.vipbg.vcu.edu/~vipbg/Tc24/GenEpiHelperFunctions.R")
> setwd("~/Dropbox/Disc_Cohorts/Twins")
> # Prepare Data
> # -----------------------------------------------------------------------
> load("All Twins Data_Wide.RData")
> attach(twinRvsL)
> describe(twinRvsL)
var n mean sd median trimmed mad min max range skew kurtosis se
Twin1or2 1 1858 1.5 0.50 1.5 1.5 0.74 1.00 2.00 1.00 0.00 -2.00 0.01
L_PC1 2 1854 0.0 0.03 0.0 0.0 0.03 -0.15 0.14 0.29 -0.25 0.65 0.00
L_PC2 3 1854 0.0 0.02 0.0 0.0 0.02 -0.11 0.09 0.20 -0.07 1.09 0.00
L_PC3 4 1854 0.0 0.01 0.0 0.0 0.01 -0.05 0.05 0.10 0.07 1.01 0.00
L_PC4 5 1854 0.0 0.01 0.0 0.0 0.01 -0.04 0.05 0.09 0.04 0.58 0.00
L_PC5 6 1854 0.0 0.01 0.0 0.0 0.01 -0.02 0.02 0.04 -0.05 0.10 0.00
L_PC6 7 1854 0.0 0.01 0.0 0.0 0.01 -0.02 0.02 0.04 0.07 0.22 0.00
R_PC1 8 1854 0.0 0.03 0.0 0.0 0.03 -0.14 0.12 0.26 -0.28 0.37 0.00
R_PC2 9 1854 0.0 0.02 0.0 0.0 0.02 -0.11 0.09 0.20 -0.15 1.21 0.00
R_PC3 10 1854 0.0 0.01 0.0 0.0 0.01 -0.04 0.06 0.10 0.28 1.21 0.00
R_PC4 11 1854 0.0 0.01 0.0 0.0 0.01 -0.05 0.04 0.09 -0.08 0.51 0.00
R_PC5 12 1854 0.0 0.01 0.0 0.0 0.01 -0.02 0.02 0.05 -0.03 0.29 0.00
R_PC6 13 1854 0.0 0.01 0.0 0.0 0.01 -0.02 0.02 0.04 -0.09 0.11 0.00
> Vars1 <- c('L_PC1','L_PC2','L_PC3','L_PC4','L_PC5','L_PC6')
> Vars2 <- c('R_PC1','R_PC2','R_PC3','R_PC4','R_PC5','R_PC6')
> nv <- 6
> ntv <- 6
> leftData <- twinRvsL[,2:7]
> rightData <- twinRvsL[,8:13]
> # DO NOT RUN! Fit Multivariate Saturated Model
> # -----------------------------------------------------------------------
> meanSV <- c(0,0,0,0,0,0)
> cholSV <- c(
+ 1,0,0,0,0,0,
+ 1,0,0,0,0,
+ 1,0,0,0,
+ 1,0,0,
+ 1,0,
+ 1)
> multiTwinSatModel <- mxModel("multiTwinSat",
+ mxModel("LEFT",
+ mxMatrix( type="Lower", nrow=ntv, ncol=ntv, free=TRUE, values=cholSV, lbound=-2, ubound=2, name="CholLEFT" ),
+ mxAlgebra( expression=CholMZ %*% t(CholMZ), name="expCovLEFT" ),
+ mxMatrix( type="Full", nrow=1, ncol=ntv, free=TRUE, values=meanSV, lbound=-0.5, ubound=0.5, name="expMeanLEFT" ),
+ mxData( observed=leftData, type="raw" ),
+ mxFIMLObjective( covariance="expCovLEFT", means="expMeanLEFT", dimnames=Vars1),
+ # Algebra's needed for equality constraints
+ mxAlgebra( expression=t(diag2vec(expCovMZ)), name="expVarLEFT"),
+ # mxAlgebra( expression=expVarMZ[1,1:nv], name="expVarMZt1"),
+ # mxAlgebra( expression=expVarMZ[1,(nv+1):ntv], name="expVarMZt2"),
+ # mxAlgebra( expression=expMeanMZ[1,1:nv], name="expMeanMZt1"),
+ # mxAlgebra( expression=expMeanMZ[1,(nv+1):ntv], name="expMeanMZt2")
+ ),
+ mxModel("RIGHT",
+ mxMatrix( type="Lower", nrow=ntv, ncol=ntv, free=TRUE, values=cholSV, lbound=-2, ubound=2, name="CholRIGHT" ),
+ mxAlgebra( expression=CholDZ %*% t(CholDZ), name="expCovRIGHT" ),
+ mxMatrix( type="Full", nrow=1, ncol=ntv, free=TRUE, values=meanSV, lbound=-0.5, ubound=0.5, name="expMeanRIGHT" ),
+ mxData( observed=rightData, type="raw" ),
+ mxFIMLObjective( covariance="expCovRIGHT", means="expMeanRIGHT", dimnames=Vars2),
+ # Algebra's needed for equality constraints
+ mxAlgebra( expression=t(diag2vec(expCovDZ)), name="expVarRIGHT"),
+ # mxAlgebra( expression=expVarDZ[1,1:nv], name="expVarDZt1"),
+ # mxAlgebra( expression=expVarDZ[1,(nv+1):ntv], name="expVarDZt2"),
+ # mxAlgebra( expression=expMeanDZ[1,1:nv], name="expMeanDZt1"),
+ # mxAlgebra( expression=expMeanDZ[1,(nv+1):ntv], name="expMeanDZt2")
+ ),
+ mxAlgebra( LEFT.objective + RIGHT.objective, name="-2sumll" ),
+ mxAlgebraObjective("-2sumll")
+ )
Error in mxModel("LEFT", mxMatrix(type = "Lower", nrow = ntv, ncol = ntv, :
argument is missing, with no default

pgseye's picture
Offline
Joined: 10/13/2009
Thanks for your help Mike.

Thanks for your help Mike. That did fix the problem.

I didn't realise about posting scripts + data - I'll do that in the future.

Thanks again.

Paul

neale's picture
Offline
Joined: 07/31/2009
Usually this error means that

Usually this error means that there's a stray comma. Probably it is because the last mxAlgebra statement you commented out does *not* have a comma at the end of the line, so it would probably be best to delete the comma at the end of the immediately preceding mxAlgebra statement:
change
mxAlgebra( expression=t(diag2vec(expCovDZ)), name="expVarRIGHT"),
into
mxAlgebra( expression=t(diag2vec(expCovDZ)), name="expVarRIGHT")

It's not all that easy to debug in this case because you have provided the error message feedback from R, rather than the R script itself and the data it uses. In the event that posting the dataset to the server is not appropriate (due to IRB, confidentiality or other concerns) then I suggest that a fake dataset be mocked up in order to help identify the source of the error.