[R-sig-ME] [R] using stimulate(model) for parametric bootstrapping in lmer repeatabilities
Søren Højsgaard
Soren.Hojsgaard at agrsci.dk
Wed Jun 8 10:34:47 CEST 2011
Dear Jenni,
In the newest version of the doBy package there is a function called PBrefdist (and PBrefdist.mer) for calculating the reference distribution of the likelihood ratio statistic for comparing nested models. Looking into this function may help you. Perhaps the functions PBmodcomp and BCmodcomp (and their .mer methods) for calculating p-values based on parametric bootstrap can also be of inspiration.
Regards
Søren
________________________________________
Fra: r-help-bounces at r-project.org [r-help-bounces at r-project.org] På vegne af Jenni Sanderson [Jennifer.Louise.Sanderson at googlemail.com]
Sendt: 8. juni 2011 09:56
Til: r-help at r-project.org; r-sig-mixed-models at r-project.org
Emne: [R] using stimulate(model) for parametric bootstrapping in lmer repeatabilities
Hi all,
I am currently doing a consistency analysis using an lmer model and
trying to use parametric bootstrapping for the confidence intervals.
My model is like this:
model<-lmer(y~A+B+(1|C/D)+(1|E),binomial)
where E is the individual level for consistency analysis, A-D are
other fixed and random effects that I have to control for.
Following Nakagawa and Scheilzeth I can work out the repeatability
estimate using the following (as it is a binomial and the residual
variance is fixed at 1).
attr(lme4::VarCorr(model)$E, "stddev")^2 /
(1*(pi^2)/3 + attr(lme4::VarCorr(model)$E, "stddev")^2 +
attr(lme4::VarCorr(model)$C, "stddev")^2 +
attr(lme4::VarCorr(model)$D, "stddev")^2 )
My question is can I use stimulate(model) to generate values that I
can then use to do parametric bootstrap analysis and generate the
confidence intervals?
Something like this:
n<-length(A)
niter<-1000
y<-matrix(nrow=n,ncol=niter*2)
for (i in 1:niter) {
y[,I(i*2-1):I(i*2)]<-simulate(model)[,1] }
rvalues<-numeric()
for (i in 1:niter) {
yboot<-cbind(y[,I(i*2-1)],y[,I(i*2)])
mboot<-lmer(y~A+B+(1|C/D)+(1|E),binomial)
rvalues[i]<- attr(lme4::VarCorr(mboot)$E, "stddev")^2/(1*(pi^2)/3 +
attr(lme4::VarCorr(mboot)$E, "stddev")^2 +
attr(lme4::VarCorr(mboot)$C, "stddev")^2 +
attr(lme4::VarCorr(mboot)$D, "stddev")^2 )}
confidence.intervals<-quantile(rvalues,c(0.05,0.95))
In the guide to lme4 it says that stimulate() "generate simulations
based on the estimated fitted models (conditional on the estimated
values of both the random and fixed effects)", which sounds like
exactly what I would need to generate values for parametric
bootstrapping but I can't find any examples of where anyone has done
this.
Any advice would be very much appreciated! Thank you very much.
Jenni
Jenni Sanderson
PhD student - Conflict and Cooperation in Vertebrate Societies
University of Exeter
______________________________________________
R-help at r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
More information about the R-sig-mixed-models
mailing list