pqstat.pl

Narzędzia użytkownika

Modele wielowymiarowe

Przygotowanie zmiennych do analizy

Interakcje

Interakcje rozważane są w modelach wielowymiarowych a ich występowanie oznacza, że wpływ zmiennej niezależnej ($X_1$) na zmienną zależną ($Y$) jest inny, w zależności od poziomu kolejnej zmiennej niezależnej ($X_2$) lub szeregu kolejnych zmiennych niezależnych. By można było rozważać interekcje w modelach wielowymiarowych należy wskazać zmienne mówiące o prawdopodobnych interakcjach, czyli iloczyny odpowiednich zmiennych. W tym celu wybieramy przycisk Interakcje w oknie wybranej analizy wielowymiarowej. W oknie ustawiania interakcji z wciśniętym przyciskiem CTRL wskazujemy zmienne, które mają tworzyć interakcje i przenosimy je do sąsiedniej listy przy pomocy strzałki. Uruchamiając przycisk OK uzyskujemy odpowiednie kolumny w arkuszu danych.

W analizie interakcji wybór odpowiedniego kodowania zmiennych dychotomicznych pozwala na uniknięcie przeparametryzowania związanego z interakcjami. Przeparametryzowanie powoduje, że efekty niższego rzędu dla zmiennych dychotomicznych są redundantne względem uwikłanych interakcji wyższego rzędu. W rezultacie uwzględnienie w modelu interakcji wyższego rzędu niweluje efekt interakcji rzędów niższych, nie pozwalając na ich prawidłową ocenę. By uniknąć przeparametryzowania w modelu w którym występują interakcje zmiennych dychotomicznych zaleca się wybierać opcję kodowanie efektów.

W modelach z interakcjami należy pamiętać o odpowiednim ich „przycinaniu”, tak by usuwając efekty główne usunąć również efekty rzędów wyższych, które są od nich zależne. To znaczy: jeśli w modelu mamy następujące zmienne (efekty główne): $X_1$, $X_2$, $X_3$ i interakcje: $X_1*X_2$, $X_1*X_3$, $X_2*X_3$, $X_1*X_2*X_3$, to usuwając z modelu zmienną $X_1$ musimy usunąć również te interakcje, w których ona występuje, czyli: $X_1*X_2$, $X_1*X_3$ oraz $X_1*X_2*X_3$.

 

Kodowanie zmiennych

Problemem w przygotowaniu danych do analizy wielowymiarowej jest odpowiednie zakodowanie zmiennych nominalnych i porządkowych. Jest to ważny element przygotowania danych do analizy, gdyż ma zasadniczy wpływ na interpretację współczynników modelu. Zmienne nominalne lub porządkowe dzielą analizowane obiekty na dwie lub więcej kategorii, przy czym zmienne dychotomiczne (o dwóch kategoriach, $k=2$) wystarczy tylko odpowiednio zakodować, a zmienne o wielu kategoriach ($k>2$) rozbić na zmienne fikcyjne (ang. dummy variable) o dwóch kategoriach oraz zakodować.

  • [$k=2$] Jeśli zmienna jest dychotomiczna, badacz sam decyduje o sposobie, w jaki wprowadzi dane ją reprezentujące, może więc wprowadzić dowolne kody liczbowe np. 0 i 1. W programie można zmienić własne kodowanie na kodowanie efektu zaznaczając tę opcję w oknie wybranej analizy wielowymiarowej. Kodowanie takie powoduje zastąpienie mniejszej wartości wartością -1 a wartości większej wartością 1.
  • [$k>2$] Jeśli zmienna ma wiele kategorii, to w oknie wybranej analizy wielowymiarowej wybieramy przycisk Zmienne fikcyjne i ustawiamy kategorię referencyjną/bazową dla tych zmiennych, które chcemy rozbić na zmienne fikcyjne. Zmienne te będą zakodowane zero-jedynkowo, chyba, że w oknie analizy zostanie wybrana opcja kodowanie efektu - wówczas kodowane będą jako -1, 0 i 1.

Kodowanie zero-jedynkowe (dummy coding) jest wykorzystywane by przy pomocy modeli wielowymiarowych odpowiedzieć na pytanie: Jak wyniki ($Y$), w każdej analizowanej kategorii, różnią się od wyników kategorii referencyjnej. Kodowanie to polega na przypisaniu wartości 0 lub 1 do każdej kategorii danej zmiennej. Kategoria zakodowana jako 0 jest wówczas kategorią referencyjną (reference).

  • [$k=2$] Gdy kodowana zmienna jest dychotomiczna, wówczas umieszczając ją w modelu regresji uzyskamy wyliczony dla niej współczynnik ($b_i$). Współczynnik ten jest odniesieniem wartości zmiennej zależnej $Y$ dla kategorii 1 do kategorii referencyjnej (w korekcji o pozostałe zmienne w modelu).
  • [$k>2$] Gdy analizowana zmienna ma więcej niż dwie kategorie, wówczas $k$ kategorii jest reprezentowanych przez $k-1$ zmiennych fikcyjnych (dummy variables) o kodowaniu zero-jedynkowym. Tworząc zmienne o kodowaniu zero-jedynkowym wybiera się kategorię, dla której nie tworzy się zmiennej fikcyjnej. Kategoria ta traktowana jest w modelach jako kategoria odniesienia (gdyż w każdej zmiennej zakodowanej w sposób zero-jedynkowy odpowiadają jej wartości 0).

Gdy tak uzyskane zmienne $X_1, X_2, ..., X_{k-1}$ o kodowaniu zero-jedynkowym zostaną umieszczone w modelu regresji, wówczas zostaną dla nich wyliczone współczynniki $b_1, b_2, ..., b_{k-1}$.

  • [$b_1$] to odniesienie wyników $Y$ (dla kodów 1 w $X_1$) do kategorii referencyjnej (w korekcji o pozostałe zmienne w modelu);
  • [$b_2$] to odniesienie wyników $Y$ (dla kodów 1 w $X_2$) do kategorii referencyjnej (w korekcji o pozostałe zmienne w modelu)
  • […]
  • [$b_{k-1}$] to odniesienie wyników $Y$ (dla kodów 1 w $X_{k-1}$) do kategorii referencyjnej (w korekcji o pozostałe zmienne w modelu);

Przykład

Zakodujemy zgodnie z kodowaniem zero-jedynkowym zmienną płeć o dwóch kategoriach (płeć męską wybierzemy jako kategorię referencyjną) i zmienną wykształcenie o 4 kategoriach (wykształcenie podstawowe wybierzemy jako referencyjne).

\begin{tabular}{|c|c|}
\hline
&\textcolor[rgb]{0.5,0,0.5}{\textbf{Zakodowna}}\\
\textbf{Płeć}&\textcolor[rgb]{0.5,0,0.5}{\textbf{płeć}}\\\hline
k&\textcolor[rgb]{0.5,0,0.5}{1}\\
k&\textcolor[rgb]{0.5,0,0.5}{1}\\
k&\textcolor[rgb]{0.5,0,0.5}{1}\\
\cellcolor[rgb]{0.88,0.88,0.88}m&\cellcolor[rgb]{0.88,0.88,0.88}0\\
\cellcolor[rgb]{0.88,0.88,0.88}m&\cellcolor[rgb]{0.88,0.88,0.88}0\\
k&\textcolor[rgb]{0.5,0,0.5}{1}\\
k&\textcolor[rgb]{0.5,0,0.5}{1}\\
\cellcolor[rgb]{0.88,0.88,0.88}m&\cellcolor[rgb]{0.88,0.88,0.88}0\\
\cellcolor[rgb]{0.88,0.88,0.88}m&\cellcolor[rgb]{0.88,0.88,0.88}0\\
k&\textcolor[rgb]{0.5,0,0.5}{1}\\
\cellcolor[rgb]{0.88,0.88,0.88}m&\cellcolor[rgb]{0.88,0.88,0.88}0\\
k&\textcolor[rgb]{0.5,0,0.5}{1}\\
\cellcolor[rgb]{0.88,0.88,0.88}m&\cellcolor[rgb]{0.88,0.88,0.88}0\\
k&\textcolor[rgb]{0.5,0,0.5}{1}\\
\cellcolor[rgb]{0.88,0.88,0.88}m&\cellcolor[rgb]{0.88,0.88,0.88}0\\
\cellcolor[rgb]{0.88,0.88,0.88}m&\cellcolor[rgb]{0.88,0.88,0.88}0\\
...&...\\\hline
\end{tabular}

\begin{tabular}{|c|ccc|}
\hline
& \multicolumn{3}{c|}{\textbf{Zakodowane wykształcenie}}\\
\textbf{Wykształcenie}&\textcolor[rgb]{0,0,1}{\textbf{zawodowe}}&\textcolor[rgb]{1,0,0}{\textbf{średnie}}&\textcolor[rgb]{0,0.58,0}{\textbf{wyższe}}\\\hline
\cellcolor[rgb]{0.88,0.88,0.88}podstawowe&\cellcolor[rgb]{0.88,0.88,0.88}0&\cellcolor[rgb]{0.88,0.88,0.88}0&\cellcolor[rgb]{0.88,0.88,0.88}0\\
\cellcolor[rgb]{0.88,0.88,0.88}podstawowe&\cellcolor[rgb]{0.88,0.88,0.88}0&\cellcolor[rgb]{0.88,0.88,0.88}0&\cellcolor[rgb]{0.88,0.88,0.88}0\\
\cellcolor[rgb]{0.88,0.88,0.88}podstawowe&\cellcolor[rgb]{0.88,0.88,0.88}0&\cellcolor[rgb]{0.88,0.88,0.88}0&\cellcolor[rgb]{0.88,0.88,0.88}0\\
\textcolor[rgb]{0,0,1}{zawodowe}&\textcolor[rgb]{0,0,1}{1}&0&0\\
\textcolor[rgb]{0,0,1}{zawodowe}&\textcolor[rgb]{0,0,1}{1}&0&0\\
\textcolor[rgb]{0,0,1}{zawodowe}&\textcolor[rgb]{0,0,1}{1}&0&0\\
\textcolor[rgb]{0,0,1}{zawodowe}&\textcolor[rgb]{0,0,1}{1}&0&0\\
\textcolor[rgb]{1,0,0}{średnie}&0&\textcolor[rgb]{1,0,0}{1}&0\\
\textcolor[rgb]{1,0,0}{średnie}&0&\textcolor[rgb]{1,0,0}{1}&0\\
\textcolor[rgb]{1,0,0}{średnie}&0&\textcolor[rgb]{1,0,0}{1}&0\\
\textcolor[rgb]{1,0,0}{średnie}&0&\textcolor[rgb]{1,0,0}{1}&0\\
\textcolor[rgb]{0,0.58,0}{wyższe}&0&0&\textcolor[rgb]{0,0.58,0}{1}\\
\textcolor[rgb]{0,0.58,0}{wyższe}&0&0&\textcolor[rgb]{0,0.58,0}{1}\\
\textcolor[rgb]{0,0.58,0}{wyższe}&0&0&\textcolor[rgb]{0,0.58,0}{1}\\\
\textcolor[rgb]{0,0.58,0}{wyższe}&0&0&\textcolor[rgb]{0,0.58,0}{1}\\
\textcolor[rgb]{0,0.58,0}{wyższe}&0&0&\textcolor[rgb]{0,0.58,0}{1}\\
...&...&...&...\\\hline
\end{tabular}

Budując na podstawie zmiennych fikcyjnych, w modelu regresji wielorakiej, moglibyśmy chcieć sprawdzić jak zmienne te wpływają na pewną zmienną zależną np. $Y$ = wysokość zarobków (wyrażoną w tysiącach złotych). W wyniku takiej analizy dla każdej zmiennej fikcyjnej uzyskamy przykładowe współczynniki:

- dla płci istotny statystycznie współczynnik $b_{i}=-0.5$ - co oznacza, że średnie zarobki kobiet są o pół tysiąca złoty niższe niż mężczyzn; przy założeniu że pozostałe zmienne w modelu pozostają na stałym poziomie;

- dla wykształcenia zawodowego istotny statystycznie współczynnik $b_{i}=0.6$ - co oznacza, że średnie zarobki osób z wykształceniem zawodowym są o 0.6 tysiąca złoty wyższe niż dla osób z wykształceniem podstawowym; przy założeniu że pozostałe zmienne w modelu pozostają na stałym poziomie;

