Glm quasibinomial, sélection de modèle - glm quasibinomial vs glm binomial

Bonjour à tous,

Grace à une régression logistique, je cherche à identifier parmi différents paramètres environnementaux, ceux qui influencent ma variable à expliquer, qui est de type binaire.

J'ai utilisé un glm quasibinomial. Pourriez-vous m'indiquer dans ce cas, le ou les moyens à privilégier afin de sélectionner le meilleur modèle possible ? On utilise généralement l'AIC me semble-t-il. Cependant, celui-ci n'est pas calculé dans le cas d'un glm quasibinomial. Doit-on dans ce cas, comparer simplement les différences de déviances entre chaque couple de modèle ?
De plus, en utilisant les mêmes termes mais dans un ordre différent je ne trouve pas les mêmes paramètres significatifs. J'ai d'autant plus de mal à choisir les paramètres à conserver dans mon modèle... Dans l'idéal, je cherche une fonction R qui permettrait cette sélection, via la méthode "pas à pas" par exemple (afin d'éviter les inconvénients des méthodes ascendantes et descendantes), mais qui, du coup, ne se baserait pas sur les AIC. Est-ce possible ?

J'ai une seconde question... Certains auteurs suggèrent que la surdispersion est quasiment la règle dans les régressions logistiques et qu'il vaut mieux systématiquement en tenir compte, plutôt que de la négliger alors qu'elle était bien présente. J'ai lu de plus, que les résultats ne seraient de toute façon pas affectés s'il n'y avait en réalité, pas de surdispersion. Est-ce gênant dans ce cas de passer spontanément à un glm quasibinomial sans avoir essayer le glm binomial et même si les données n'indiquaient pas particulièrement une surdispersion à la base ?

Merci beaucoup pour votre aide,

Ginger

Alors je vais essayer de répondre à une partie de tes questions.

Pour commencer le glm quasibinomial est en faite un glm binomiale qui ne fixe pas le paramètre dispersion, ainsi si ton modèle est surdispersé ce paramètre est calculé en conséquence, mais si il n'y a pas de surdispersion cela ne change a priori rien. Je serai à ta place j'utiliserai donc le modèle quasibinomial directement (ça t'économise du temps de calcul).

Je ne sais pas combien tu as de variables et observations à étudier dans ton modèle mais la fonction step() de R permet de sélectionner le meilleur modèle en fonction du critère AIC, et tu as le choix de faire la méthode descendante, montante, ou les deux à la fois. Attention : si tu as un nombre de données conséquente et que tu fais tourner ton programme en local, les temps de calcul peuvent être très long, en effet, la fonction teste tout les modèles possibles...

Par ailleurs le critère AIC est le critère de sélection de modèle généralement utilisé. Mais il en existe d'autres : BIC, cp de Mallows, etc...

Ensuite si tu ne veux pas utiliser AIC, tu peux utiliser la fonction anova() sur le résultat de ta fonction glm() qui te permet de tester la significativité des coefficients associés à tes variables explicatives selon différents tests (Fisher, chi-2,...).
Par contre, il n'existe pas, à ma connaissance, de fonction retirant progressivement les variables dont le coefficient n'est pas significatif (au seuil 5% par défaut), mais tu peux coder une fonction qui le ferait.
Je pense que cette solution est moins gourmande en temps de calcul vu que tu recalcule juste un glm à chaque itération.

Pour finir, tu n'en as pas parlé, mais n'oublie pas que tu peux prendre en compte, ou pas, les effets croisés de tes variables explicatives, ce qui complique le problème (tu passe de par exemple 10 variables explicatives à 100 si tu prends en compte les effets croisés) mais peut, parfois, être intéressant.

Dernier conseil : évite de mélanger les différentes méthodes pour sélectionner le meilleur modèle, cela complique les interprétations en général et n'a pas forcément de "logique statistique"...

Bon courage ! ;)

Bonjour,

