Comprendre la création de variables factices (manuelles ou automatisées) dans GLM

13

Si une variable de facteur (par exemple, le sexe avec les niveaux M et F) est utilisée dans la formule glm, des variables fictives sont créées et peuvent être trouvées dans le résumé du modèle glm avec leurs coefficients associés (par exemple, genderM)

Si, au lieu de compter sur R pour diviser le facteur de cette manière, le facteur est codé dans une série de variables numériques 0/1 (par exemple, genderM (1 pour M, 0 pour F), genderF (1 pour F, 0 pour M) et ces variables sont ensuite utilisées comme variables numériques dans la formule glm, le résultat du coefficient serait-il différent?

Fondamentalement, la question est: R utilise-t-il un calcul de coefficient différent lorsque vous travaillez avec des variables factorielles par rapport à des variables numériques?

Question de suivi (éventuellement répondue par ce qui précède): outre l'efficacité de laisser R créer des variables fictives, y a-t-il un problème avec les facteurs de recodage en tant que série de variables numériques 0,1 et en utilisant plutôt ceux du modèle?

Bryan
la source
2
Les coefficients seront les mêmes, le codage par défaut de R pour les facteurs est exactement comme vous l'avez décrit (c'est ce qu'on appelle le codage "factice"). Si vous effectuez cette opération manuellement, vous devrez connaître le «piège des variables factices» - vous ne pouvez pas inclure toutes les variables créées, mais seulement (alternativement, exclure l'interception); sinon votre modèle est sur-identifié. Il existe également d'autres méthodes de codage des variables de facteur. N - 1NN1
Affine
@ Affine ma supposition est que si je me nourrissais à la fois de genderM et de genderF, l'un d'eux retournerait NA pour les coefficients (avec le message qu'une variable a été exclue à cause des singularités). Cela a du sens car ils sont parfaitement linéaires dans ce cas. Mais vous dites que je ne peux pas inclure tous les N; cela signifie-t-il que, même si genderF est réglé sur NA, cela entraînerait des différences / problèmes pour le coefficient genderM? Ou, plus simplement, si GLM / LM exclut les variables en raison des singularités, l'utilisation d'un modèle sur-identifié est-elle un problème? (Je suis d'accord avec votre point de vue - juste en interrogeant les ramifications pratiques)
Bryan

Réponses:

22

Les variables catégorielles (appelées « facteurs » dans R) doivent être représentées par des codes numériques dans plusieurs modèles de régression. Il existe de très nombreuses façons de construire des codes numériques de manière appropriée (voir cette grande liste sur le site d'aide des statistiques de l'UCLA). Par défaut, R utilise le codage au niveau de référence (que R appelle "contr.treatment"), et qui est à peu près la statistique par défaut à l'échelle. Cela peut être modifié pour tous les contrastes de votre session R entière en utilisant les options? Ou pour des analyses / variables spécifiques en utilisant les contrastes ou le "C (notez la majuscule). Si vous avez besoin de plus d'informations sur le codage au niveau de référence, je l'explique ici: Régression basée par exemple sur les jours de la semaine.

Certaines personnes trouvent le codage de niveau de référence déroutant, et vous n'avez pas besoin de l'utiliser. Si vous le souhaitez, vous pouvez avoir deux variables pour les hommes et les femmes; cela s'appelle le niveau signifie le codage. Cependant, si vous faites cela, vous devrez supprimer l'ordonnée à l'origine ou la matrice du modèle sera singulière et la régression ne peut pas être ajustée comme les notes @Affine ci-dessus et comme je l'explique ici: le codage des variables qualitatives conduit à des singularités . Pour supprimer l'interception, vous modifiez votre formule en ajoutant -1ou +0comme ceci: y~... -1ou y~... +0.

