[R-sig-ME] lme4 and mitools interface

[Paul Johnson]

On Sat, Mar 5, 2011 at 6:10 PM, Jessica Salvatore <salv0041 at umn.edu> wrote:
> Hello,
>
> I am attempting to combine parameter estimates for an lmer model run on 5
> multiply imputed datasets. I have tried this two ways: (1) Using
> with.imputationList from {mitools} and (2) Using a for loop. Both of these
> strategies failed, although for different reasons. Below I have created a
> toy example that reproduces the problems. I seem to get the most traction
> with the for-loop strategy; however, I am having trouble in the last
> combination step. Any help would be much appreciated.
>
> Thank you,
>
> Jessica Salvatore


Hello, Jessica

I didn't see anybody had answered your question.  Here is a working
example from a project I did 2 years ago. I generated the imputations
with the Amelia package, then this puts them in an imputationLIst and
runs through lmer.   The thing "myimpdat" is Amelia output, and I
first have to take the mX object out of each one and give to
imputationList.

I promise this did work, but I've not re-run it today to make sure it
still works

library(mitools)

myimplist <- imputationList(list(myimpdat$m1, myimpdat$m2,
myimpdat$m3, myimpdat$m4, myimpdat$m5, myimpdat$m6, myimpdat$m7,
myimpdat$m8, myimpdat$m9, myimpdat$m10))


###Do for whites, then blacks
mymodels1 <- with(myimplist, lmer ( (TR2COPN==0) ~     (AGERANGEN +
RMARITAL + (round(POLINTN > 0)) + EDUCN + GENDERN +  OWNN +  TRNEIN +
TROWNN + GRPETHN + QOLN + (round(EFFCOMN) > 1) + COPDIFF +(COPKILLB >
kdividers[2])  * BLKCHIEF + COPRATE + (CRIME98 > qdividers[2]) +
(BADCRIME > wdividers[2]) + (MYCTYBLK > 50)  + scale1)  + ( 1 |
PLACE33) ,  subset=RACE_ALL %in% c("White"),  family=binomial, x=T))



mymodels3 <- with(myimplist, lmer ( (TR2COPN==0) ~  -1 + RACE_ALL /
(AGERANGEN + RMARITAL + (round(POLINTN > 0)) + EDUCN + GENDERN +  OWNN
+  TRNEIN + TROWNN + GRPETHN + QOLN + (round(EFFCOMN) > 1) + COPDIFF
+(COPKILLB > kdividers[2])  * BLKCHIEF + COPRATE + (CRIME98 >
qdividers[2]) + (BADCRIME > wdividers[2]) + (MYCTYBLK > 50)  + scale1)
 + ( -1 + RACE_ALL | PLACE33) ,  family=binomial, x=T))



mymodels <- mymodels3
###mymodels <- blonlynocontext

betas <- MIextract(mymodels, fun=fixef)
vars <- MIextract(mymodels, fun=vcov)

vars2 <- list()
for ( i in 1:10) { vars2[[i]] <- as.matrix(vars[[i]])  }

mycombine <- MIcombine(betas, vars2)

mycombine

summary(mycombine)


>
> #load needed libraries
> library(lme4)
> library(mitools)
> library(mice)
>
> #Dataset has 4 variables (id, female, y.1, y.2, y.3 (i.e., the dependent
> #variable measured at three timepoints)
> id<-c(1, 2, 3, 4, 5, 6, 7, 8, 9, 10)
> female<-c(0, 0, 0, 0, 0, 1, 1, 1, 1, 1)
> y.1<-c(19, 7, 29, 16, 15, 30, 2, 20, 0, 26)
> y.2<-c(3, 28, NA, 1, 25, 2, NA, 22, 17, 23)
> y.3<-c(6, NA, 13, 20, 13, 2, 4, 30, NA, 22)
> data.wide<-cbind(id, female, y.1, y.2, y.3)
> data.wide<-as.data.frame(data.wide)
>
> #run model without imputing missing data (reshape, order, and eliminate
> missing data first)
> test.long<-reshape(data.wide, varying = c("y.1", "y.2", "y.3"), idvar =
> "id",
> v.names = "y", timevar = "year", times = c(1, 2, 3), direction = "long")
> test.long<-test.long[order(test.long$id, test.long$year),]
> test.long<-na.omit(test.long)
> test.model<-lmer(y ~ 1+ female + year + (1+year|id), data = test.long)
> summary(test.model)
>
> #multiply impute missing data 5 times
> nrep<-5  #number of imputations
> myseed<-1234  #set seed so that results can be reproduced
> set.seed (myseed)
> imp<-mice(data.wide[2:5], m = nrep, im = c(rep("", 2), rep("norm", 2)))
> #imputation procedure
>
> #reshape imputed datasets into long format. Change .id to numeric,
> rename.imp
> comp.data.wide<-complete(imp, action = "long")
> comp.data.wide$myid<-as.integer(comp.data.wide$.id)
> comp.data.wide$myimp<-comp.data.wide$.imp
> #order by imputation, id, and year
> data.long<-reshape(comp.data.wide, varying = c("y.1", "y.2", "y.3"), idvar =
> "id", v.names = "y", timevar = "year", times = c(1, 2, 3), direction =
> "long")
> data.long<-data.long[order(data.long$myimp, data.long$myid,
> data.long$year),]   #order by imputation, id, and year
>
> ###############################
> #METHOD 1: WITH.IMPUTATIONLIST
> ##############################
> #create imputationList
> impList<-imputationList(data.long[1:5])
> #run model using  with.ImputationList
> models<-with(impList, lmer(y ~ female + year + (1+year|myid)))
>    #When I run this, I get the following error message:
>    #Error in eval(expr, dataset, enclos = pf) :
>    #  numeric 'envir' arg not of length one
>
>
> ####################
> #METHOD 2: FOR LOOP
> ####################
>        LMM<- list(nrep)
>        LMM.beta<- list(nrep)
>        cov.beta<- list(nrep)
>        #Loop over replications
>        for(i in 1:nrep){
>            LMM[[i]]<-lmer(y~ 1 + female +year+ (1+year|myid), data =
> data.long, subset = myimp==i)
>            LMM.beta[[i]]<-fixef(LMM[[i]])
>            cov.beta[[i]]<- vcov(LMM[[i]])}
>        #Access results (i.e., show that the for loop works)
>        summary(LMM[[1]])
>        #Combine results
>        myresults <- MIcombine(LMM.beta,cov.beta)
>        smyresults <- summary(myresults)
>        #When I do this, I get the following error message
>        #Multiple imputation results:
>        #      MIcombine.default(LMM.beta, cov.beta)
>        #Error in y[1L + 0L:(m - 1L) * (n + 1L)] <- x :
>        #  incompatible types (from S4 to double) in subassignment type fix
>
>        #I get the same error when I use MIextract:
>        LMMbetas2<-MIextract(LMM, fun = fixef)
>        cov.betas2<-MIextract(LMM, fun = vcov)
>        summary(MIcombine(LMMbetas2, cov.betas2))
>
>
> --
> Jessica Salvatore
> Doctoral Candidate
> Institute of Child Development || University of Minnesota
> 51 East River Road || Minneapolis MN || 55455
> Office 166 || 612.624.5096
>
>        [[alternative HTML version deleted]]
>
> _______________________________________________
> R-sig-mixed-models at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>



-- 
Paul E. Johnson
Professor, Political Science
1541 Lilac Lane, Room 504
University of Kansas