Comment extraire / calculer l'effet de levier et les distances de Cook pour les modèles d'effets mixtes linéaires

10

Est-ce que quelqu'un sait comment calculer (ou extraire) l'effet de levier et les distances de Cook pour un merobjet de classe (obtenu via le lme4package)? Je voudrais les représenter pour une analyse des résidus.

r mixed-model linear-model residuals leverage Roey Angel
la source

Réponses:

13

Vous devriez jeter un oeil à l'ensemble de R influence.ME. Il vous permet de calculer des mesures de données influentes pour les modèles d'effets mixtes générés par lme4.

Un exemple de modèle:

library(lme4)
model <- lmer(mpg ~ disp + (1 | cyl), mtcars)

La fonction influenceest la base de toutes les étapes suivantes:

library(influence.ME)
infl <- influence(model, obs = TRUE)

Calculez la distance de Cook:

cooks.distance(infl)

Tracer la distance de Cook:

plot(infl, which = "cook")

entrez la description de l'image ici

Sven Hohenstein
la source
Merci! Cela aide certainement. Que diriez-vous de calculer l'effet de levier pour la distance d'un cuisinier par rapport à l'effet de levier?
Roey Angel
@RoeyAngel Je suppose que ce n'est pas possible avec le influence.MEpaquet. Malheureusement, je n'ai pas de solution pour cette tâche.
Sven Hohenstein
Ne devrait-il pas en être ainsi infl <- influence(model, group = "cyl"), car vous avez spécifié l'effet aléatoire comme (1|cyl)? Je ne sais pas, je ne comprends pas du tout, j'ai juste installé de l'influence ... mais je ne sais pas vraiment quand l'utiliser obs = TRUEet quand l'utiliser group...
Curieux
Je voudrais ajouter ce qui suit: Si vous souhaitez obtenir le numéro de ligne sur lequel les distances D de Cook se produisent - le même nombre apparaissant dans le tracé sans tracer, vous pouvez utiliser la formule r suivante sur les nombres de distances D de Cooks avec une coupe par exemple 0,1cooksD_data<-as.data.frame(cooks.distance(ft1)) cooksD_data_select<-cooksd[cooksD_data>0.1,drop=FALSE,] cooksD_oultiers<-as.numeric(rownames(cooksD_data_select))]
Elias Estatistics
Est-ce mieux que la hatvalues()fonction recommandée ici ?
Curieux