Utiliser le codage de niveau signifie au lieu du codage de niveau de référence changera les coefficients estimés et la signification des tests d'hypothèse qui sont imprimés avec votre sortie. Lorsque vous avez un facteur à deux niveaux (par exemple, homme contre femme) et que vous utilisez le codage de niveau de référence, vous verrez l'interception appelée (constant)et une seule variable répertoriée dans la sortie (peut-être sexM). L'ordonnée à l'origine est la moyenne du groupe de référence (peut-être des femmes) et sexMla différence entre la moyenne des hommes et la moyenne des femmes. La valeur de p associée à l'ordonnée à l'origine est un test un échantillon pour déterminer si le niveau de référence a une moyenne de et la valeur de p associée àt0t 0sexMvous indique si les sexes diffèrent dans votre réponse. Mais si vous utilisez le codage de niveau signifie à la place, vous aurez deux variables répertoriées et chaque valeur p correspondra à un test à un échantillon pour savoir si la moyenne de ce niveau est . Autrement dit, aucune des valeurs de p ne sera un test pour savoir si les sexes diffèrent. t0

set.seed(1)
y    = c(    rnorm(30), rnorm(30, mean=1)         )
sex  = rep(c("Female",  "Male"          ), each=30)
fem  = ifelse(sex=="Female", 1, 0)
male = ifelse(sex=="Male", 1, 0)

ref.level.coding.model   = lm(y~sex)
level.means.coding.model = lm(y~fem+male+0)