- dla wykształcenia średniego istotny statystycznie współczynnik $b_{i}=1$ - oznacza, że średnie zarobki osób z wykształceniem średnim są o tysiąc złoty wyższe niż dla osób z wykształceniem podstawowym; przy założeniu że pozostałe zmienne w modelu pozostają na stałym poziomie;

- dla wykształcenia wyższego istotny statystycznie współczynnik $b_{i}=1.5$ - co oznacza, że średnie zarobki osób z wykształceniem wyższym są o 1.5 tysiąca wyższe niż dla osób z wykształceniem podstawowym; przy założeniu że pozostałe zmienne w modelu pozostają na stałym poziomie.

Kodowanie efektów (effect coding) jest wykorzystywane, by przy pomocy modeli wielowymiarowych odpowiedzieć na pytanie: Jak wyniki ($Y$), w każdej analizowanej kategorii, różnią się od wyników średniej (nieważonej) uzyskanej z próby. Kodowanie to polega na przypisaniu wartości -1 lub 1 do każdej kategorii danej zmiennej. Kategoria zakodowana jako -1 jest wówczas kategorią bazową (base)

  • [$k=2$] Gdy kodowana zmienna jest dychotomiczna, wówczas umieszczając ją w modelu regresji uzyskamy wyliczony dla niej współczynnik ($b_i$). Współczynnik ten jest odniesieniem $Y$ dla kategorii 1 do nieważonej średniej ogólnej (w korekcji o pozostałe zmienne w modelu).

Gdy analizowana zmienna ma więcej niż dwie kategorie, wówczas $k$ kategorii jest reprezentowanych przez $k-1$ zmiennych fikcyjnych o kodowaniu efektu. Tworząc zmienne o kodowaniu efektu wybiera się kategorię dla której nie tworzy się oddzielnej zmiennej. Kategoria ta traktowana jest w modelach jako kategoria bazowa (gdyż w każdej zmiennej zapisanej poprzez kodowanie efektu odpowiadają jej wartości -1).

Gdy tak uzyskane zmienne $X_1, X_2, ..., X_{k-1}$ o kodowaniu efektu zostaną umieszczone w modelu regresji, wówczas zostaną dla nich wyliczone współczynniki $b_1, b_2, ..., b_{k-1}$.

  • [$b_1$] to odniesienie wyników $Y$ (dla kodów 1 w $X_1$) do nieważonej średniej ogólnej (w korekcji o pozostałe zmienne w modelu);
  • [$b_2$] to odniesienie wyników $Y$ (dla kodów 1 w $X_2$) do nieważonej średniej ogólnej (w korekcji o pozostałe zmienne w modelu);
  • […]
  • [$b_{k-1}$] to odniesienie wyników $Y$ (dla kodów 1 w $X_{k-1}$) do nieważonej średniej ogólnej (w korekcji o pozostałe zmienne w modelu);

Przykład

Zakodujemy przy pomocy kodowania efektu zmienną płeć o dwóch kategoriach (płeć męską wybierzemy jako kategorię bazową) i zmienną wskazującą region zamieszkania na terenie analizowanego kraju. Wyróżniono 5 regionów: północny, południowy, wschodni, zachodni i centralny - region centralny wybierzemy jako bazowy.

\begin{tabular}{|c|c|}
\hline
&\textcolor[rgb]{0.5,0,0.5}{\textbf{Zakodowona}}\\
\textbf{Płeć}&\textcolor[rgb]{0.5,0,0.5}{\textbf{płeć}}\\\hline
k&\textcolor[rgb]{0.5,0,0.5}{1}\\
k&\textcolor[rgb]{0.5,0,0.5}{1}\\
k&\textcolor[rgb]{0.5,0,0.5}{1}\\
\cellcolor[rgb]{0.88,0.88,0.88}m&\cellcolor[rgb]{0.88,0.88,0.88}-1\\
\cellcolor[rgb]{0.88,0.88,0.88}m&\cellcolor[rgb]{0.88,0.88,0.88}-1\\
k&\textcolor[rgb]{0.5,0,0.5}{1}\\
k&\textcolor[rgb]{0.5,0,0.5}{1}\\
\cellcolor[rgb]{0.88,0.88,0.88}m&\cellcolor[rgb]{0.88,0.88,0.88}-1\\
\cellcolor[rgb]{0.88,0.88,0.88}m&\cellcolor[rgb]{0.88,0.88,0.88}-1\\
k&\textcolor[rgb]{0.5,0,0.5}{1}\\
\cellcolor[rgb]{0.88,0.88,0.88}m&\cellcolor[rgb]{0.88,0.88,0.88}-1\\
k&\textcolor[rgb]{0.5,0,0.5}{1}\\
\cellcolor[rgb]{0.88,0.88,0.88}m&\cellcolor[rgb]{0.88,0.88,0.88}-1\\
k&\textcolor[rgb]{0.5,0,0.5}{1}\\
\cellcolor[rgb]{0.88,0.88,0.88}m&\cellcolor[rgb]{0.88,0.88,0.88}-1\\
\cellcolor[rgb]{0.88,0.88,0.88}m&\cellcolor[rgb]{0.88,0.88,0.88}-1\\
...&...\\\hline
\end{tabular}

\begin{tabular}{|c|cccc|}
\hline
\textbf{Region}& \multicolumn{4}{c|}{\textbf{Zakodowany region}}\\
\textbf{zamieszkania}&\textcolor[rgb]{0,0,1}{\textbf{zachodni}}&\textcolor[rgb]{1,0,0}{\textbf{wschodni}}&\textcolor[rgb]{0,0.58,0}{\textbf{północny}}&\textcolor[rgb]{0.55,0,0}{\textbf{południowy}}\\\hline
\cellcolor[rgb]{0.88,0.88,0.88}centralny&\cellcolor[rgb]{0.88,0.88,0.88}-1&\cellcolor[rgb]{0.88,0.88,0.88}-1&\cellcolor[rgb]{0.88,0.88,0.88}-1&\cellcolor[rgb]{0.88,0.88,0.88}-1\\
\cellcolor[rgb]{0.88,0.88,0.88}centralny&\cellcolor[rgb]{0.88,0.88,0.88}-1&\cellcolor[rgb]{0.88,0.88,0.88}-1&\cellcolor[rgb]{0.88,0.88,0.88}-1&\cellcolor[rgb]{0.88,0.88,0.88}-1\\
\cellcolor[rgb]{0.88,0.88,0.88}centralny&\cellcolor[rgb]{0.88,0.88,0.88}-1&\cellcolor[rgb]{0.88,0.88,0.88}-1&\cellcolor[rgb]{0.88,0.88,0.88}-1&\cellcolor[rgb]{0.88,0.88,0.88}-1\\
\textcolor[rgb]{0,0,1}{zachodni}&\textcolor[rgb]{0,0,1}{1}&0&0&0\\
\textcolor[rgb]{0,0,1}{zachodni}&\textcolor[rgb]{0,0,1}{1}&0&0&0\\
\textcolor[rgb]{0,0,1}{zachodni}&\textcolor[rgb]{0,0,1}{1}&0&0&0\\
\textcolor[rgb]{0,0,1}{zachodni}&\textcolor[rgb]{0,0,1}{1}&0&0&0\\
\textcolor[rgb]{1,0,0}{wschodni}&0&\textcolor[rgb]{1,0,0}{1}&0&0\\
\textcolor[rgb]{1,0,0}{wschodni}&0&\textcolor[rgb]{1,0,0}{1}&0&0\\
\textcolor[rgb]{1,0,0}{wschodni}&0&\textcolor[rgb]{1,0,0}{1}&0&0\\
\textcolor[rgb]{1,0,0}{wschodni}&0&\textcolor[rgb]{1,0,0}{1}&0&0\\
\textcolor[rgb]{0,0.58,0}{północny}&0&0&\textcolor[rgb]{0,0.58,0}{1}&0\\
\textcolor[rgb]{0,0.58,0}{północny}&0&0&\textcolor[rgb]{0,0.58,0}{1}&0\\
\textcolor[rgb]{0.55,0,0}{południowy}&0&0&0&\textcolor[rgb]{0.55,0,0}{1}\\
\textcolor[rgb]{0.55,0,0}{południowy}&0&0&0&\textcolor[rgb]{0.55,0,0}{1}\\
\textcolor[rgb]{0.55,0,0}{południowy}&0&0&0&\textcolor[rgb]{0.55,0,0}{1}\\
...&...&...&...&...\\\hline
\end{tabular}

Budując na podstawie zmiennych fikcyjnych, w modelu regresji wielorakiej, moglibyśmy chcieć sprawdzić jak zmienne te wpływają na pewną zmienną zależną np. $Y$ = wysokość zarobków (wyrażoną w tysiącach złotych). W wyniku takiej analizy dla każdej zmiennej fikcyjnej uzyskamy przykładowe współczynnik:

- dla płci istotny statystycznie współczynnik $b_{i}=-0.5$ - co oznacza, że średnie zarobki kobiet są o pół tysiąca złoty niższe niż średnie zarobki w kraju; przy założeniu że pozostałe zmienne w modelu pozostają na stałym poziomie;

- dla regionu zachodniego istotny statystycznie współczynnik $b_{i}=0.6$ - co oznacza, że średnie zarobki osób zamieszkujących na zachodzie kraju są o 0.6 tysiąca złoty wyższe niż średnie zarobki w kraju; przy założeniu że pozostałe zmienne w modelu pozostają na stałym poziomie;

- dla regionu wschodniego istotny statystycznie współczynnik $b_{i}=-1$ - oznacza, że średnie zarobki osób zamieszkujących na wschodzie kraju są o tysiąc złoty niższe niż średnie zarobki w kraju; przy założeniu że pozostałe zmienne w modelu pozostają na stałym poziomie;

- dla regionu północnego istotny statystycznie współczynnik $b_{i}=0.4$ - co oznacza, że średnie zarobki osób zamieszkujących na północy są o 0.4 tysiąca wyższe niż średnie zarobki w kraju; przy założeniu że pozostałe zmienne w modelu pozostają na stałym poziomie;

- dla regionu południowego nieistotny statystycznie współczynnik $b_{i}=0.1$ - co oznacza, że średnie zarobki osób zamieszkujących na południu nie różnią się istotnie od średnich zarobków w kraju; przy założeniu że pozostałe zmienne w modelu pozostają na stałym poziomie.

 
 

Liniowa regresja wieloraka

Okno z ustawieniami opcji Regresji wielorakiej wywołujemy poprzez menu StatystykaModele wielowymiaroweRegresja wieloraka

Budowany model regresji wielorakiej pozwala na zbadanie wpływu wielu zmiennych niezależnych ($X_1$, $X_2$, $\ldots$, $X_k$) na jedną zmienną zależną ($Y$). Najczęściej wykorzystywaną odmianą regresji wielorakiej jest Liniowa Regresja Wieloraka. Jest ona rozszerzeniem modeli regresji liniowej opartej o współczynnik korelacji liniowej Pearsona. Zakłada ona występowanie liniowego związku pomiędzy badanymi zmiennymi. Liniowy model regresji wielorakiej przyjmuje postać: \begin{displaymath}
Y=\beta_0+\beta_1X_1+\beta_2X_2+\ldots+\beta_kX_k+\epsilon.
\end{displaymath}

gdzie:

$Y$ - zmienna zależna, objaśniana przez model,

$X_1,X_2,\ldots X_k$ - zmienne niezależne, objaśniające,

$\beta_0,\beta_1,\beta_2,\ldots \beta_k$ - parametry,

$\epsilon$ - składnik losowy (reszta modelu).

Jeśli model został stworzony w oparciu o próbę o liczności $n$ powyższe równanie można przedstawić w postaci macierzowej: \begin{displaymath}
Y=X\beta+\epsilon.
\end{displaymath}

gdzie: 
$
Y=\left( \begin{array}{ccc}
y_1\\
y_2\\
\vdots\\
y_n
\end{array}\right),
$

$
X=\left( \begin{array}{ccccc}
1 & x_{11} & x_{21} & \ldots & x_{k1}\\
1 & x_{12} & x_{22} & \ldots & x_{k2}\\
\vdots & \vdots & \vdots & \ddots & \vdots\\
1 & x_{1n} & x_{2n} & \ldots & x_{kn}
\end{array} \right),
$

