[R-sig-ME] looping with multiple response variables lme4

Antoine Tremblay trea26 at gmail.com
Sun Apr 8 03:47:46 CEST 2012


Hi Justus,
Maybe the following functions from package "LMERConvenienceFunctions" 
can be of help:

bfFixefLMER_F.fnc       Back-fits an LMER model on p-values from ANOVA
                         and, optionally, on LLRT.
bfFixefLMER_t.fnc       Back-fits an LMER model on absolute t-value
                         and, optionally, on LLRT.
fitLMER.fnc             Back-fit fixed effects and forward-fit random
                         effects of an LMER mo

Basically you can have a model with say a 10-way interaction and all 
lower-level interactions and main effects) and these fucntions will 
iteratively remove any non-significant model term starting with the 
highest-order interaction and going down to the main effects. If a 
lower-order interaction or main effect is part of a higher-order one, 
the term won't be pruned (the marginality principle). It also prints a 
log of the steps it went through. Here's an example:

# fit initial model
m1 = lmer(Fz ~ FreqB * LengthB * WMC + (1 | Subject) +
        (1 | Item), data = eeg)

# backfit fixed effects on F
m2 = bfFixefLMER_F.fnc(mB, eeg, log.file = FALSE,
        llrt = FALSE)

processing model terms of interaction level 3
     iteration 1
         p-value for term "FreqB:LengthB:WMC" = 0.1152 > 0.05
         not part of higher-order interaction
         removing term
processing model terms of interaction level 2
     iteration 2
         p-value for term "FreqB:WMC" = 0.9552 > 0.05
         not part of higher-order interaction
         removing term
     iteration 3
         p-value for term "FreqB:LengthB" = 0.6564 > 0.05
         not part of higher-order interaction
         removing term
processing model terms of interaction level 1
     iteration 4
         p-value for term "WMC" = 0.7796 > 0.05
         part of higher-order interaction
         skipping term
     iteration 5
         p-value for term "FreqB" = 0.7635 > 0.05
         not part of higher-order interaction
         removing term
     iteration 6
         p-value for term "LengthB" = 0.0924 > 0.05
         part of higher-order interaction
         skipping term
pruning random effects structure ...
     nothing to prune


The output will simply be a pruned lme model. You can certainly look at 
the code and tweak it so it'll return SD, SE, and AIC.

Antoine
NeuroCognitive Imaging Laboratory
Dalhousie University
Halifax, Nova Scotia
Canada


 > Hi Justus:
 >
 > I think you forgot to surround the formula with quotes, like this:
 > stmLM <- as.formula(paste(Response[i],"~ccov + x5k_c + ltr + x5k_fr + 
 > elev + (1|P_ID) + (1|mtx) + (1|site)"))
 >
 > The error is basically saying you are trying to paste a string into a 
 > formula, which does not work.
 >
 > Hope this helps,
 > Roby
 >
 >
 > On Apr 5, 2012, at 9:29 PM, Justus Deikumah wrote:
 >
 >> Dear lme4 professionals
 >>
 >> Thank you for any help on the subject matter.
 >>
 >> I am doing a linear mixed effects modelling with 10 response 
variables and 5 predictors using lme4. I have 126 model combinations if 
I run my model with one response variable at a time. This seems very 
manual and will take the rest of my PhD. I have tried haven read some of 
your recent papers on loops and iterations so I have tried to write a 
sript for this with help from some friends as follows:
 > What I want the loop to do for me is to model each response with all 
the candidate predictors, drop one predictor and run the models again 
for all the 126 combinations. Then I want to produce Estimated 
coefficients of fixed effects, SD, SE, AICs etc as output that can be 
sent into an excel template for model averaging?
 >> attach(var1)
 >> Response <- c("gnst", "spst", "fvt", "Oph", "cniv", "frgv", "gniv", 
"inct","omnv","Sobs", "ACE", "Chao1","Chao2", "Abund")
 >
 > ## Output storage step up
 > results <- vector("list", length(Response))
 > ##Start the loop
 > #Create the formula string
 > #stmLM <- as.formula(paste(Response[i]~ccov + x5k_c + ltr + x5k_fr + 
elev + (1|P_ID) + (1|mtx) + (1|site)))
 >
 >
 >
 > lmer <- lmer(stmLM)
 >
 > # Save the results{
 > results[[i]] <- summary(lmer)
 > }
 > {
 > print(results[[i]])
 > }
 > My problem is anytime I reach summary and results portion of the 
analysis then i receive the error message below;
 >
 > Error: length(formula <- as.formula(formula)) == 3 is not TRUE
 > this is making my R-training uneasy for me. Can someone pls tell me 
what is wrong as I believe there may be a problem with the for loop of 
which i am new to anyway. What does the error message say? Has anyone 
any scripts for this??I am using the lmer function? What kind of object 
doest lmer produce?
 > Hope someone understands my question and can help, Please!!
 > Thank you
 >
 >
 >
 > Deikumah Justus
 > PhD Candidate (Confirmed) - Landscape Ecology & Conservation
 > Landscape Ecology and Conservation Group
 > School of Geography, Planning and Environmental Management
 > The University of Queensland
 > Brisbane Qld 4072
 > Australia
 >
 > Chamberlain Building (35)
 > Campbell Rd - St Lucia Campus
 >
 > W: http://www.gpem.uq.edu.au/cser-jdeikumah
 > E: j.deikumah2 at uq.edu.au
 > M: +61 424481796, +61 479094200
 > If you don't like something, Change it. If you can't change it, 
change your Attitude. DON'T COMPLAIN (Maya Angelou)
 >
 >     [[alternative HTML version deleted]]
 >
 > _______________________________________________
 > R-sig-mixed-models at r-project.org mailing list
 > https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models


-- 
Antoine Tremblay, PhD
NeuroCognitive Imaging Laboratory
Dalhousie University
Halifax, NS B3H 3J5,
Canada

Tel.: (902) 494-1911
eom




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