summary(ref.level.coding.model)
# ...
# Coefficients:
#             Estimate Std. Error t value Pr(>|t|)    
# (Intercept)  0.08246    0.15740   0.524    0.602    
# sexMale      1.05032    0.22260   4.718 1.54e-05 ***
#   ---
#   Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
# ...
summary(level.means.coding.model)
# ...
# Coefficients:
#      Estimate Std. Error t value Pr(>|t|)    
# fem   0.08246    0.15740   0.524    0.602    
# male  1.13277    0.15740   7.197 1.37e-09 ***
#   ---
#   Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
# ...
gung - Réintégrer Monica
la source
1
Merci pour l'ajout du code: cela démontre clairement que l'interception dans le codage de la cellule de référence = sexFemale in cell signifie le codage (j'ai également été confondu par "ce qui est arrivé" à la variable N-1 ...) Peut-être a-t-il besoin d'une autre question, mais c'est facile à comprendre avec une variable, qu'en est-il de 2 ou plus? ex: âge: "vieux" "jeune". Dans le codage des cellules ref, les coefficients s'afficheraient-ils pour sexMale, ageYoung (par exemple) et le compte d'interception pour BOTH sexFemale et ageOld?
Bryan
1
Si vous avez un facteur avec 3 niveaux, l'ordonnée à l'origine est la moyenne du niveau de référence, et les 2 autres seront représentés dans la sortie. Leurs coefficients seront à la fois la différence entre eux et le niveau de référence et leur ps sera la signification de ces difs. Si vous avez 2 facteurs, les deux auront des niveaux de référence, et l'interception sera la moyenne de ces personnes qui sont dans les deux groupes de référence (par exemple, young F) et les autres niveaux seront dif b / t le niveau donné de facteur 1 w / le niveau de référence de l'autre facteur et le groupe des deux niveaux de référence. Par exemple oldest old F- «jeune F , & M» est young M- young F.
gung - Rétablir Monica
1
J'ai joué un peu avec cela et j'ai connu une R^2différence substantielle entre les deux approches. Je sais que ce n'est qu'un R^2, mais y a-t-il une explication à cela?
hans0l0
@ hans0l0, je n'en ai aucune idée. Il ne devrait pas y avoir de différence.
gung - Réintègre Monica
1
@confused, vous pouvez le trouver dans la documentation ,? glm .
gung - Rétablir Monica
2

Les coefficients estimés seraient les mêmes sous réserve que vous créiez vos variables fictives (c'est-à-dire les variables numériques) cohérentes avec R. Par exemple: permet de créer de fausses données et d'ajuster une glm de Poisson à l'aide d'un facteur. Notez que la glfonction crée une variable de facteur.

> counts <- c(18,17,15,20,10,20,25,13,12)
> outcome <- gl(3,1,9)
> outcome
[1] 1 2 3 1 2 3 1 2 3
Levels: 1 2 3
> class(outcome)
[1] "factor"
> glm.1<- glm(counts ~ outcome, family = poisson())
> summary(glm.1)

Call:
glm(formula = counts ~ outcome, family = poisson())

Deviance Residuals: 
    Min       1Q   Median       3Q      Max  
-0.9666  -0.6713  -0.1696   0.8471   1.0494  

Coefficients:
            Estimate Std. Error z value Pr(>|z|)    
(Intercept)   3.0445     0.1260  24.165   <2e-16 ***
outcome2     -0.4543     0.2022  -2.247   0.0246 *  
outcome3     -0.2930     0.1927  -1.520   0.1285    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1   1

(Dispersion parameter for poisson family taken to be 1)

    Null deviance: 10.5814  on 8  degrees of freedom
Residual deviance:  5.1291  on 6  degrees of freedom
AIC: 52.761

Number of Fisher Scoring iterations: 4

Étant donné que le résultat a trois niveaux, je crée deux variables muettes (dummy.1 = 0 si le résultat = 2 et dummy.2 = 1 si le résultat = 3) et le réajustent en utilisant ces valeurs numériques:

> dummy.1=rep(0,9)
> dummy.2=rep(0,9)
> dummy.1[outcome==2]=1
> dummy.2[outcome==3]=1
> glm.2<- glm(counts ~ dummy.1+dummy.2, family = poisson())
> summary(glm.2)

Call:
glm(formula = counts ~ dummy.1 + dummy.2, family = poisson())

Deviance Residuals: 
    Min       1Q   Median       3Q      Max  
-0.9666  -0.6713  -0.1696   0.8471   1.0494  

Coefficients:
            Estimate Std. Error z value Pr(>|z|)    
(Intercept)   3.0445     0.1260  24.165   <2e-16 ***
dummy.1      -0.4543     0.2022  -2.247   0.0246 *  
dummy.2      -0.2930     0.1927  -1.520   0.1285    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1   1

(Dispersion parameter for poisson family taken to be 1)

    Null deviance: 10.5814  on 8  degrees of freedom
Residual deviance:  5.1291  on 6  degrees of freedom
AIC: 52.761

Number of Fisher Scoring iterations: 4

Comme vous pouvez le voir, les coefficients estimés sont les mêmes. Mais vous devez être prudent lors de la création de vos variables factices si vous souhaitez obtenir le même résultat. Par exemple, si je crée deux variables muettes comme (dummy.1 = 0 si le résultat = 1 et dummy.2 = 1 si le résultat = 2), les résultats estimés sont différents comme suit:

> dummy.1=rep(0,9)
> dummy.2=rep(0,9)
> dummy.1[outcome==1]=1
> dummy.2[outcome==2]=1
> glm.3<- glm(counts ~ dummy.1+dummy.2, family = poisson())
> summary(glm.3)

Call:
glm(formula = counts ~ dummy.1 + dummy.2, family = poisson())

Deviance Residuals: 
    Min       1Q   Median       3Q      Max  
-0.9666  -0.6713  -0.1696   0.8471   1.0494  

Coefficients:
            Estimate Std. Error z value Pr(>|z|)    
(Intercept)   2.7515     0.1459   18.86   <2e-16 ***
dummy.1       0.2930     0.1927    1.52    0.128    
dummy.2      -0.1613     0.2151   -0.75    0.453    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1   1

(Dispersion parameter for poisson family taken to be 1)

    Null deviance: 10.5814  on 8  degrees of freedom
Residual deviance:  5.1291  on 6  degrees of freedom
AIC: 52.761

Number of Fisher Scoring iterations: 4

En effet, lorsque vous ajoutez une outcomevariable dans glm.1, R crée par défaut deux variables factices, à savoir outcome2et outcome3et les définit de manière similaire à dummy.1et dummy.2dans glm.2, c'est-à-dire que le premier niveau de résultat est lorsque toutes les autres variables factices ( outcome2et outcome3) sont définies pour être zéro.

Stat
la source
Merci pour la preuve de code des coefficients estimés étant les mêmes. L'avertissement sur la création de la mienne est également utile: je voulais créer la mienne car une variable de modèle serait liée directement à une colonne de base de données par son nom (ce qui pourrait être utile en aval), mais il semble que je doive comprendre les problèmes de Je vais faire ça.
Bryan