2017-03-28 2 views
3

Je courais un modèle mixte en utilisant lme4 dans R:bootstrapping pour lmer avec terme d'interaction

full_mod3=lmer(logcptplus1 ~ logdepth*logcobb + (1|fyear) + (1 |flocation), 
data=cpt, REML=TRUE) 

Résumé:

Formula: logcptplus1 ~ logdepth * logcobb + (1 | fyear) + (1 | flocation) 
    Data: cpt 

REML criterion at convergence: 577.5 

Scaled residuals: 
    Min  1Q Median  3Q  Max 
-2.7797 -0.5431 0.0248 0.6562 2.1733 

Random effects: 
Groups Name  Variance Std.Dev. 
fyear  (Intercept) 0.2254 0.4748 
flocation (Intercept) 0.1557 0.3946 
Residual    0.9663 0.9830 
Number of obs: 193, groups: fyear, 16; flocation, 16 

Fixed effects: 
        Estimate Std. Error t value 
(Intercept)  4.3949  1.2319 3.568 
logdepth   0.2681  0.4293 0.625 
logcobb   -0.7189  0.5955 -1.207 
logdepth:logcobb 0.3791  0.2071 1.831 

J'ai utilisé le paquet effects et fonction R pour calculer 95 % intervalles de confiance pour la sortie du modèle. J'ai calculé et extrait l'IC 95% et l'erreur standard en utilisant le package effects afin que je puisse examiner la relation entre la variable prédictive d'importance et la variable de réponse en maintenant la variable prédictive secondaire (logdepth) constante à la médiane (2.5) en l'ensemble de données:

gm=4.3949 + 0.2681*depth_median + -0.7189*logcobb_range + 0.3791* 
(depth_median*logcobb_range) 

ef2=effect("logdepth*logcobb",full_mod3, 
    xlevels=list(logcobb=seq(log(0.03268),log(0.37980),,200))) 

Grand Mean model output with 95% CI

nous avons tenté d'amorcer les 95% d'IC ​​en utilisant le code de here. Cependant, j'ai besoin de calculer les IC à 95% pour seulement la profondeur médiane (2,5). Existe-t-il un moyen de spécifier dans le code confint() afin que je puisse calculer les CI nécessaires pour visualiser les résultats bootstrap comme dans le graphique ci-dessus?

confint(full_mod3,method="boot",nsim=200,boot.type="perc") 

Répondre

4

Vous pouvez le faire en spécifiant une fonction personnalisée:

library(lme4) 
?confint.merMod 

FUN: fonction d'amorçage; Si 'NULL', une fonction interne qui renvoie les paramètres à effet fixe ainsi que les paramètres à effet aléatoire sur l'échelle d'écart-type/de corrélation seront utilisés. Voir 'bootMer' pour plus de détails.

Alors FUN peut être une fonction de prédiction (?predict.merMod) qui utilise un argument newdata qui varie et fixe les variables prédictives appropriées.

Un exemple avec données intégrée (pas tout à fait aussi intéressant que le vôtre car il y a une seule variable explicative continue, mais je pense qu'il devrait illustrer l'approche assez clairement):

fm1 <- lmer(Reaction ~ Days + (Days | Subject), sleepstudy) 
pframe <- data.frame(Days=seq(0,20,by=0.5)) 
## predicted values at population level (re.form=NA) 
pfun <- function(fit) { 
    predict(fit,newdata=pframe,re.form=NA) 
} 
set.seed(101) 
cc <- confint(fm1,method="boot",FUN=pfun) 

Image:

par(las=1,bty="l") 
matplot(pframe$Days,cc,lty=2,col=1,type="l", 
    xlab="Days",ylab="Reaction") 

enter image description here