$
\beta=\left( \begin{array}{ccc}
\beta_0\\
\beta_1\\
\beta_2\\
\vdots\\
\beta_k
\end{array} \right),
$

$
\epsilon=\left( \begin{array}{ccc}
\epsilon_1\\
\epsilon_2\\
\vdots\\
\epsilon_n
\end{array}\right).
$

Rozwiązaniem równania jest wówczas wektor ocen parametrów $\beta_0,\beta_1,\ldots,\beta_k$ nazywanych współczynnikami regresji:


$
b=\left( \begin{array}{ccc}
b_0\\
b_1\\
b_2\\
\vdots\\
b_k
\end{array}\right). $

Współczynniki te szacowane są poprzez klasyczną metodę najmniejszych kwadratów. Na podstawie tych wartości możemy wnioskować o wielkości wpływu zmiennej niezależnej (dla której ten współczynnik został oszacowany) na zmienną zależną. Podają o ile jednostek zmieni się zmienna zależna, gdy zmienną niezależną zmienimy o 1 jednostkę. Każdy współczynnik obarczony jest pewnym błędem szacunku. Wielkość tego błędu wyliczana jest ze wzoru:

\begin{displaymath}
SE_b=\sqrt{\frac{1}{n-(k+1)}e^Te(X^TX)^{-1}},
\end{displaymath}

gdzie:

$e=Y-\widehat{Y}$ to wektor reszt modelu (różnica pomiędzy rzeczywistymi wartościami zmiennej zależnej Y a wartościami $\widehat{Y}$ przewidywanymi na podstawie modelu).

Zmienne fikcyjne i interakcje w modelu

Omówienie przygotowania zmiennych fikcyjnych i interakcji przedstawiono w rozdziale Przygotowanie zmiennych do analizy w modelach wielowymiarowych.

Uwaga! Budując model należy pamiętać, że liczba obserwacji musi być większa lub równa liczbie szacowanych parametrów modelu ($n\ge k+1$).

Weryfikacja modelu

  • Istotność statystyczna poszczególnych zmiennych w modelu.

Na podstawie współczynnika oraz jego błędu szacunku możemy wnioskować czy zmienna niezależna, dla której ten współczynnik został oszacowany wywiera istotny wpływ na zmienną zależną. W tym celu posługujemy się testem t-Studenta.

Hipotezy:

\begin{array}{cc}
\mathcal{H}_0: & \beta_i=0,\\
\mathcal{H}_1: & \beta_i\ne 0.
\end{array}

Wyliczmy statystykę testową według wzoru:

\begin{displaymath}
t=\frac{b_i}{SE_{b_i}}
\end{displaymath}

Statystyka testowa ma rozkład t-Studenta z $n-k$ stopniami swobody.

Wyznaczoną na podstawie statystyki testowej wartość $p$ porównujemy z poziomem istotności $\alpha$:

\begin{array}{ccl}
$ jeżeli $ p \le \alpha & \Longrightarrow & $ odrzucamy $ \mathcal{H}_0 $ przyjmując $ 	\mathcal{H}_1, \\
$ jeżeli $ p > \alpha & \Longrightarrow & $ nie ma podstaw, aby odrzucić $ \mathcal{H}_0. \\
\end{array}

  • Jakość zbudowanego modelu liniowej regresji wielorakiej możemy ocenić kilkoma miarami.
  • Błąd standardowy estymacji - jest miarą dopasowania modelu:

\begin{displaymath}
SE_e=\sqrt{\frac{\sum_{i=1}^ne_i^2}{n-(k+1)}}.
\end{displaymath}

Miara ta opiera się na resztach modelu $e_i=y_i-\widehat{y}_i$, czyli rozbieżności pomiędzy rzeczywistymi wartościami zmiennej zależnej $y_i$ w próbie a wartościami zmiennej zależnej $\widehat{y}_i$ wyliczonej na podstawie zbudowanego modelu. Najlepiej byłoby, gdyby różnica ta była jak najbliższa zeru dla wszystkich badanych obiektów próby. Zatem, aby model był dobrze dopasowany, błąd standardowy estymacji ($SE_e$) wyrażony jako wariancja $e_i$, powinien być jak najmniejszy.

  • Współczynnik korelacji wielorakiej $R=\sqrt{R^2} \in <0; 1>$ - określa siłę oddziaływania zespołu zmiennych $X_1,X_2,\ldots X_k$ na zmienną zależną $Y$.
  • Współczynnik determinacji wielorakiej $R^2$ - jest miarą dopasowania modelu.

Wartość tego współczynnika mieści się w przedziale $<0; 1>$, gdzie 1 oznacza doskonałe dopasowanie modelu, 0 - zupełny bark dopasowania. W jego wyznaczeniu posługujemy się następującą równością:

\begin{displaymath}
T_{SS}=E_{SS}+R_{SS},
\end{displaymath}

gdzie:

$T_{SS}$ - całkowita suma kwadratów,

$E_{SS}$ - suma kwadratów wyjaśniona przez model,

$R_{SS}$ - resztowa suma kwadratów.

Współczynnik determinacji wyliczamy z wzoru:

\begin{displaymath}
R^2=\frac{T_{SS}}{E_{SS}}.
\end{displaymath}

Wyraża on procent zmienności zmiennej zależnej tłumaczony przez model.

Ponieważ wartość współczynnika $R^2$ zależy od dopasowania modelu, ale jest również wrażliwa na ilość zmiennych w modelu i liczność próby, bywają sytuacje, w których może być obarczona pewnym błędem. Dalego też wyznacza się poprawianą wartość tego parametru:

\begin{displaymath}
R^2_{adj}=R^2-\frac{k(1-R^2)}{n-(k+1)}.
\end{displaymath}

  • Istotność statystyczna wszystkich zmiennych w modelu

Podstawowym narzędziem szacującym istotność wszystkich zmiennych w modelu jest test analizy wariancji (test F). Test ten weryfikuje jednocześnie 3 równoważne hipotezy:

\begin{array}{cc}
\mathcal{H}_0: & \textrm{wszystkie } \beta_i=0,\\
\mathcal{H}_0: & R^2=0,\\
\mathcal{H}_0: & $liniowość związku$,
\end{array} \begin{array}{cc}
\mathcal{H}_1: & \textrm{istnieje }\beta_i\neq0;\\
\mathcal{H}_1: & R^2\neq0;\\
\mathcal{H}_1: & $brak związku liniowego$.
\end{array}

Statystyka testowa ma postać:

\begin{displaymath}
F=\frac{E_{MS}}{R_{MS}}
\end{displaymath}

gdzie:

$\displaystyle E_{MS}=\frac{E_{SS}}{df_{E}}$ - średnia kwadratów wyjaśniona przez model,

$\displaystyle R_{MS}=\frac{R_{SS}}{df_{R}}$ - resztowa średnia kwadratów,

$df_E=k$, $df_R=n-(k+1)$ - odpowiednie stopnie swobody.

Statystyka ta podlega rozkładowi F Snedecora z $df_E$ i $df_R$ stopniami swobody.

Wyznaczoną na podstawie statystyki testowej wartość $p$ porównujemy z poziomem istotności $\alpha$:

\begin{array}{ccl}
$ jeżeli $ p \le \alpha & \Longrightarrow & $ odrzucamy $ \mathcal{H}_0 $ przyjmując $ 	\mathcal{H}_1, \\
$ jeżeli $ p > \alpha & \Longrightarrow & $ nie ma podstaw, aby odrzucić $ \mathcal{H}_0. \\
\end{array}

Przykład c.d. (plik wydawca.pqs)

 

Więcej informacji o zmiennych w modelu

  • Standaryzowane $b_1,b_2,\ldots,b_k$ - w odróżnieniu od parametrów surowych (które w zależności od opisywanej zmiennej są wyrażone w różnych jednostkach miary i nie mogą być bezpośrednio porównywane) standaryzowane oceny parametrów modelu pozwalają porównywać wkład poszczególnych zmiennych w wyjaśnienie zmienności zmiennej zależnej $Y$.
  • Macierz korelacji - zawiera informacje o sile związku pomiędzy poszczególnymi zmiennymi, czyli współczynnik korelacji Pearsona $r_p \in <-1; 1>$. Wspólczynnikiem tym badamy korelację dla każdej pary zmiennych, nie uwzględniając wpływu pozostałych zmiennych w modelu.
  • Macierz kowariancji - podobnie jak macierz korelacji, zawiera informacje o związku liniowym pomiędzy poszczególnymi zmiennymi. Przy czym wartość ta nie jest wystandaryzowana.
  • Współczynnik korelacji cząstkowej - należy do przedziału $<-1; 1>$ i jest miarą korelacji pomiędzy konkretną zmienną niezależną $X_i$ (uwzględniając jej skorelowanie z pozostałymi zmiennymi w modelu) a zmienną zależną $Y$ (uwzględniając jej skorelowanie z pozostałymi zmiennymi w modelu).

Kwadrat tego współczynnika to współczynnik determinacji cząstkowej - należy do przedziału $<0; 1>$ i oznacza stosunek wyłącznej zmienności danej zmiennej niezależnej $X_i$ do tej zmienności zmiennej zależnej $Y$, która nie została wyjaśniona przez pozostałe zmienne w modelu.

Im wartość tych współczynników znajduje się bliżej 0, tym bardziej bezużyteczną informację niesie badana zmienna, czyli jest ona nadmiarowa.

  • Współczynnik korelacji semicząstkowej - należy do przedziału $<-1; 1>$ i jest miarą korelacji pomiędzy konkretną zmienną niezależną $X_i$ (uwzględniając jej skorelowanie z pozostałymi zmiennymi w modelu) a zmienną zależną $Y$ (NIE uwzględniając jej skorelowanie z pozostałymi zmiennymi w modelu).

Kwadrat tego współczynnika to współczynnik determinacji semicząstkowej - należy do przedziału $<0; 1>$ i oznacza stosunek wyłącznej zmienności danej zmiennej niezależnej $X_i$ do całkowitej zmienności zmiennej zależnej $Y$.

Im wartość tych współczynników znajduje się bliżej zera, tym bardziej bezużyteczną informację niesie badana zmienna, czyli jest ona nadmiarowa.

  • R-kwadrat ($R^2 \in <0; 1>$) - wyraża on procent zmienności danej zmiennej niezależnej $X_i$ tłumaczony przez pozostałe zmienne niezależne. Im bliżej wartości 1 tym silniej badana zmienna związana jest liniowo z pozostałymi zmiennymi niezależnymi, co może oznaczać, że jest ona zmienną nadmiarową.
  • Tolerancja = $1-R^2 \in<0; 1>$- wyraża on procent zmienności danej zmiennej niezależnej $X_i$ NIE tłumaczony przez pozostałe zmienne niezależne. Im wartość tolerancji jest bliższa 0 tym silniej badana zmienna związana jest liniowo z pozostałymi zmiennymi niezależnymi, co może oznaczać, że jest ona zmienną nadmiarową.
  • Porównanie modelu pełnego z modelem po usunięciu danej zmiennej

Porównanie tych dwóch modeli dokonujemy:

  • testem F, w sytuacji gdy z modelu usuwamy jedną zmienną lub wiecej niż jedną zmienną (patrz porównywanie modeli),
  • testem t-Studenta, gdy z modelu usuwamy tylko jedną zmienną. Jest to ten sam test, którym badamy istotność poszczególnych zmiennych w modelu.

W przypadku usunięcia tylko jednej zmiennej wyniki obu tych testów są tożsame.

Jeśli różnica pomiędzy porównywanymi modelami jest istotna statystycznie (wartość $p \le \alpha$), wówczas model pełny jest istotnie lepszy niż model zredukowany. To oznacza, że badana zmienna nie jest nadmiarowa, wywiera ona istotny wpływ na dany model i nie powinna być z niego usuwana.

  • Wykresy rozrzutu

Wykresy te pozwalają dokonać subiektywnej oceny liniowości związku pomiędzy zmiennymi i zidentyfikować punkty odstające. Dodatkowo wykresami rozrzutu możemy posłużyć się w analizie reszt modelu.

 

