Étant donné trois variables, y
et x
qui sont positives continues et z
qui sont catégoriques, j'ai deux modèles candidats donnés par:
fit.me <- lmer( y ~ 1 + x + ( 1 + x | factor(z) ) )
et
fit.fe <- lm( y ~ 1 + x )
J'espère comparer ces modèles pour déterminer quel modèle est le plus approprié. Il me semble que dans un certain sens fit.fe
est imbriqué à l'intérieur fit.me
. En règle générale, lorsque ce scénario général se vérifie, un test du chi carré peut être effectué. Dans R
, nous pouvons effectuer ce test avec la commande suivante,
anova(fit.fe,fit.me)
Lorsque les deux modèles contiennent des effets aléatoires (générés par à lmer
partir du lme4
package), la anova()
commande fonctionne correctement. En raison des paramètres aux limites, il est normalement conseillé de tester la statistique du chi carré résultant par simulation, néanmoins, nous pouvons toujours utiliser la statistique dans la procédure de simulation.
Lorsque les deux modèles ne contiennent que des effets fixes, cette approche --- et, la anova()
commande associée --- fonctionnent correctement.
Cependant, lorsqu'un modèle contient des effets aléatoires et que le modèle réduit ne contient que des effets fixes, comme dans le scénario ci-dessus, la anova()
commande ne fonctionne pas.
Plus précisément, j'obtiens l'erreur suivante:
> anova(fit.fe, fit.me)
Error: $ operator not defined for this S4 class
Y a-t-il quelque chose de mal à utiliser l'approche Chi-Square d'en haut (avec simulation)? Ou est-ce simplement un problème de anova()
ne pas savoir comment gérer les modèles linéaires générés par différentes fonctions?
En d'autres termes, serait-il approprié de générer manuellement la statistique du chi carré dérivée des modèles? Si oui, quels sont les degrés de liberté appropriés pour comparer ces modèles? Par mes calculs:
Nous estimons deux paramètres dans le modèle à effets fixes (pente et interception) et deux autres paramètres (paramètres de variance pour la pente aléatoire et interception aléatoire) dans le modèle à effets mixtes. En règle générale, le paramètre d'interception n'est pas compté dans le calcul des degrés de liberté, ce qui implique que et ; cela dit, je ne sais pas si les paramètres de variance pour les paramètres à effets aléatoires doivent être inclus dans le calcul des degrés de liberté; les estimations de la variance pour les paramètres à effet fixe ne sont pas prises en compte , mais je pense que c'est parce que les estimations des paramètres pour les effets fixes sont supposées être des constantes inconnues alors qu'elles sont considérées comme des variables aléatoires inconnaissablespour des effets mixtes. J'apprécierais une aide sur cette question.
Enfin, quelqu'un a-t-il une R
solution plus appropriée ( basée sur) pour comparer ces modèles?
la source
lm()
par àgls()
partir dunlme
package, etlmer()
parlme()
(à nouveau à partir dunlme
package), tout fonctionnera correctement. Mais notez que vous obtiendrez un test conservateur ( valeurs de p trop grandes ), car les paramètres du modèle plus simple se trouvent à la limite de l'espace des paramètres. Et vraiment, le choix d'inclure ou non les effets aléatoires devrait être basé sur la théorie (par exemple, le plan d'échantillonnage), pas sur un test statistique.Réponses:
Techniquement, vous pouvez le faire fonctionner en changeant simplement l'ordre des paramètres:
Fonctionnera très bien. Si vous passez un objet généré par
lmer
first, leanova.merMod
sera appelé à la place deanova.lm
(qui ne sait pas comment gérer leslmer
objets). Voir:Bien que le choix d'un modèle mixte ou d'un modèle fixe soit un choix de modélisation qui doit prendre en compte la conception expérimentale, pas un problème de sélection de modèle. Voir https://bbolker.github.io/mixedmodels-misc/glmmFAQ.html#testing-significance-of-random-effects @ BenBolker pour plus de détails:
la source