[R-sig-ME] lme4: constrain sigma to 0 [SEC=Unclassified]

David Duffy David.Duffy at qimr.edu.au
Wed Feb 12 09:37:16 CET 2014


For those who would like to fit these models, you could muck around with 
the OpenMx package, which allows you to fix or constrain any parameters 
you like ;)  To fix the residuals to zero, change free to FALSE in the
"residual variances" bit.  Parameters with the same "labels" are 
automatically constrained equal.

#
# Based on an example in the OpenMX documentation
#
require(lme4)
data(sleepstudy)
sleep2 <- reshape(sleepstudy, direction="wide", idvar="Subject", 
timevar="Days")
rownames(sleep2) <- sleep2$Subject
sleep2 <- sleep2[, -1]
names(sleep2) <- gsub("\\.","",names(sleep2)) # Mx no like "."

require(OpenMx)
growthCurveModel <- mxModel("Sleepstudy as Linear Growth Curve Model",
     type="RAM",
     mxData(
         observed=sleep2,
         type="raw"
     ),
     manifestVars=names(sleep2),
     latentVars=c("intercept","slope"),
     # residual variances
     mxPath(
         from=names(sleep2),
         arrows=2,
         free=TRUE,
         values = rep(1, ncol(sleep2)),
         labels=rep("residual", ncol(sleep2))
     ),
     # latent variances and covariance
     mxPath(
         from=c("intercept","slope"),
         arrows=2,
                 connect="unique.pairs",
         free=TRUE,
         values=c(1, 1, 1),
         labels=c("vari", "cov", "vars")
     ),
     # intercept loadings
     mxPath(
         from="intercept",
         to=names(sleep2),
         arrows=1,
         free=FALSE,
         values = rep(1, ncol(sleep2))
     ),
     # slope loadings
     mxPath(
         from="slope",
         to=names(sleep2),
         arrows=1,
         free=FALSE,
         values=seq(0,9)
                              ),
     # manifest means
     mxPath(from="one",
         to=names(sleep2),
         arrows=1,
         free=FALSE,
         values = rep(0, ncol(sleep2))
     ),
     # latent means
     mxPath(from="one",
         to=c("intercept", "slope"),
         arrows=1,
         free=TRUE,
         values=c(1, 1),
         labels=c("meani", "means")
     )
) # close model

growthCurveFit <- mxRun(growthCurveModel)
print(summary(growthCurveFit))
print(growthCurveFit at output$estimate)



| David Duffy (MBBS PhD)
| email: David.Duffy at qimrberghofer.edu.au  ph: INT+61+7+3362-0217 fax: -0101
| Genetic Epidemiology, QIMR Berghofer Institute of Medical Research
| 300 Herston Rd, Brisbane, Queensland 4006, Australia  GPG 4D0B994A



More information about the R-sig-mixed-models mailing list