Analiza reszt modelu

By otrzymać poprawny model regresji, powinniśmy sprawdzić podstawowe założenia dotyczące reszt modelu.

  • Obserwacje odstające

Badając reszty modelu szybko można uzyskać wiedzę na temat wartości odstających. Obserwacje takie mogą bardzo zaburzyć równanie regresji, ponieważ mają duży wpływ na wartości współczynników tego równania. Jeśli dana reszta $e_i$ jest oddalona o więcej niż 3 odchylenia standardowe od wartości średniej, wówczas obserwacje taką można uznać za obserwacje odstającą. Usunięcie obserwacji odstającej może w znaczącym stopniu przyczynić się do poprawy modelu.

  • Normalność rozkładu reszt modelu

Założenie to sprawdzamy przy pomocy testu normalności Lillieforsa. Duża różnica między rozkładem reszt a rozkładem normalnym (wartość $p \le \alpha$) może zaburzać ocenę istotności współczynników poszczególnych zmiennych modelu.

  • Homoskedastyczność (stałość wariancji)

By sprawdzić czy istnieją obszary, gdzie wariancja reszt modelu jest zwiększona lub zmniejszona posługujemy się wykresami:

  • reszty względem wartości przewidywanych
  • kwadrat reszty względem wartości przewidywanych
  • reszty względem wartości obserwowanych
  • kwadrat reszty względem wartości obserwowanych
  • Autokorelacja reszt modelu

Aby zbudowany model można było uznać za poprawny, wartości reszt nie powinny być ze sobą skorelowane (dla wszystkich par $e_i, e_j$). Założenie to możemy sprawdzić wyliczając statystykę testu Durbina-Watsona

\begin{displaymath}
d=\frac{\sum_{t=2}^n\left(e_t-e_{t-1}\right)^2}{\sum_{t=1}^ne_t^2},
\end{displaymath}

Aby sprawdzić dodatnią autokorelację na poziomie istotności $\alpha$, sprawdzamy położenie statystyki $d$ w stosunku do górnej ($d_{U,\alpha}$) i dolnej ($d_{L,\alpha}$) wartości krytycznej:

  • Jeżeli $d <d_{L,\alpha}$ - błędy są dodatnio skorelowane;
  • Jeśli $d> d_{U,\alpha}$ - błędy nie są dodatnio skorelowane;
  • Jeśli $d_{L,\alpha}<d <d_{U,\alpha}$ - wynik testu jest niejednoznaczny.

Aby sprawdzić ujemną autokorelację na poziomie istotności $\alpha$, sprawdzamy położenie wartości $4-d$ w stosunku do górnej ($d_{U,\alpha}$) i dolnej ($d_{L,\alpha}$) wartości krytycznej:

  • Jeżeli $4-d <d_{L,\alpha}$ - błędy są ujemnie skorelowane;
  • Jeśli $4-d> d_{U,\alpha}$ - błędy nie są ujemnie skorelowane;
  • Jeśli $d_{L,\alpha}<4-d <d_{U,\alpha}$ - wynik testu jest niejednoznaczny.

Wartości krytyczne testu Durbina-Watsona dla poziomu istotności $\alpha=0.05$ znajdują się na stronie internetowej (pqstat) - źródło tablic: Savina i White (1977)1)

Przykład c.d. (plik wydawca.pqs)

 

Predykcja na podstawie modelu

Najczęściej ostatnim etapem analizy regresji jest wykorzystanie zbudowanego i uprzednio zweryfikowanego modelu do predykcji. Przewidywanie wartości zmiennej zależnej jest możliwe dla zadanych wartości zmiennych niezależnych. Oszacowana wartość wyliczana jest z pewnym błędem. Dlatego też dodatkowo dla wyliczonej wartości wyznaczane są granice wynikające z błędu:

  • dla wartości oczekiwanej wyznaczane są granice ufności,
  • dla pojedynczego punktu wyznaczane są granice predykcji.

Przykład (plik wydawca.pqs)

Pewien wydawca książek chciał się dowiedzieć, jaki wpływ na zysk brutto ze sprzedaży mają takie zmienne jak: koszty produkcji, koszty reklamy, koszty promocji bezpośredniej, suma udzielonych rabatów, popularność autora. W tym celu przeanalizował 40 pozycji wydanych w ciągu ostatniego roku. Fragment danych przedstawia poniższy rysunek:

Pięć pierwszych zmiennych wyrażonych jest w tysiącach dolarów - są to więc zmienne zebrane na skali interwałowej. Natomiast ostatnia zmienna: popularność autora $-$ to zmienna dychotomiczna, gdzie 1 oznacza autora znanego, 0 oznacza autora nieznanego.

Na podstawie uzyskanej wiedzy wydawca planuje przewidzieć zysk brutto z kolejnej wydawanej książki znanego autora. Koszty, jakie zamierza ponieść to: koszty produkcji $\approx 11$, koszty reklamy $\approx 13$, koszty promocji bezpośredniej $\approx 0.5$, suma udzielonych rabatów $\approx 0.5$.

Budujemy model liniowej regresji wielorakiej wybierając: zysk brutto $-$ jako zmienną zależną $Y$, koszty produkcji, koszty reklamy, koszty promocji bezpośredniej, suma udzielonych rabatów, popularność autora $-$ jako zmienne niezależne $X_1, X_2, X_3, X_4, X_5$. W rezultacie wyliczone zostaną współczynniki równania regresji oraz miary pozwalające ocenić jakość modelu.

Na podstawie oszacowanej wartości współczynnika $b$, związek pomiędzy zyskiem brutto a wszystkimi zmiennymi niezależnymi możemy opisać równaniem: \begin{displaymath}
zysk_{brutto}=4.18+2.56(k_{prod})+2(k_{rekl})+4.67(k_{prom})+1.42(rabaty)+10.15(popul_{autora})+[8.09]
\end{displaymath} Uzyskane współczynniki interpretujemy następująco:

  • Jeśli koszt produkcji wzrośnie o 1 tysiąc dolarów, to zysk brutto wzrośnie o około 2.56 tysiące dolarów, przy złożeniu, że pozostałe zmienne się nie zmienią;
  • Jeśli koszt reklamy wzrośnie o 1 tysiąc dolarów, to zysk brutto wzrośnie o około 2 tysiące dolarów, przy złożeniu, że pozostałe zmienne się nie zmienią;
  • Jeśli koszt promocji bezpośrenidej wzrośnie o 1 tysiąc dolarów, to zysk brutto wzrośnie o około 4.67 tysiące dolarów, przy złożeniu, że pozostałe zmienne się nie zmienią;
  • Jeśli suma udzielonych rabatów wzrośnie o 1 tysiąc dolarów, to zysk brutto wzrośnie o około 1.42 tysiące dolarów, przy złożeniu, że pozostałe zmienne się nie zmienią;
  • Jeśli książka została napisana przez autora znanego (oznaczonego przez 1), to w modelu popularność autora przyjmujemy jako wartość 1 i otrzymujemy równanie:

\begin{displaymath}
zysk_{brutto}=14.33+2.56(k_{prod})+2(k_{rekl})+4.67(k_{prom})+1.42(rabaty)
\end{displaymath} Jeśli natomiast książka została napisana przez autora nieznanego (oznaczonego przez 0), to w modelu popularność autora przyjmujemy jako wartość 0 i otrzymujemy równanie: \begin{displaymath}
zysk_{brutto}=4.18+2.56(k_{prod})+2(k_{rekl})+4.67(k_{prom})+1.42(rabaty)
\end{displaymath} Wynik testu t-Studenta uzyskany dla każdej zmiennej wskazuje, że tylko koszt produkcji, koszt reklamy oraz popularność autora wywiera istotny wpływ na otrzymany zysk. Jednocześnie, dla tych zmiennych standaryzowane współczynniki $b$ są największe.

Dodatkowo, model jest dobrze dopasowany o czym świadczy: mały błąd standardowy estymacji $SE_e=8.086501$, wysoka wartość współczynnika determinacji wielorakiej $R^2=0.850974$ i poprawionego współczynnika determinacji wielorakiej $R_{adj}^2= 0.829059$ oraz wynik testu F analizy wariancji: $p<0.000001$.

Na podstawie interpretacji dotychczasowych wyników możemy przypuszczać, że część zmiennych nie wywiera istotnego wpływu na zysk i może być zbyteczna. Aby model był dobrze sformułowany interwałowe zmienne niezależne powinny być silnie skorelowane ze zmienną zależną i stosunkowo słabo pomiędzy sobą. Możemy to sprawdzić wyliczając macierz korelacji i macierz kowariancji:

Najbardziej spójną informację, pozwalającą znaleźć te zmienne w modelu, które są zbędne (nadmiarowe) niesie analiza korelacji cząstkowej i semicząstkowej i nadmiarowości:

Wartości współczynników korelacji cząstkowej i semicząstkowej wskazują, że najmniejszy wkład w budowany model mają: koszt promocji bezpośredniej i suma udzielonych rabatów. Jednak, są to zmienne najmniej skorelowane z pozostałymi w modelu, o czym świadczy niska wartość $R^2$ i wysoka wartość tolerancji. Ostatecznie, ze statystycznego punktu widzenia, modele bez tych zmiennych nie były by modelami gorszymi niż model obecny (patrz wynik testu t-Studenta dla porównywania modeli). To od decyzji badacza zależy, czy pozostawi ten model, czy zbuduje nowy model pozbawiony kosztów promocji bezpośredniej i sumy udzielonych rabatów. My pozostawiamy model obecny.

Na koniec przeprowadzimy analizę reszt. Fragment tej analizy znajduje się poniżej:

Możemy zauważyć, że jedna z reszt modelu jest obserwacją odstającą $-$ jest oddalona o więcej niż 3 odchylenia standardowe od wartości średniej. Jest to obserwacja o numerze 16. Obserwację te możemy łatwo znaleźć kreśląc wykres resz względem obserwowanych lub przewidywanych wartości zmiennej $Y$.

Ten odstający punkt zaburza założenie dotyczące homoskedastyczności. Założenie homoskedastyczności było by spełnione (tzn. wariancja reszt opisana na osi $Y$ byłaby podobna, gdy przechodzimy wzdłuż osi $X$), gdybyśmy ten punkt odrzucili. Dodatkowo, rozkład reszt nieco odbiega od rozkładu normalnego (wartość $p$ testu Lilieforsa wynosi $p=0.016415$):

Przyglądając się dokładniej punktowi odstającemu (pozycja 16 w danych do zadania) widzimy, że książka ta jako jedyna wykazuje wyższe koszty niż zysk brutto (zysk brutto = 4 tysiące dolarów, suma kosztów = (8+6+0.33+1.6) = 15.93 tysiące dolarów).

Uzyskany model możemy poprawić usuwając z niego punkt odstający. Wymaga to ponownego przeprowadzenia analizy z włączonym filtrem wykluczającym punkt odstający.

W rezultacie uzyskaliśmy bardzo podobny model, ale obarczony mniejszym błędem i lepiej dopasowany:

\begin{displaymath}
zysk_{brutto}=6.89+2.68(k_{prod})+2.08(k_{rekl})+1.92(k_{prom})+1.33(rabaty)+7.38(popul_{autora})+[4.86]
\end{displaymath} Ostatecznie zbudowany model wykorzystamy do predykcji. Na podstawie przewidywanych nakładów w wysokości: koszty produkcji $\approx 11$ tysięcy dolarów, koszty reklamy $\approx 13$ tysięcy dolarów, koszty promocji bezpośredniej $\approx 0.5$ tysiąca dolarów, suma udzielonych rabatów $\approx 0.5$ tysiąca dolarów,\\oraz faktu, że jest to autor znany (popularność autora $\approx 1$) wyliczamy przewidywany zysk brutto wraz z przedziałem ufności:

Przewidziany zysk wynosi 72 tysiące dolarów.

 
 

Porównywanie modeli liniowej regresji wielorakiej

Okno z ustawieniami opcji porównywania modeli wywołujemy poprzez menu StatystykaModele wielowymiaroweRegresja wieloraka - porównywanie modeli