merci pour cette réponse qui m'éclaire un peu plus...

Pour la fonction step, par contre, j'ai un problème pratique. En faisant quelques recherches, j'ai pensé à utiliser cette fonction mais dans mon cas, le message d'erreur suivant s'affiche :

"> step(modele_d, direction="both")
Erreur dans step(modele_d, direction = "both") : 
 AIC n'est pas défini pour ce modèle, donc 'step' ne peut poursuivre"

Cela correspond bien à ce que j'ai remarqué et lu quelque part également... En utilisant un modèle quasibinomial, il n'y aurait pas de calcul d'AIC !?? (R indique "AIC : NA"). Du coup, je ne sais pas sur quoi me baser à la place ? Et quelle autre fonction réaliserait cette même démarche pour ce cas particulier ?

Alors je ne sais pas s'il existe des fonctions toutes prêtes pour ton problème alors (je pense que oui mais je n'en connais pas), mais tu peux faire une fonction à la main comme je t'ai indiqué précédemment à l'aide de la fonction anova() et en sélectionnant les variables dont les coefficients sont significatifs (au seuil alpha que tu fixeras toi même, en général 5%). Tu t'arrêtes lorsque tout les coefficients sont significatifs.

Merci, je vais essayer cette méthode !

Reviens nous dire ici si cela a marché, à mon avis cela doit être bon, mais je suis curieux :D

Bonjour,

J'ai utilisé la méthode anova comme recommandé... D'abord pour tester la significativité des paramètres intégrés dans le modèle : anova(glm, test="F"), puis pour comparer les modèles successifs : anova (glm1, glm2, test="F" ou test= Chi2) afin de savoir si la différence de déviance était significative ou non. Lorsque la différence n'était pas significative, j'ai choisi le modèle le plus simple (un paramètre en moins) pour respecter le principe de parcimonie, et lorsqu'elle était significative j'ai choisi le modèle avec la meilleure déviance.

Pour tester la significativité des paramètres, j'ai utilisé le test F (plutôt que le Chi2) suite à la lecture d'un cours trouvé sur le net :

Code:

---

1 2 3

anova(z, test = "Chi") # Analysis of deviance - sequential anova(z, test = "F") # Use F test if using gaussian, quasibinomial or quasipoisson link# anova(z1, z2, test = "Chi") # compare fits of 2 models, "reduced" vs "full"

---

Cependant je n'ai pas encore pu conclure car selon l'ordre de mes variables dans le modèle, les variables considérées comme significatives changent légèrement. Ainsi, dans un cas, 3 de mes variables (var2, var3, var4) sont sélectionnées dans le modèle final tandis que dans un autre, seules 2 d'entre elles sont sélectionnées (var2, var3).
Pour essayer de m'affranchir de ce problème, j'ai essayé quelques packages dont j'ai entendu parler. Ils permettraient d'automatiser la sélection (packages MuMIn et glmulti)...mais sans succès. J'ai eu des messages d'erreur et n'ai pas réussi à obtenir de résultat final, même après de nombreux essais. Je me cantonnerai donc plutôt à l'anova, plus facile à gérer pour un niveau moyen en stats et programmation comme le mien. Mais alors, comment trancher sur mon modèle final entre les 2 options que j'ai trouvées ?

J'ai quelques autres questions secondaires :
- le glm semble bien fonctionner (je n'ai pas non plus de message d'erreur type : glmfit : l'algorithme n'a pas convergé), pourtant mes graphs des résidus sont assez "moches". Dans ce cas, mon modèle reste-t-il valable ?
- Dois-je intégrer les interactions entre variables dès le début dans la sélection du modèle ou puis-je les tester une fois que mes variables significatives sont sélectionnées ?

Merci beaucoup pour ton aide parce que à ce stade de détail je rame carrément !

Ginger

Bonjour,

