2013-09-25 1 views
3

J'adapte un tableau de contingence en 3 dimensions (non fourni ici mais je peux si cela peut aider) avec un modèle loglinear, à la fois avec loglm et avec glm. Les deux résultats que je reçois en termes de coefficients sont les suivants:coefficients dans GLM vs coefficients dans loglm

> coefficients(nodnox_loglm_model) 
$`(Intercept)` 
[1] 10.18939 

$w 
     0.05   0.1  0.15   0.2  0.25   0.3  0.35   0.4  0.45 
-1.04596513 -0.41193617 -0.08840858 0.06407334 -0.06862606 0.02999039 0.17084795 0.45838071 0.35307375 
     0.5 
0.53856982 

$s 
      2   3   4   5 
0.36697307 0.15164360 -0.48264571 -0.03597096 

et

> coefficients(nodnox_glm_model) 
(Intercept)   s3   s4   s5  w0.1  w0.15  w0.2  w0.25  w0.3 
    9.5104005 -0.2153295 -0.8496188 -0.4029440 0.6340290 0.9575566 1.1100385 0.9773391 1.0759555 
     w0.35  w0.4  w0.45  w0.5 
    1.2168131 1.5043458 1.3990389 1.5845350 

Je sais que ces deux méthodes ont différentes procédure numérique - Je ne me soucie pas de cela - tout ce que je veux savoir est comment puis-je relier les coefficients glm aux coefficients loglm?

Tout ce que je trouve sur Internet et la documentation que je cherchais avant de venir à stackoverflow est cette note:

Le tableau de coefficient de GLM fonctionne exactement comme le résumé pour ANOVA produit par lm: le niveau premier par ordre alphabétique (s2, w0.5) est utilisé comme interception, et tous les niveaux suivants sont testés par rapport au premier (donc les coefficients restants sont des différences par rapport à la moyenne, pas signifie eux-mêmes).

Pour moi, cependant, cela ne suffit pas pour comprendre comment obtenir les coefficients de la sortie de glm sous la forme de loglm. Maintenant, votre question pourrait être: "pourquoi ne pas utiliser loglm directement?" Loglm ne fonctionnerait pas dans mon cas (ce n'est pas celui que je compare ici, mais il a une table en 5 dimensions avec des zéros.) Si j'utilise loglm sur la table d'origine, ça me donne tous les coefficients comme NaNs) . Donc je suis coincé avec glm et je veux vraiment obtenir les coefficients comme dans loglm.

Merci beaucoup!

Répondre

4

Il semble que vous ayez une table croisée bidirectionnelle avec 10 niveaux de facteur w et 5 niveaux de facteur s sans interactions dans le modèle. Avec glm(), le schéma de codage par défaut pour les variables catégorielles est treatment coding où le premier groupe d'un facteur est le niveau de référence et le paramètre respectif de chaque groupe restant est sa différence avec cette référence. L'estimation (Intercept) est pour la cellule avec tous les groupes = niveau de référence pour leur facteur.

Avec loglm(), les paramètres sont pour le codage de déviation, ce qui signifie que chaque groupe a son propre paramètre et les paramètres pour un facteur somme à zéro. (Intercept) est le grand moyen qui s'ajoute à tous les effets de groupe.

Dans votre exemple, vous pouvez dire glm() d'utiliser le codage de déviation pour obtenir les mêmes estimations des paramètres comme avec loglm() (voir exemple ci-dessous), ou convertir les estimations des paramètres de codage de traitement comme suit:

  • w = 0,05 et s = 2 est la cellule de référence: glm() 9,5104005 = loglm() 10,18939 + -1,04596513 + 0,36697307
  • w = 0,1 et s = 2 est le niveau de référence pour s mais a besoin de la différence de w = 0,1 à la référence w = 0,05: glm() 9,5104005 + 0,6340290 = loglm() 10,188939 + -0,41193617 + 0,36697307
  • w = 0.1 et s = 3 mais nécessite la différence de w = 0,1 à la référence w = 0,05 et la différence de s = 3 à la référence s = 2: glm() 9,5104005 0,6340290 + + -0,2153295 = loglm() 10,18939 + -0,41193617 + 0,15164360, et etc.

Exemple avec glm() en utilisant la déviation de codage (UCBAdmissions est un tableau croisé dont les fréquences absolues intégrées dans la base R):

> library(MASS)        # for loglm() 
> llmFit <- loglm(~ Admit + Gender + Dept, data=UCBAdmissions) 
> coef(llmFit) 
$`(Intercept)` 
[1] 5.177567 

$Admit 
    Admitted Rejected 
-0.2283697 0.2283697 

$Gender 
     Male  Female 
0.1914342 -0.1914342 

$Dept 
      A   B   C   D   E   F 
0.23047857 -0.23631478 0.21427076 0.06663476 -0.23802565 -0.03704367 

> UCBdf <- as.data.frame(UCBAdmissions) # convert to data frame for glm() 
> glmFit <- glm(Freq ~ Admit + Gender + Dept, family=poisson(link="log"), 
+    contrasts=list(Admit=contr.sum, Gender=contr.sum, Dept=contr.sum), 
+    data=UCBdf) 
> coef(glmFit) 
(Intercept)  Admit1  Gender1  Dept1  Dept2  Dept3  Dept4 
5.17756677 -0.22836970 0.19143420 0.23047857 -0.23631478 0.21427076 0.06663476 
     Dept5 
-0.23802565 

Notez que glm() ne mentionne pas ces paramètres es les dates qui sont entièrement déterminées (aliasées) à travers la contrainte de somme à zéro pour les paramètres d'un facteur.

Questions connexes