Liniowa regresja wieloraka daje możliwość jednoczesnej analizy wielu zmiennych niezależnych. Pojawia się więc problem wyboru optymalnego modelu. W natłoku informacji jakie niesie zbyt duży model istnieje możliwość zagubienia ważnych informacji. Zbyt mały może pominąć te cechy, które w wiarygodny sposób mogłyby opisać badane zjawisko. Bowiem nie liczba zmiennych w modelu, ale ich jakość decyduje o jakości modelu. W wyborze zmiennych niezależnych niezbędna jest wiedza i doświadczenie związane z badanym zjawiskiem. Należy pamiętać, by w modelu znajdowały się zmienne silnie skorelowane ze zmienną zależną i słabo skorelowane między sobą.

Nie istnieje jedna prosta reguła statystyczna, która decydowałaby o liczbie zmiennych niezbędnych w modelu. Najczęściej w porównaniu posługujemy się miarami dopasowania modelu takimi jak: $R_{adj}^2$ - poprawiona wartość współczynnika determinacji wielorakiej (im wyższa wartość tym lepiej dopasowany model), $SE_e$ - błąd standardowy estymacji (im niższa wartość tym lepiej dopasowany model). W tym celu można również wykorzystać test F oparty o współczynnik determinacji wielorakiej $R^2$. Test ten służy do weryfikacji hipotezy, że dopasowanie obu porównywanych modeli jest tak samo dobre.

Hipotezy:

\begin{array}{cc}
\mathcal{H}_0: & R_{FM}^2=R_{RM}^2,\\
\mathcal{H}_1: & R_{FM}^2\ne R_{RM}^2,
\end{array}

gdzie:

$R_{FM}^2, R_{RM}^2$ $-$ współczynniki determinacji wielorakiej w porównywanych modelach (pełnym i zredukowanym).

Statystyka testowa ma postać:

\begin{displaymath}
F=\frac{R_{FM}^2-R_{RM}^2}{k_{FM}-k_{RM}}\cdot\frac{n-k_{FM}-1}{1-R_{FM}^2},
\end{displaymath}

Statystyka ta podlega rozkładowi F Snedecora z $df_1=k_{FM}-k_{RM}$ i $df_2=n-k_{FM}-1$ stopniami swobody.

Wyznaczoną na podstawie statystyki testowej wartość $p$ porównujemy z poziomem istotności $\alpha$:

\begin{array}{ccl}
$ jeżeli $ p \le \alpha & \Longrightarrow & $ odrzucamy $ \mathcal{H}_0 $ przyjmując $ 	\mathcal{H}_1, \\
$ jeżeli $ p > \alpha & \Longrightarrow & $ nie ma podstaw, aby odrzucić $ \mathcal{H}_0. \\
\end{array}

Jeśli porównywane modele nie różnią się istotnie, to powinniśmy wybrać ten z mniejszą liczbą zmiennych. Brak różnicy oznacza bowiem, że zmienne które są w modelu pełnym, a nie ma ich w modelu zredukowanym, nie wnoszą istotnej informacji. Jeśli natomiast różnica w jakości dopasowania modeli jest istotna statystycznie oznacza to, że jeden z nich (ten z większą liczbą zmiennych, o większym $R^2$) jest istotnie lepszy niż drugi.

W programie PQStat porównywanie modeli możemy przeprowadzić ręcznie lub automatycznie.

  • Ręczne porównywanie modeli- polega na zbudowaniu 2 modeli:

pełnego - modelu z większą liczbą zmiennych,

zredukowanego - modelu z mniejszą liczbą zmiennych $-$ model taki powstaje z modelu pełnego po usunięciu zmiennych, które z punktu widzenia badanego zjawiska są zbędne.

Wybór zmiennych niezależnych w porównywanych modelach a następnie wybór lepszego modelu, na podstawie uzyskanych wyników porównania, należy do badacza.

  • Automatyczne porównywanie modeli jest wykonywane w kilku krokach:

[krok 1] Zbudowanie modelu z wszystkich zmiennych.

[krok 2] Usunięcie jednej zmiennej z modelu. Usuwana zmienna to ta, która ze statystycznego punktu widzenia wnosi do aktualnego modelu najmniej informacji.

[krok 3] Porównanie modelu pełnego i zredukowanego.

[krok 4] Usunięcie kolejnej zmiennej z modelu. Usuwana zmienna to ta, która ze statystycznego punktu widzenia wnosi do aktualnego modelu najmniej informacji.

[krok 5] Porównanie modelu wcześniejszego i nowo zredukowanego.

[…]

W ten sposób powstaje wiele, coraz mniejszych modeli. Ostatni model zawiera tylko 1 zmienną niezależną.

W rezultacie, każdy model jest opisany miarami dopasowania ($R_{adj}^2$, $SE_e$), a kolejno powstające (sąsiednie) modele są porównywane testem F. Model, który zostanie ostatecznie zaznaczony jako statystycznie optymalny, to model o największym $R_{adj}^2$ i najmniejszym $SE_e$. Ponieważ jednak żadna z metod statystycznych nie potrafi w pełni odpowiedzieć na pytanie który z modeli jest najlepszy, to badacz, na podstawie uzyskanych wyników, powinien wybrać model.

Przykład c.d. (plik wydawca.pqs)

Do przewidywania zysku brutto ze sprzedaży książek wydawca planuje brać pod uwagę takie zmienne jak: koszty produkcji, koszty reklamy, koszty promocji bezpośredniej, suma udzielonych rabatów, popularność autora. Nie wszystkie te zmienne muszą jednak wpływać znacząco na zysk. Spróbujemy wybrać taki model regresji liniowej, który będzie zawierał optymalną (ze statystycznego punktu widzenia) liczbę zmiennych.

  • Ręczne porównywanie modeli.

Na podstawie zbudowanego wcześniej modelu pełnego możemy podejrzewać, że zmienne: koszty promocji bezpośredniej, suma udzielonych rabatów mają niewielki wpływ na budowany model (tzn. te zmienne nie pomagają przewidzieć wielkości zysku). Sprawdzimy czy, ze statystycznego punktu widzenia, model pełny jest lepszy niż model po usunięciu tych dwóch zmiennych.

Okazuje się, że nie ma podstaw by uważać, że model pełny jest lepszy niż model zredukowany (wartość $p$ testu F służącego porównywaniu modeli wynosi $p=0.401345$). Dodatkowo, model zredukowany jest nieco lepiej dopasowany niż model pełny (dla modelu zredukowanego $R_{adj}^2=0.82964880$, dla modelu pełnego $R_{adj}^2=0.82905898$).

  • Automatyczne porównywanie modeli.

W przypadku automatycznego porównywania modeli uzyskaliśmy bardzo podobne wyniki. Najlepszym modelem jest model o największym współczynniku $R^2_{adj}$ i najmniejszym błędzie standardowym estymacji $SE_e$. Sugerowanym najlepszym modelem jest tu model zawierający tylko 3 zmienne niezależne: koszty produkcji, koszty reklamy, popularność autora. Na podstawie powyższych analiz, ze statystycznego punktu widzenia, optymalnym modelem jest model zawierający 3 najważniejsze zmienne niezależne: koszty produkcji, koszty reklamy, popularność autora. Jednak ostateczna decyzja, który model wybrać należy do osoby posiadającą specjalistyczną widzę z zakresu badania - w tym przypadku wydawcy. Należy pamiętać, że wybrany model powinien zostać ponownie zbudowany a jego założenia zweryfikowane w oknie Regresja wieloraka.

 

Regresja logistyczna

Okno z ustawieniami opcji Regresji logistycznej wywołujemy poprzez menu StatystykaModele wielowymiaroweRegresja logistyczna

Budowany model regresji logistycznej (podobnie jak liniowej regresji wielorakiej) pozwala na zbadanie wpływu wielu zmiennych niezależnych ($X_1, X_2,..., X_k$) na jedną zmienną zależną ($Y$). Tym razem jednak zmienna zależna przyjmuje jedynie dwie wartości, np. chory/zdrowy, niewypłacalny/wypłacalny itp.

Owe dwie wartości kodowane są jako (1)/(0) gdzie:

(1) wartość wyróżniona - posiadanie danej cechy

(0) brak danej cechy.

Funkcja, na której oparty jest model regresji logistycznej wylicza nie dwupoziomową zmienną $Y$, a prawdopodobieństwo przyjęcia przez tą zmienną wyróżnionej wartości:

\begin{displaymath}
P(Y=1|X_1,X_2,...,X_k)=\frac{e^Z}{1+e^Z},
\end{displaymath} gdzie:

$P(Y=1|X_1,X_2,...,X_k)$\textendash prawdopodobieństwo przyjęcia wartości wyróżnionej (1) pod warunkiem uzyskania konkretnych wartości zmiennych niezależnych, tzw. prawdopodobieństwo przewidywane dla 1.

$Z$ najczęściej wyrażone jest zależnością liniową:

$Z=\beta_0+\sum_{i=1}^k\beta_iX_i$,

$X_1,X_2,\ldots X_k$ - zmienne niezależne, objaśniające,

$\beta_0,\beta_1,\beta_2,\ldots \beta_k$ - parametry.

Zmienne fikcyjne i interakcje w modelu

Omówienie przygotowania zmiennych fikcyjnych i interakcji przedstawiono w rozdziale Przygotowanie zmiennych do analizy w modelach wielowymiarowych.

Uwaga! Funkcja Z może być również opisana zależnością wyższego stopnia np. kwadratową - do modelu wprowadzamy wówczas zmienną zawierającą kwadrat danej zmiennej niezależnej $X_i^2$.

Logitem nazywamy przekształcenie tego modelu do postaci:

\begin{displaymath}
\ln\left(\frac{P}{1-P}\right)=Z.
\end{displaymath}

Macierze biorące udział w równaniu, dla próby o liczności $n$, zapisujemy następująco:


$
Y=\left( \begin{array}{ccc}
y_1\\
y_2\\
\vdots\\
y_n
\end{array} \right),
$

$
X=\left( \begin{array}{ccccc}
1 & x_{11} & x_{21} & \ldots & x_{k1}\\
1 & x_{12} & x_{22} & \ldots & x_{k2}\\
\vdots & \vdots & \vdots & \ddots & \vdots\\
1 & x_{1n} & x_{2n} & \ldots & x_{kn}
\end{array} \right),
$

$
\beta=\left( \begin{array}{ccc}
\beta_0\\
\beta_1\\
\beta_2\\
\vdots\\
\beta_k
\end{array}\right),
$

Rozwiązaniem równania jest wówczas wektor ocen parametrów $\beta_0,\beta_1,\ldots,\beta_k$ nazywanych współczynnikami regresji:


$
b=\left( \begin{array}{ccc}
b_0\\
b_1\\
b_2\\
\vdots\\
b_k
\end{array} \right). $

Współczynniki te szacowane są poprzez metodę największej wiarygodności czyli poprzez poszukiwanie maksimum funkcji wiarygodności $L$ (w programie użyto algorytm iteracyjny Newton-Raphson) . Na podstawie tych wartości możemy wnioskować o wielkości wpływu zmiennej niezależnej (dla której ten współczynnik został oszacowany) na zmienną zależną.

Każdy współczynnik obarczony jest pewnym błędem szacunku. Wielkość tego błędu wyliczana jest ze wzoru:

\begin{displaymath}
SE_b=\sqrt{diag(H^{-1})_b}
\end{displaymath} gdzie:

$diag(H^{-1})$ to główna przekątna macierzy kowariancji.

Uwaga! Budując model należy pamiętać, że liczba obserwacji powinna być dziesięciokrotnie większa lub równa liczbie szacowanych parametrów modelu ($n\ge 10(k+1)$).

Uwaga! Budując model należy pamiętać, że zmienne niezależne nie powinny być współliniowe. W przypadku gdy występuje współliniowość, estymacja może być niepewna a uzyskane wartości błędów bardzo wysokie. Zmienne współliniowe należy usunąć z modelu bądź zbudować z nich jedną zmienna niezależną np. zamiast współliniowych zmiennych: wiek matki i wiek ojca można zbudować zmienną wiek rodziców.

Uwaga! Kryterium zbieżności funkcji algorytmu iteracyjnego Newtona-Raphsona można kontrolować przy pomocy dwóch parametrów: limitu iteracji zbieżności (podaje maksymalną ilość iteracji w jakiej algorytm powinien osiągnąć zbieżność) i kryterium zbieżności (podaje wartość poniżej której uzyskana poprawa estymacji uznana będzie za nieznaczną i algorytm zakończy działanie).

