Zrozumienie tworzenia zmiennych fikcyjnych (ręcznych lub automatycznych) w GLM

Jeśli zmienna czynnikowa (np. Płeć o poziomach M i F) jest używana we wzorze glm, tworzone są zmienne zmienne, które można znaleźć w podsumowaniu modelu glm wraz z powiązanymi z nimi współczynnikami (np. Płeć)

Jeśli zamiast polegać na R w celu podzielenia współczynnika w ten sposób, czynnik jest zakodowany w szeregu liczbowych zmiennych 0/1 (np. GenderM (1 dla M, 0 dla F), genderF (1 dla F, 0 dla M) i te zmienne są następnie używane jako zmienne numeryczne we wzorze glm, czy wynik współczynnika byłby inny?

Zasadniczo pytanie brzmi: czy R stosuje inne obliczenia współczynnika podczas pracy ze zmiennymi czynnikowymi w porównaniu ze zmiennymi numerycznymi?

Pytanie uzupełniające (na które może odpowiedzieć powyższe pytanie): czy oprócz wydajności pozwalającej R tworzyć zmienne zastępcze, czy jest jakiś problem z ponownym kodowaniem czynników jako szeregu liczbowych zmiennych 0,1 i zastosowaniem tych w modelu?

Zmienne kategorialne (zwane „ czynnikami ” w R) muszą być reprezentowane przez kody numeryczne w modelach regresji wielokrotnej. Istnieje wiele możliwych sposobów prawidłowego konstruowania kodów numerycznych (zobacz tę wspaniałą listę na stronie pomocy statystyk UCLA). Domyślnie R używa kodowania na poziomie odniesienia (który R nazywa „contr.treatment”) i który jest w zasadzie domyślny dla całej statystyki. Można to zmienić dla wszystkich kontrastów dla całej sesji R za pomocą opcji lub dla określonych analiz / zmiennych przy użyciu kontrastów lub C (zwróć uwagę na kapitał). Jeśli potrzebujesz więcej informacji na temat kodowania na poziomie odniesienia, wyjaśnię to tutaj: Regresja oparta na przykład na dniach tygodnia.

Niektóre osoby uważają kodowanie poziomu odniesienia za mylące i nie trzeba go używać. Jeśli chcesz, możesz mieć dwie zmienne dla mężczyzn i kobiet; to się nazywa poziom oznacza kodowanie. Jednak jeśli to zrobisz, będziesz musiał powstrzymać przechwytywanie, w przeciwnym razie matryca modelu będzie pojedyncza, a regresja nie będzie pasować, jak powyższe nuty @Affine i jak wyjaśnię tutaj: Kodowanie zmiennych jakościowych prowadzi do osobliwości . Aby wyłączyć przechwytywanie, zmodyfikuj formułę, dodając -1lub w ten +0sposób: y~... -1lub y~... +0.

Zastosowanie kodowania poziomu zamiast kodowania poziomu odniesienia zmieni szacowane współczynniki i znaczenie testów hipotez wydrukowanych wraz z wynikami. Kiedy masz współczynnik dwupoziomowy (np. Męski vs. żeński) i używasz kodowania poziomu referencyjnego, zobaczysz wywołany intercept (constant)i tylko jedną zmienną wymienioną na wyjściu (być może sexM). Punkt przecięcia jest średnią grupy odniesienia (być może kobiet) i sexMstanowi różnicę między średnią mężczyzn a średnią kobiet. Wartość p powiązana z punktem przecięcia jest testem jednej próby na to, czy poziom odniesienia ma średnią a wartość p powiązana z$t$$0$sexMinformuje cię, czy płcie różnią się w twojej odpowiedzi. Ale jeśli zamiast tego użyjesz kodowania oznaczającego poziom, będziesz mieć na liście dwie zmienne, a każda wartość p będzie odpowiadać testowi próby z jedną próbą, czy średnia tego poziomu wynosi . Oznacza to, że żadna z wartości p nie będzie testem na to, czy płcie się różnią. $t$$0$

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
# ...

Oszacowane współczynniki byłyby takie same pod warunkiem, że utworzysz zmienne fikcyjne (tj. Numeryczne) zgodne z R. Na przykład: stwórzmy fałszywe dane i dopasuj współczynnik Poissona glm za pomocą współczynnika. Zauważ, że glfunkcja tworzy zmienną czynnikową.

> 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

Ponieważ wynik ma trzy poziomy, tworzę dwie zmienne fikcyjne (fikcyjne 1 = 0, jeśli wynik = 2 i fikcyjne 2 = 1, jeśli wynik = 3) i ponownie używam tych wartości liczbowych:

> 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

Jak widać szacunkowe współczynniki są takie same. Ale musisz zachować ostrożność podczas tworzenia zmiennych zastępczych, jeśli chcesz uzyskać ten sam wynik. Na przykład, jeśli utworzę dwie zmienne fikcyjne jako (fikcyjny 1 = 0, jeśli wynik = 1 i fikcyjny 2 = 1, jeśli wynik = 2), wówczas szacunkowe wyniki są różne w następujący sposób:

> 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

Wynika to z tego, że gdy dodajesz outcomezmienną w glm.1, R domyślnie tworzy dwie zmienne fikcyjne mianowicie outcome2i outcome3definiuje je podobnie do dummy.1i dummy.2w glm.2, tzn. Pierwszym poziomem wyniku jest sytuacja, gdy wszystkie pozostałe zmienne fikcyjne ( outcome2i outcome3) są ustawione na zero.