Tu as bien fais d'utiliser la fonction anova() si tu n'avais pas d'autre moyen.
Par contre, il est très étrange que tu ne retrouves pas le même modèle suivant l'ordre des variables. As- tu bien enlevé à chaque itération la variable la moins significative, jusqu'à que celles-ci soient toutes significatives au seuil désiré ?

Pour ce qui est du warning sur la convergence de l'algorithme c'est déjà un bon point, sinon tes p-valeurs risquent de ne pas avoir de sens.
Cependant, ce n'est pas parce que ton glm converge que tu as un "bon" modèle. En effet, as-tu vérifié la distribution de tes variables avant de lancer ton modèle linéaire généralisé ? Est-ce que ta fonction lien correspond bien à ce que tu as pu voir dans une première étude exploratoire des données ?

Il est vrai que tes résidus semblent particulièrement moches ! Mais ceci n'est certainement pas un problème. En effet, le graphe des résidus est utile lorsque l'on fait un glm classique (c'est à dire que l'on considère une loi normale), ici tu utilise une modèle quasibinomial, les résidus ne suivent donc pas une lois normale centrée réduite mais ceci doit être normale.
Par contre, je ne connais pas les méthodes pour valider (ou invalider) ton modèle dans ce cas ci.

Enfin, il peut tout de même être utile (si ce n'est pas trop gourmand en temps de calcul) de regarder les interactions, tu as parfois un meilleur modèle en prenant en compte les interactions. Par contre si tu prends en compte le modèle avec interaction il faut que tu effectues ta méthode dès le début, et ne prend pas en compte que les variables sélectionnées dans ton modèle simple. L'interaction de deux variables peut être significative alors qu'une des deux variables ne l'est pas.

Après tout dépend si tu cherches un modèle précis ou simple, un modèle d'ajustement ou de prédiction.

Bon courage !

Pas de doute concernant ma variable y, elle est binaire, d'où le lien logit de la famille quasibinomiale. Par contre je n'ai pas particulièrement vérifié la distribution de mes variables x. Certaines sont discrètes, d'autres continues... Aurais-je dû ?

Concernant les interactions, modèle d'ajustement ou de prédiction ? Je dois dire que les termes ne m'évoquent pas grand chose. Je cherche un modèle simple mais qui expliquera au mieux ma variable y afin de connaître les paramètres explicatifs en étant au plus près de la "réalité" . Cependant lorsque j'ai entré dans R le modèle glm(y~x1*x1......*x6), le logiciel n'a pas beaucoup apprécié. J'ai abandonné le test au bout de quelques heures du coup...

Citation:

---

Je cherche un modèle simple mais qui expliquera au mieux ma variable y afin de connaître les paramètres explicatifs en étant au plus près de la "réalité"

---

Si je comprends bien, tu cherches à faire en sorte que ton modèle s'ajuste au mieux à tes données, en comparaison avec un modèle de prédiction où tu peux tolérer que ton modèle s'ajuste moins bien à tes données si il prédit bien.

En général, toute étude statistique est précédé par une étape exploratoire de tes données. Il faut donc que tu étudies ta variable à expliquer MAIS aussi tes variables explicatives. Par exemple tu peux avoir deux variables explicatives très fortement corrélés entre elles qui apporte donc la même information en double, tu peux donc en enlever une des deux dans ton modèle.

Ce n'est pas obligatoire mais il me semble qu'essayer de transformer tes variables explicatives pour les rendre plus "normales" peut aussi apporter une sorte de stabilité à ton modèle. Là il faudrait chercher dans des cours sur le net. Mais même sans ça tu peux tout de même faire un modèle linéaire généralisé.

Pour finir je ne comprends pas pourquoi tu as utilisé la fonction glm(y~x1*x1......*x6) ?
Normalement tu utilises la fonction glm() en précisant la fonction lien à l'aide du paramètre family de cette fonction.
Enfin, si tu veux prendre en compte les interactions il faut que tu fasses glm(y~(x1*x1......*x6)^2, family = ...)

Bon courage ! :D