Iloraz Szans

Jednostkowy Iloraz Szans

Na podstawie współczynników, dla każdej zmiennej niezależnej w modelu, wylicza się łatwą w interpretacji miarę jaką jest jednostkowy Iloraz Szans:

\begin{displaymath}
OR_i=e^{\beta_i}.
\end{displaymath}

Otrzymany Iloraz Szans wyraża zmianę szansy na wystąpienie wyróżnionej wartości (1), gdy zmienna niezależna rośnie o 1 jednostkę. Wynik ten jest skorygowany o pozostałe zmienne niezależne znajdujące się w modelu w ten sposób, że zakłada iż pozostają one na stałym poziomie podczas, gdy badana zmienna niezależna rośnie o jednostkę.

Wartość OR interpretujemy następująco:

  • $OR >1$ oznacza stymulujący wpływ badanej zmiennej niezależnej na uzyskanie wyróżnionej wartości (1), tj. mówi o ile wzrasta szansa na wystąpienie wyróżnionej wartości (1), gdy zmienna niezależna wzrasta o jeden poziom.
  • $OR <1$ oznacza destymulujący wpływ badanej zmiennej niezależnej na uzyskanie wyróżnionej wartości (1), tj. mówi o ile spada szansa na wystąpienie wyróżnionej wartości (1), gdy zmienna niezależna wzrasta o jeden poziom.
  • $OR\approx1$ oznacza, że badana zmienna niezależna nie ma wpływu na uzyskanie wyróżnionej wartości (1).

[Iloraz Szans - wzór ogólny]

Program PQStat wylicza jednostkowy Iloraz Szans. Jego modyfikacja, na podstawie ogólnego wzoru, umożliwia zmianę interpretacji uzyskanego wyniku.

Iloraz szans na wystąpienie stanu wyróżnionego w ogólnym przypadku jest wyliczany jako iloraz dwóch szans. Zatem dla zmiennej niezależnej $X_1$ dla $Z$ wyrażonego zależnością liniową wyliczamy:

szansę dla kategorii pierwszej:

\begin{displaymath}
Szansa(1)=\frac{P(1)}{1-P(1)}=e^Z(1)=e^{\beta_0+\beta_1X_1(1)+\beta_2X_2+...+\beta_kX_k},
\end{displaymath}

szansę dla kategorii drugiej:

\begin{displaymath}
Szansa(2)=\frac{P(2)}{1-P(2)}=e^Z(2)=e^{\beta_0+\beta_1X_1(2)+\beta_2X_2+...+\beta_kX_k}.
\end{displaymath}

Iloraz Szans dla zmiennej $X_1$ wyraża się wówczas wzorem:

\begin{displaymath}
\begin{array}{lll}
OR_1(2)/(1) &=&\frac{Szansa(2)}{Szansa(1)}=\frac{e^{\beta_0+\beta_1X_1(2)+\beta_2X_2+...+\beta_kX_k}}{e^{\beta_0+\beta_1X_1(1)+\beta_2X_2+...+\beta_kX_k}}\\
&=& e^{\beta_0+\beta_1X_1(2)+\beta_2X_2+...+\beta_kX_k-[\beta_0+\beta_1X_1(1)+\beta_2X_2+...+\beta_kX_k]}\\
&=& e^{\beta_1X_1(2)-\beta_1X_1(1)}=e^{\beta_1[X_1(2)-X_1(1)]}=\\
&=& \left(e^{\beta_1}\right)^{[X_1(2)-X_1(1)]}.
\end{array}
\end{displaymath}

Przykład

Jeśli zmienną niezależną jest wiek wyrażony w latach, to różnica pomiędzy sąsiadującymi kategoriami wieku np. 25 lat i 26 lat wynosi 1 rok $\left(X_1(2)-X_1(1)=26-25=1\right)$. Wówczas otrzymamy jednostkowy Iloraz Szans: \begin{displaymath}OR=\left(e^{\beta_1}\right)^1,\end{displaymath} który mówi o ile zmieni się szansa na wystąpienie wyróżnionej wartości gdy wiek zmieni się o 1 rok.

Iloraz szans wyliczony dla niesąsiadujących kategorii zmiennej wiek np. 25 lat i 30 lat będzie pięcioletnim Ilorazem Szans, ponieważ różnica $X_1(2)-X_1(1)=30-25=5$. Wówczas otrzymamy pięcioletni Iloraz Szans: \begin{displaymath}OR=\left(e^{\beta_1}\right)^5,\end{displaymath} który mówi o ile zmieni się szansa na wystąpienie wyróżnionej wartości gdy wiek zmieni się o 5 lat.

Uwaga!

Jeśli analizę przeprowadzamy dla modelu innego niż liniowy, lub uwzględniamy interakcję, wówczas na podstawie ogólnego wzoru możemy wyliczyć odpowiedni Ilorazu Szans zmieniając formułę wyrażającą $Z$.

Przykład c.d. (plik zadanie.pqs)

Przykład c.d. (wada.pqs)

 

Weryfikacja modelu

Istotność statystyczna poszczególnych zmiennych w modelu (istotność ilorazu szans)

Na podstawie współczynnika oraz jego błędu szacunku możemy wnioskować czy zmienna niezależna, dla której ten współczynnik został oszacowany wywiera istotny wpływ na zmienną zależną. W tym celu posługujemy się testem Walda.

Hipotezy:

\begin{array}{cc}
\mathcal{H}_0: & \beta_i=0,\\
\mathcal{H}_1: & \beta_i\ne 0.
\end{array} lub równoważnie: \begin{array}{cc}
\mathcal{H}_0: & OR_i=1,\\
\mathcal{H}_1: & OR_i\ne 1.
\end{array}

Statystykę testową testu Walda wyliczamy według wzoru: \begin{displaymath}
\chi^2=\left(\frac{b_i}{SE_{b_i}}\right)^2
\end{displaymath} Statystyka ta ma asymptotycznie (dla dużych liczności) rozkład chi-kwadrat z $1$ stopniem swobody .

Wyznaczoną na podstawie statystyki testowej wartość $p$ porównujemy z poziomem istotności $\alpha$:

\begin{array}{ccl}
$ jeżeli $ p \le \alpha & \Longrightarrow & $ odrzucamy $ \mathcal{H}_0 $ przyjmując $ 	\mathcal{H}_1, \\
$ jeżeli $ p > \alpha & \Longrightarrow & $ nie ma podstaw, aby odrzucić $ \mathcal{H}_0. \\
\end{array}

Jakość zbudowanego modelu

Dobry model powinien spełniać dwa podstawowe warunki: powinien być dobrze dopasowany i możliwie jak najprostszy. Jakość modelu regresji logistycznej możemy ocenić kilkoma miarami, które opierają się na:

$L_{FM}$ - maksimum funkcji wiarygodności modelu pełnego (z wszystkimi zmiennymi),

$L_0$ - maksimum funkcji wiarygodności modelu zawierającego jedynie wyraz wolny,

$n$ - liczności próby.

  • Kryteria informacyjne opierają się na entropii informacji niesionej przez model (niepewności modelu) tzn. szacują utraconą informację, gdy dany model jest używany do opisu badanego zjawiska. Powinniśmy zatem wybierać model o minimalnej wartości danego kryterium informacyjnego.

$AIC$, $AICc$ i $BIC$ jest rodzajem kompromisu pomiędzy dobrocią dopasowania i złożonością. Drugi element sumy we wzorach na kryteria informacyjne (tzw. funkcja straty lub kary) mierzy prostotę modelu. Zależy on od liczby zmiennych w modelu ($k$) i liczności próby ($n$). W obu przypadkach element ten rośnie wraz ze wzrostem liczby zmiennych i wzrost ten jest tym szybszy im mniejsza jest liczba obserwacji.

Kryterium informacyjne nie jest jednak miarą absolutną, tzn. jeśli wszystkie porównywane modele źle opisują rzeczywistość w kryterium informacyjnym nie ma sensu szukać ostrzeżenia.

  • Kryterium informacyjne Akaikego (ang. Akaike information criterion)

\begin{displaymath}
AIC=-2\ln L_{FM}+2k,
\end{displaymath}

Jest to kryterium asymptotyczne - odpowiednie dla dużych prób.

  • Poprawione kryterium informacyjne Akaikego

\begin{displaymath}
AICc=AIC+\frac{2k(k+1)}{n-k-1},
\end{displaymath}

Poprawka kryterium Akaikego dotyczy wielkości próby, przez co jest to miara rekomendowana również dla prób o małych licznościach.

  • Bayesowskie kryterium informacyjne Schwartza (ang. Bayes Information Criterion lub Schwarz criterion)

\begin{displaymath}
BIC=-2\ln L_{FM}+k\ln(n),
\end{displaymath}

Podobnie jak poprawione kryterium Akaikego uwzględnia wielkość próby.

  • Pseudo R$^2$ - tzw. McFadden R$^2$ jest miarą dopasowania modelu (odpowiednikiem współczynnika determinacji wielorakiej $R^2$ wyznaczanego dla liniowej regresji wielorakiej) .

Wartość tego współczynnika mieści się w przedziale $<0; 1)$, gdzie wartości bliskie 1 oznaczają doskonałe dopasowanie modelu, $0$ - zupełny bark dopasowania. Współczynnik $R^2_{Pseudo}$ wyliczamy z wzoru:

\begin{displaymath}
R^2_{Pseudo}=1-\frac{\ln L_{FM}}{\ln L_0}.
\end{displaymath}

Ponieważ współczynnik $R^2_{Pseudo}$ nie przyjmuje wartości 1 i jest wrażliwy na ilość zmiennych w modelu, wyznacza się jego poprawioną wartość:

\begin{displaymath}
R^2_{Nagelkerke}=\frac{1-e^{-(2/n)(\ln L_{FM}-\ln L_0)}}{1-e^{(2/n)\ln L_0}} \quad \textrm{lub}\quad R^2_{Cox-Snell}=1-e^{\frac{(-2\ln L_0)-(-2\ln L_{FM})}{n}}.
\end{displaymath}

  • Istotność statystyczna wszystkich zmiennych w modelu

Podstawowym narzędziem szacującym istotność wszystkich zmiennych w modelu jest test ilorazu wiarygodności. Test ten weryfikuje hipotezę:

\begin{array}{cc}
\mathcal{H}_0: & \textrm{wszystkie }\beta_i=0,\\
\mathcal{H}_1: & \textrm{istnieje }\beta_i\neq0.
\end{array}

Statystyka testowa ma postać:

\begin{displaymath}
\chi^2=-2\ln(L_0/L_{FM})=-2\ln(L_0)-(-2\ln(L_{FM})).
\end{displaymath}

Statystyka ta ma asymptotycznie (dla dużych liczności) rozkład chi-kwadrat z $k$ stopniami swobody.

Wyznaczoną na podstawie statystyki testowej wartość $p$ porównujemy z poziomem istotności $\alpha$:

\begin{array}{ccl}
$ jeżeli $ p \le \alpha & \Longrightarrow & $ odrzucamy $ \mathcal{H}_0 $ przyjmując $ 	\mathcal{H}_1, \\
$ jeżeli $ p > \alpha & \Longrightarrow & $ nie ma podstaw, aby odrzucić $ \mathcal{H}_0. \\
\end{array}

  • Test Hosmera-Lemeshowa - Test ten, dla różnych podgrup danych, porównuje obserwowane liczności występowania wartości wyróżnionej $O_g$ i przewidywane prawdopodobieństwo $E_g$. Jeśli $O_g$ i $E_g$ są wystarczająco bliskie, wówczas można założyć, że zbudowano dobrze dopasowany model.

Do obliczeń najpierw obserwacje są dzielone na $G$ podgrup - zwykle na decyle ($G=10$).

Hipotezy:


$
\begin{array}{cl}
\mathcal{H}_0: & O_g=E_g $ dla wszystkich kategorii,$\\
\mathcal{H}_1: & O_g \neq E_g $ dla przynajmniej jednej kategorii.$
\end{array}
$

Statystyka testowa ma postać:

\begin{displaymath}
H=\sum_{g=1}^G\frac{(O_g-E_g)^2}{E_g(1-\frac{E_g}{N_g})},
\end{displaymath} gdzie:

$N_g$ - liczba obserwacji w grupie $g$.

Statystyka ta ma asymptotycznie (dla dużych liczności) rozkład chi-kwadrat z $G-2$ stopniami swobody.

Wyznaczoną na podstawie statystyki testowej wartość $p$ porównujemy z poziomem istotności $\alpha$:

\begin{array}{ccl}
$ jeżeli $ p \le \alpha & \Longrightarrow & $ odrzucamy $ \mathcal{H}_0 $ przyjmując $ 	\mathcal{H}_1, \\
$ jeżeli $ p > \alpha & \Longrightarrow & $ nie ma podstaw, aby odrzucić $ \mathcal{H}_0. \\
\end{array}

  • AUC - pole pod krzywą ROC - Krzywa ROC, zbudowana w oparciu o wartość zmiennej zależnej oraz przewidywane prawdopodobieństwo zmiennej zależnej $P$, pozwala na ocenę zdolności zbudowanego modelu regresji logistycznej do klasyfikacji przypadków do dwóch grup: (1) i (0). Powstała w ten sposób krzywa, a w szczególności pole pod nią, obrazuje jakość klasyfikacyjną modelu. Gdy krzywa ROC pokrywa się z przekątną $y = x$, to decyzja o przyporządkowaniu przypadku do wybranej klasy (1) lub (0) podejmowana na podstawie modelu jest tak samo dobra jak losowy podział badanych przypadków do tych grup. Jakość klasyfikacyjna modelu jest dobra, gdy krzywa znajduje się znacznie powyżej przekątnej $y=x$, czyli gdy pole pod krzywą ROC jest znacznie większe niż pole pod prostą $y=x$, zatem większe niż $0.5$

Hipotezy:

\begin{array}{cl}
\mathcal{H}_0: & AUC=0.5, \\
\mathcal{H}_1: & AUC\neq 0.5.
\end{array}

Statystyka testowa ma postać:

\begin{displaymath}
Z=\frac{AUC-0.5}{SE_{0.5}},
\end{displaymath} gdzie:

$SE_{0.5}$ - błąd pola.

Statystyka $Z$ ma asymptotycznie (dla dużych liczności) rozkład normalny.

Wyznaczoną na podstawie statystyki testowej wartość $p$ porównujemy z poziomem istotności $\alpha$:

\begin{array}{ccl}
$ jeżeli $ p \le \alpha & \Longrightarrow & $ odrzucamy $ \mathcal{H}_0 $ przyjmując $ 	\mathcal{H}_1, \\
$ jeżeli $ p > \alpha & \Longrightarrow & $ nie ma podstaw, aby odrzucić $ \mathcal{H}_0. \\
\end{array}

Dodatkowo, dla krzywej ROC podawana jest proponowana wartość punktu odcięcia prawdopodobieństwa przewidywanego, oraz tabela podająca wielkość czułości i swoistości dla każdego możliwego punktu odcięcia.

Uwaga! Więcej możliwości w wyliczeniu punktu odcięcia daje moduł **Krzywa ROC**. Analizę przeprowadzamy na podstawie wartości obserwowanych i prawdopodobieństwa przewidywanego, które uzyskujemy w analizie regresji logistycznej.

  • Klasyfikacja

Na podstawie wybranego punktu odcięcia prawdopodobieństwa przewidywanego można sprawdzić jakość klasyfikacji. Punkt odcięcia, to domyślnie wartość 0.5. Użytkownik może zmienić tę wartość na dowolną wartość z przedziału $(0,1)$ np. wartość sugerowaną przez krzywą ROC.

W wyniku uzyskamy tabelę klasyfikacji oraz procent poprawnie zaklasyfikowanych przypadków, procent poprawnie zaklasyfikowanych (0) - swoistość oraz procent poprawnie zaklasyfikowanych (1) - czułość.

Predykcja na podstawie modelu

Na podstawie wybranego punktu odcięcia prawdopodobieństwa przewidywanego oraz zadanych wartości zmiennych niezależnych, można wyliczyć przewidywaną wartość zmiennej zależnej (0) lub (1). Punkt odcięcia, to domyślnie wartość $0.5$. Użytkownik może zmienić tę wartość na dowolną wartość z przedziału $(0,1)$ np. wartość sugerowaną przez krzywą ROC.

Przykład (plik wada.pqs)

Przeprowadzono badanie mające na celu identyfikację czynników ryzyka pewnej rzadko występującej wady wrodzonej u dzieci. W badaniu wzięło udział 395 matek dzieci z ta wadą oraz 375 matek dzieci zdrowych. Zebrane dane to: miejsce zamieszkania, płeć dziecka, masa urodzeniowa dziecka, wiek matki, kolejność ciąży, przebyte poronienia samoistne, infekcje oddechowe, palenie tytoniu, wykształcenie matki.

Budujemy model regresji logistycznej by sprawdzić które zmienne mogą wywierać istotny wpływ na występowanie wady. Jako zmienną zależną ustawiamy kolumnę GRUPA, wartością wyróżnioną w tej zmiennej jako $1$ jest grupa badana, czyli matki dzieci z wadą wrodzoną. Kolejne $9$ zmiennych, to zmienne niezależne:

MiejsceZam (2=miasto/1=wieś),

Płeć (1=mężczyzna/0=kobieta),

MasaUr (w kilogramach z dokładnością do 0.5kg),

WiekM (w latach),

KolCiąży (dziecko z której ciąży),

PoronSamo (1=tak/0=nie),

InfOddech (1=tak/0=nie),

Palenie (1=tak/0=nie),

WyksztM (1=podstawowe lub niżej/2=zawodowe/3=średnie/4=wyższe).

Jakość dopasowania modelu nie jest wysoka ($R^2_{Pseudo}=0.11$, $R^2_{Nagelkerke}=0.19$ i $R^2_{Cox-Snell}=0.14$). Jednocześnie model jest istotny statystycznie (wartość $p<0.000001$ testu ilorazu wiarygodności), a zatem część zmiennych niezależnych znajdujących się w modelu jest istotna statystycznie. Wynik testu Hosmera-Lemeshowa wskazuje na brak istotności ($p=0.2753$). Przy czym, w przypadku testu Hosmera-Lemeshowa pamiętamy o tym, że brak istotności jest pożądany, bo wskazuje na podobieństwo liczności obserwowanych i prawdopodobieństwa przewidywanego.

Interpretacja poszczególnych zmiennych w modelu zaczyna się od sprawdzenia ich istotności. W tym przypadku zmienne, które w istotny sposób są związane z występowaniem wady to:

Płeć: $p=0.0063$,

MasaUr: $p=0.0188$,

KolCiąży: $p=0.0035$,

InfOddech: $p<0.000001$,

Palenie: $p=0.0003$.

Badana wada wrodzona jest wadą rzadką, ale szansa na jej wystąpienie zależy od wymienionych zmiennych w sposób opisany poprzez iloraz szans:

  • zmienna Płeć: $OR[95\%CI]=1.60[1.14;2.22]$ - szansa wystąpienia wady u chłopca jest $1.6$ krotnie większa niż u dziewczynki;
  • zmienna MasaUr: $OR[95\%CI]=0.74[0.57;0.95]$ - im wyższa masa urodzeniowa, tym szansa wystąpienia wady u dziecka jest mniejsza;
  • zmienna KolCiąży: $OR[95\%CI]=1.34[1.10;1.63]$ - szansa wystąpienia wady u dziecka wzrasta wraz z każdą kolejną ciążą $1.34$ krotnie;
  • zmienna InfOddech: $OR[95\%CI]=4.46[2.59;7.69]$ - szansa wystąpienia wady u dziecka, gdy matka w czasie ciąży przechodziła infekcje oddechową jest $4.46$ krotnie większa niż gdyby jej nie przechodziła;
  • zmienna Palenie: $OR[95\%CI]=4.44[1.98;9.96]$ - matka paląca w czasie ciąży zwiększa $4.44$ krotnie szansę na wystąpienia wady u dziecka

W przypadku zmiennych nieistotnych statystycznie przedział ufności dla Ilorazu Szans zawiera jedynkę co oznacza, że zmienne te nie zwiększają ani nie zmniejszają szansy na wystąpienie badanej wady. Nie można więc interpretować uzyskanego ilorazu w podobny sposób jak dla zmiennych istotnych statystycznie.

Wpływ poszczególnych zmiennych niezależnych na występowanie wady możemy równiez opisać przy pomocy wykresu dotyczącego ilorazu szans:

Przykład c.d. (wada.pqs)

Zbudujemy raz jeszcze model regresji logistycznej, ale tym razem zmienną wykształcenie rozbijemy na zmienne fikcyjne (kodowanie zero-jedynkowe). Tracimy tym samym informację o uporządkowaniu kategorii wykształcenia, ale zyskujemy możliwość wnikliwszej analizy poszczególnych kategorii. Rozbicia na zmienne fikcyjne dokonano tworząc 3 zmienne dotyczące wykształcenia:

WZawodowe (1=tak/0=nie),

WŚrednie(1=tak/0=nie),

WWyższe(1=tak/0=nie).

Brak jest zmiennej dotyczącej wykształcenia podstawowego, ponieważ będzie ono stanowić kategorię odniesienia.

W rezultacie zmienne opisujące wykształcenie stają się istotne statystycznie. Dopasowanie modelu nie ulega znacznej zmianie, ale zmienia się sposób interpretacji ilorazu szans dla wykształcenia:

\begin{tabular}{|l|l|}
\hline
\textbf{Zmienna}& $OR[95\%CI]$ \\\hline
Wykształcenie podstawowe& kategoria referencyjna\\
Wykształcenie zawodowe& $0.51[0.26;0.99]$\\
Wykształcenie średnie& $0.42[0.22;0.80]$\\
Wykształcenie wyższe& $0.45[0.22;0.92]$\\\hline
\end{tabular}

Szansa na wystąpienie badanej wady w każdej kategorii wykształcenia odnoszona jest zawsze do szansy wystąpienia wady przy wykształceniu podstawowym. Widzimy, że dla bardziej wykształconych matek, iloraz szans jest niższy. Dla matki z wykształceniem:

  • zawodowym szansa wystąpienia wady u dziecka stanowi 0.51 część szansy na urodzenie dziecka z wadą jaką ma matka z wykształceniem podstawowym;
  • średnim szansa wystąpienia wady u dziecka stanowi 0.42 część szansy na urodzenie dziecka z wadą jaką ma matka z wykształceniem podstawowym;
  • wyższym szansa wystąpienia wady u dziecka stanowi 0.45 część szansy na urodzenie dziecka z wadą jaką ma matka z wykształceniem podstawowym.

Przykład (plik zadanie.pqs)

Przeprowadzono eksperyment mający na celu zbadanie umiejętność koncentracji grupy dorosłych podczas sytuacji niekomfortowych. W eksperymencie wzięło udział 130 osób. Każda badana osoba dostała pewne zadanie, którego rozwiązanie wymagało skupienia uwagi. Podczas eksperymentu niektóre osoby zostały poddane działaniu czynnika zakłócającego jakim była podwyższona temperatura powietrza do 32 stopni Celsiusza. Osoby biorące udział w eksperymencie zapytano dodatkowo o ich miejsce zamieszkania, płeć, wiek i wykształcenie. Czas na rozwiązanie zadania ograniczono do 45 minut. Dla osób, które skończyły przed czasem odnotowano rzeczywisty czas poświęcony na rozwiązanie.

Zmienna ROZWIĄZANIE (tak/nie) zawiera wynik eksperymentu, czyli informację o tym, czy zadanie zostało rozwiązane poprawnie czy też nie. Pozostałe zmienne, które mogły wpływać na wynik eksperymentu to:

MIEJSCEZAM (1=miasto/0=wieś),

PŁEĆ (1=kobieta/0=mężczyzna),

WIEK (w latach),

WYKSZTAŁCENIE (1=podstawowe, 2=zawodowe, 3=średnie, 4=wyższe),

CZAS rozwiązywania (w minutach),

ZAKŁÓCENIA (1=tak/0=nie).

Na bazie wszystkich zmiennych zbudowano model regresji logistycznej, gdzie jako stan wyróżniony zmiennej ROZWIĄZANIE wybrano tak.

Jakość jego dopasowania opisują współczynniki: $R^2_{Pseudo}=0.27$, $R^2_{Nagelkerke}=0.41$ i $R^2_{Cox-Snell}=0.30$. Na wystarczającą jakość dopasowania wskazuje również wynik testu Hosmera-Lemeshowa $(p=0.1725)$. Cały model jest istotny statystycznie o czym mówi wynik testu ilorazu wiarygodności $(p<0.000001)$.

Wartości obserwowane i prawdopodobieństwo przewidywane możemy zobaczyć na wykresie:

W modelu zmienne, które w sposób istotny wpływają na wynik to:

WIEK: $p=0.0014$,

CZAS: $p=0.0012$,

ZAKŁÓCENIA: $p=0.0001$.

Przy czym, im osoba rozwiązująca jest młodsza, czas rozwiązywania krótszy i brak jest czynnika zakłócającego, tym większe prawdopodobieństwo poprawnego rozwiązania:

WIEK: $OR[95\%CI]=0.90[0.85;0.96]$,

CZAS: $OR[95\%CI]=0.91[0.87;0.97]$,

ZAKŁÓCENIA: $OR[95\%CI]=0.15[0.06;0.37]$.

Uzyskane wyniki Ilorazu Szans przedstawiono na poniższym wykresie:

Jeśli model miałby zostać użyty do prognozowania, to należy przyjrzeć się jakości klasyfikacji. Wyliczamy w tym celu krzywe ROC.

Rezultat wydaje się zadowalający. Pole pod krzywą wynosi $AUC=0.83$ i jest istotnie większe niż $0.5$ $(p<0.000001)$, więc na podstawie zbudowanego modelu można klasyfikować. Proponowany punkt odcięcia dla krzywej ROC wynosi $0.60$ i jest nieco wyższy niż standardowo używany w regresji poziom $0.5$. Klasyfikacja wyznaczona na bazie tego punktu odcięcia daje $78.46\%$ przypadków zaklasyfikowanych poprawnie, z czego poprawnie zaklasyfikowanych wartości tak jest $77.92\%$ (czułość$[95\%CI]=77.92\% [67.02\%;86.58\%]$), wartości nie jest $79.25\%$ (swoistość$[95\%CI]=79.25\% [65.89\%;89.16\%]$).

Na tym etapie możemy zakończyć analizę klasyfikacji, lub jeśli wynik nie jest wystarczający bardziej wnikliwą analizę krzywej ROC możemy przeprowadzić w module Krzywa ROC.

Ponieważ uznaliśmy, że klasyfikacja na podstawie modelu jest zadowalająca, możemy wyliczyć prognozowaną wartość zmiennej zależnej dla dowolnie zadanych warunków. Sprawdźmy jakie szanse na rozwiązanie zadania ma osoba dla której:

MIEJSCEZAM (1=miasto),

PŁEĆ (1=kobieta),

WIEK (50 lat),

WYKSZTAŁCENIE (1=podstawowe),

CZAS rozwiązywania (20 minut),

ZAKŁÓCENIA (1=tak).

W tym celu na podstawie wartości współczynnika $b$ wyliczane jest prawdopodobieństwo przewidywane (prawdopodobieństwo uzyskania odpowiedzi tak pod warunkiem określenia wartości zmiennych zależnych):

\begin{displaymath}
\begin{array}{l}
P(Y=tak|MIEJSCEZAM,PŁEĆ,WIEK,WYKSZTAŁCENIE,CZAS,ZAKŁÓCENIA)=\\[0.2cm]
=\frac{e^{7.23-0.45\textrm{\scriptsize \textit{MIEJSCEZAM}}-0.45\textrm{\scriptsize\textit{PŁEĆ}}-0.1\textrm{\scriptsize\textit{WIEK}}+0.46\textrm{\scriptsize\textit{WYKSZTAŁCENIE}}-0.09\textrm{\scriptsize\textit{CZAS}}-1.92\textrm{\scriptsize\textit{ZAKŁÓCENIA}}}}{1+e^{7.23-0.45\textrm{\scriptsize\textit{MIEJSCEZAM}}-0.45\textrm{\scriptsize\textit{PŁEĆ}}-0.1\textrm{\scriptsize\textit{WIEK}}+0.46\textrm{\scriptsize\textit{WYKSZTAŁCENIE}}-0.09\textrm{\scriptsize\textit{CZAS}}-1.92\textrm{\scriptsize\textit{ZAKŁÓCENIA}}}}=\\[0.2cm]
=\frac{e^{7.231-0.453\cdot1-0.455\cdot1-0.101\cdot50+0.456\cdot1-0.089\cdot20-1.924\cdot1}}{1+e^{7.231-0.453\cdot1-0.455\cdot1-0.101\cdot50+0.456\cdot1-0.089\cdot20-1.924\cdot1}}
\end{array}
\end{displaymath}

W rezultacie tych obliczeń program zwróci wynik:

Uzyskane prawdopodobieństwo rozwiązania zadania wynosi $0.1215$, więc na podstawie punktu odcięcia $0.60$ przewidziany wynik to $0$ - czyli zadanie nie rozwiązane poprawnie.

 
 

Porównywanie modeli regresji logistycznej

Okno z ustawieniami opcji porównywania modeli wywołujemy poprzez menu StatystykaModele wielowymiaroweRegresja logistyczna - porównywanie modeli

Ze względu na możliwość jednoczesnej analizy wielu zmiennych niezależnych w jednym modelu regresji logistycznej, podobnie jak w liniowej regresji wielorakiej, istnieje problem wyboru optymalnego modelu. Wybierając zmienne niezależne należy pamiętać, by w modelu znajdowały się zmienne silnie skorelowane ze zmienną zależną i słabo skorelowane między sobą.

Porównując modele z różną liczbą zmiennych niezależnych zwracamy uwagę na dopasowanie modelu ($R_{Pseudo}^2$, $R^2_{Nagelkerke}$, $R^2_{Cox-Snell}$). Dla każdego modelu wyliczamy również maksimum funkcji wiarygodności, które następnie porównujemy przy użyciu testu ilorazu wiarygodności.

Hipotezy:

\begin{array}{cc}
\mathcal{H}_0: & L_{FM}=L_{RM},\\
\mathcal{H}_1: & L_{FM}\ne L_{RM},
\end{array}

gdzie:

$L_{FM}, L_{RM}$ - maksimum funkcji wiarygodności w porównywanych modelach (pełnym i zredukowanym).

Statystyka testowa ma postać:

\begin{displaymath}
\chi^2=-2\ln(L_{RM}/L_{FM})=-2\ln(L_{RM})-(-2\ln(L_{FM}))
\end{displaymath}

Statystyka ta ma asymptotycznie (dla dużych liczności) rozkład chi-kwadrat z $df=k_{FM}-k_{RM}$ stopniami swobody, gdzie $k_{FM}$ i $k_{RM}$ to ilość szacowanych parametrów w porównywanych modelach.

Wyznaczoną na podstawie statystyki testowej wartość $p$ porównujemy z poziomem istotności $\alpha$:

\begin{array}{ccl}
$ jeżeli $ p \le \alpha & \Longrightarrow & $ odrzucamy $ \mathcal{H}_0 $ przyjmując $ 	\mathcal{H}_1, \\
$ jeżeli $ p > \alpha & \Longrightarrow & $ nie ma podstaw, aby odrzucić $ \mathcal{H}_0. \\
\end{array}

Decyzję o tym, który model wybrać podejmujemy na podstawie wielkości $R_{Pseudo}^2$, $R^2_{Nagelkerke}$, $R^2_{Cox-Snell}$ oraz wyniku testu ilorazu wiarygodności porównującego kolejno powstające (sąsiednie) modele. Jeśli porównywane modele nie różnią się istotnie, to powinniśmy wybrać ten z mniejszą liczbą zmiennych. Brak różnicy oznacza bowiem, że zmienne które są w modelu pełnym, a nie ma ich w modelu zredukowanym, nie wnoszą istotnej informacji. Jeśli natomiast różnica jest istotna statystycznie oznacza to, że jeden z nich (ten z większą liczbą zmiennych, o większym $R^2$) jest istotnie lepszy niż drugi.

W programie PQStat porównywanie modeli możemy przeprowadzić ręcznie lub automatycznie.

Ręczne porównywanie modeli - polega na zbudowaniu 2 modeli:

pełnego - modelu z większą liczbą zmiennych,

zredukowanego - modelu z mniejszą liczbą zmiennych - model taki powstaje z modelu pełnego po usunięciu zmiennych, które z punktu widzenia badanego zjawiska są zbędne.

Wybór zmiennych niezależnych w porównywanych modelach a następnie wybór lepszego modelu, na podstawie uzyskanych wyników porównania, należy do badacza.

Automatyczne porównywanie modeli jest wykonywane w kilku krokach:

[krok 1] Zbudowanie modelu z wszystkich zmiennych.

[krok 2] Usunięcie jednej zmiennej z modelu. Usuwana zmienna to ta, która ze statystycznego punktu widzenia wnosi do aktualnego modelu najmniej informacji.

[krok 3] Porównanie modelu pełnego i zredukowanego.

[krok 4] Usunięcie kolejnej zmiennej z modelu. Usuwana zmienna to ta, która ze statystycznego punktu widzenia wnosi do aktualnego modelu najmniej informacji.

[krok 5] Porównanie modelu wcześniejszego i nowo zredukowanego.

[…]

W ten sposób powstaje wiele, coraz mniejszych modeli. Ostatni model zawiera tylko 1 zmienną niezależną.

Przykład c.d. (plik zadanie.pqs)

W eksperymencie badającym umiejętność koncentracji zbudowano model regresji logistycznej w oparciu o następujące zmienne:

zmienna zależna: ROZWIĄZANIE (tak/nie) - informacja o tym, czy zadanie zostało rozwiązane poprawnie czy też nie;

zmienne niezależne:

MIEJSCEZAM (1=miasto/0=wieś),

PŁEĆ (1=kobieta/0=mężczyzna),

WIEK (w latach),

WYKSZTAŁCENIE (1=podstawowe, 2=zawodowe, 3=średnie, 4=wyższe),

CZAS rozwiązywania (w minutach),

ZAKŁÓCENIA (1=tak/0=nie).

Sprawdzimy, czy wszystkie zmienne niezależne są w modelu niezbędne.

  • Ręczne porównywanie modeli.

Na podstawie zbudowanego wcześniej modelu pełnego możemy podejrzewać, że zmienne: MIEJSCEZAM i PŁEĆ mają niewielki wpływ na budowany model (tzn. na podstawie tych zmiennych nie możemy z sukcesem dokonywać klasyfikacji). Sprawdzimy czy, ze statystycznego punktu widzenia, model pełny jest lepszy niż model po usunięciu tych dwóch zmiennych.

Wynik testu ilorazu wiarygodności ($p=0.3051$) wskazuje , że nie ma podstaw by uważać, że model pełny jest lepszy niż model zredukowany. Zatem, przy nieznacznej utracie jakości modelu, miejsce zamieszkania i płeć mogą zostać pominięte.

Uwaga!

Porównania obu modeli pod względem zdolności do klasyfikacji możemy dokonać porównując krzywe ROC dla tych modeli. W tym celu wykorzystujemy moduł Zależne Krzywe ROC - porównanie opisanego w rozdziale.

  • Automatyczne porównywanie modeli.

W przypadku automatycznego porównywania modeli uzyskaliśmy bardzo podobne wyniki. Najlepszym modelem jest model zbudowany na podstawie zmiennych niezależnych: WIEK, WYKSZTAŁCENIE, CZAS rozwiązywania, ZAKŁUCENIA.

Na podstawie powyższych analiz, ze statystycznego punktu widzenia, optymalnym modelem jest model zawierający 4 najważniejsze zmienne niezależne: WIEK, WYKSZTAŁCENIE, CZAS rozwiązywania, ZAKŁUCENIA. Dokładną jego analizę możemy przeprowadzić w module Regresja Logistyczna. Jednak ostateczna decyzja, który model wybrać należy do eksperymentatora.

 
1) Savin N.E. and White K.J. (1977), The Durbin-Watson Test for Serial Correlation with Extreme Sample Sizes or Many Regressors. Econometrica 45, 1989-1996

Narzędzia strony