Analiza współzależności pomiędzy zmiennymi

Pomiędzy zjawiskami występują związki (zależności.) Nauki formułują te związki w postaci praw. Jak takie prawo naukowe powstaje? Typowo w dwu etapach, najpierw za pomocą dedukcji stawia się hipotezę, potem konfrontuje się hipotezę z danymi (podejście hipotetyczno-dedukcyjne). Na tym drugim etapie używa się statystyki (lub matematyki jeżeli prawo ma charakter deterministyczny)

Upraszczając metoda hypodedukcji sprowadza się do dedukcyjnego sformułowania hipotezy, która następnie jest empirycznie falsyfikowana, tj. próbuje się wykazać, że jest ona nieprawdziwa. Konsekwencje: nie można dowieść prawdziwości żadnej hipotezy, można natomiast wykazać, że hipoteza jest fałszywa.

Związki między cechami mogą być: funkcyjne (nauki przyrodnicze) – wartościom jednej zmiennej odpowiada tylko jedna wartość drugiej zmiennej lub stochastyczne – wartościom jednej zmiennej odpowiadają z pewnym przybliżeniem wartości innej zmiennej.

Problem: czy istnieje związek (zależność) pomiędzy cechami? Dla uproszczenie pomiędzy dwoma cechami, np. czy istnieje związek pomiędzy paleniem a chorobą nowotworową, wiekiem a prawdopodobieństwem zgonu z powodu COVID19 itd

Jaki jest charakter zależności? Jaka jest siła zależności?

Rodzaj metod zależy od sposobu pomiaru danych (nominalne, porządkowe, liczbowe)

Dwie zmienne nominalne

Tabele korelacyjne

Łączny rozkład dwóch (dwudzielna) lub większej (wielodzielna) liczby zmiennych można przedstawić w tabeli (zwanej także korelacyjną). W języku angielskim tego typu tabele noszą nazwę two-way tables albo contingency tables.

Ograniczmy się do tabel dwudzielnych, których ,,teorię’’ przedstawimy na prostym przykładzie: dla pewnej grupy osób odnotowujemy ich status-społeczno-ekonomiczny (wysoki/high, średni/middle, niski/low) oraz status-względem-palenia (wartości: pali/current, palił-nie-pali/former, nigdy-nie-palił/never).

Uwaga: status-społeczno-ekonomiczny to powiedzmy miara prestiżu używana w socjologii (można na Wikipedii doczytać co to dokładnie jest)

Obie zmienne są nominalne, obie mają po trzy wartości. Można poklasyfikować wszystkich badanych w następujący sposób:

High Low Middle Sum
current 51 43 22 116
former 92 28 21 141
never 68 22 9 99
Sum 211 93 52 356

Nieświadomie skonstruowaliśmy pierwszą tabelę korelacyjną. Taka tabela składa się z wierszy i kolumn. Dolny wiersz (Sum czyli Razem po polsku) zawiera łączną liczebność dla wszystkich wierszy w danej kolumnie. Podobnie prawa skrajna kolumna zawiera łączną liczebność dla wszystkich kolumn dla danego wiersza. Dolny wiersz/Prawą kolumnę nazywamy rozkładami brzegowymi. Pozostałe kolumny/wiersze (ale bez wartości łącznych) nazywane są rozkładami warunkowymi.

Przy warunku że osoba miała status SES równy high, 51 takich osób paliło, 92 kiedyś paliły a 68 nigdy nie paliły. Albo: Przy warunku że osoba nigdy nie paliła (never), 68 takich osób miało status high, 22 low a 9 middle.

Tabelę korelacyjną można także przedstawiać jako udziały (%):

High Low Middle Sum
current 14.32584 12.078652 6.179775 32.58427
former 25.84270 7.865169 5.898876 39.60674
never 19.10112 6.179775 2.528090 27.80899
Sum 59.26966 26.123596 14.606742 100.00000

Linus Pauling i witamina C

Linus Pauling (laureat nagrody Nobla/witamina C). Podawanie witaminy C a przeziębienie/brak przeziębienia.

Eksperyment polegał na tym, że podzielił 280 narciarzy na dwie grupy po 140 osób; przez 5–7 dni podawał witaminę C jednej grupie oraz placebo drugiej grupie; obserwował zachorowania na przeziębienie przez następne dwa tygodnie.

Eksperyment można przedstawić w postaci tablicy korelacyjnej (P/C oznacza czy narciarz zażywał witaminę czy placebo; cold/nocold czy zachorował czy nie zachorował na katar):

nocold cold razem
C 122 17 139
P 109 31 140
Sum 231 48 279

Jeden narciarz nie dokończył eksperymentu. Historia milczy dlaczego :-)

Rozkład brzegowy: Zachorował/nie zachorował (48 vs 231 albo 48/231=20,8% spośród 279 narciarzy zachorowało na katar). Takie coś przypominam nazywamy szansą. Szansa zachorowania na katar wynosi 0,208 albo jeden do 4,8

Nasze rozumowanie jest takie że jeżeli branie witaminy C nie ma wpływu na katar to zarówno ci narciarze co brali C jak i ci narciarze co brali P powinni mieć jednakowe szanse na zachorowanie. Te szanse można obliczyć na podstawie rozkładów warunkowych

Dwa rozkłady warunkowe:

Brał P→zachorował: 31 vs 109 albo 31/(31+109)=0,2214. Interpretacja: 22,14% tych którzy brali placebo zachorowało.

Brał C→zachorował: 17 vs 122 albo 17/(17+122.0) = 0.1223. Interpretacja: 12,23% tych którzy brali witaminę C zachorowało.

Zatem mamy rozbieżność powinno być 20,8% a jest 22,14% dla tych co brali P (czyli więcej) oraz 12,23% dla tych co brali C.

Na oko księgowego witamina C działa (bo jest różnica) ale dla statystyka liczy się czy ta różnica jest na tyle duża że (z założonym prawdopodobieństwem) można wykluczyć działanie przypadku.

Rozumowanie jest następujące: jeżeli prawdopodobieństwo wystąpienia tak dużej różnicy jest małe, to cechy nie są niezależne Jest to istota i jedyny wniosek z czegoś co się nazywa testem istotności-chi-kwadrat. Test chi-kwadrat porównuje liczebności tablicy korelacyjnej z idealną-tablicą-korelacyjną która zakłada niezależność jednej zmiennej od drugiej.

Można udowodnić że taka tablica powstanie przez przemnożenie dla każdego elementu tablicy odpowiadających mu wartości brzegowych a następnie podzieleniu tego przez łączną liczebność:

N Y Sum
C 115.086 23.91398 139
P 115.914 24.08602 140
Sum 231.000 48.00000 279

czyli dla pierwszego elementu było to: 139 * 231 / 279 = 115,086. Proszę zwrócić uwagę że rozkłady brzegowe są identyczne, identyczna jest też łączna liczebność. Różnią się tylko rozkłady warunkowe.

W każdym arkuszu kalkulacyjnym taki test jest dostępny w postaci stosownej funkcji. Argumentem tej funkcji jest obszar zawierający tabelę korelacyjną (bez rozkładów brzegowych).

W przypadku tablicy Paulinga wykonanie testu da w wyniku:

## [1] "0.041864"

powyższe 0.0418644 oznacz że wystąpienie tak dużych różnic pomiędzy oczekiwanymi przy założeniu o niezależności zmiennych liczebnościami a obserwowanymi liczebnościami zdarza się około 4 razy na 100. Formalnie 0.0418644 to prawdopodobieństwo oczywiście.

Przypominam ideę testu: jeżeli prawdopodobieństwo zaobserwowanych różnic jest małe to albo zakładamy że

  • albo mamy pecha i pięć razy podrzucając monetą zawsze nam spadła reszka (prawdopodobieństwo około 0,03), albo

  • albo że założenie co do niezależności jest fałszywe.

Statystyk zawsze wybierze drugie. Pozostaje tylko ustalenie co to znaczy małe.

Małe to takie które jest mniejsze od arbitralnie przyjętego przez statystyka. Zwykle jest to 0,05 lub 0,01 (czasami 0,1) co oznacza że odrzucając założenie o braku związku pomiędzy katarem a braniem witaminy C pomylimy się pięć lub raz na 100.

W statystyce nie ma wyników ze 100% pewnością Statystyk zawsze musi przyjąć prawdopodobieństwo popełnienia błędu zwane poziomem istotności (w sensie, że zaobserwowane wielkości – w tym przypadku odchyleń – są istotnie różne od przypadkowych).

Zatem w tym konkretnym przypadku: na poziomie istotności 0,05 hipotezę o niezależności cech należy odrzucić ale na poziomie 0,01 już nie. (Dlaczego?)

Dwie zmienne liczbowe

W tym przypadku dobrze jest rozpocząć analizę od wykresu

Korelacyjny wykres rozrzutu (korelogram, wykres XY w Excelu, scatter plot)

W układzie kartezjańskim każdej obserwacji odpowiada kropka o współrzędnych XY.

O występowaniu związku świadczy układanie się kropek według jakiegoś kształtu (krzywej). O braku związku świadczy chmura punktów niepodobna do żadnej krzywej.

Punkty układające się według prostej świadczą o zależności liniowej (wyjątek: linia pozioma lub pionowa) Punkty układające się według krzywej świadczą o zależności nieliniowej.

Przykład: Czy istnieje zależność pomiędzy zamożnością a spożyciem mięsa? Dysponujemy danymi dotyczącymi spożycia mięsa per capita oraz GDP per capita (w 2013 roku).

Pomiar siły zależności: współczynnik korelacji liniowej Pearsona

Kowariancja to suma iloczynów odchyleń wartości zmiennych \(XY\) od ich wartości średnich:

\[cov (xy) = \frac{1}{n} \sum_{i=1}^N (x - \bar x) (y - \bar y)\]

Kowariancja zależy od rozproszenia (im większe tym większa), ma też dziwną jednostkę (jednostkaX · jednostkaY) oraz zależy od wybranych skal (tony vs gramy na przykład)

https://tinystats.github.io/teacups-giraffes-and-statistics/05_correlation.html

Dlatego do pomiaru związku pomiędzy cechami nie używa się kowariancji, ale współczynnika korelacji (liniowej, Pearson linear correlation coefficient):

\[r_xy = cov(xy) / (S_x \cdot S_y)\]

Współczynnik jest miarą niemianowaną, o wartościach ze zbioru \([-1;1]\); Skrajne wartości \(\pm 1\) świadczą o związku funkcyjnym (wszystkie punkty układają się na linii prostej); wartość zero świadczy o braku związku (linia pozioma/pionowa)

Interpretacja opisowa: wartości powyżej 0,9 świadczą o silnej zależności

Przykład: korelacja między spożyciem mięsa a GDP

Współczynnik korelacji liniowej wynosi 0.6823158 (umiarkowana korelacja).

Czy ta wartość jest istotnie różna od zera? Jest na to stosowny test statystyczny, który sprowadza się do określenia jakie jest prawdopodobieństwo otrzymania r = 0.6823158 przy założeniu że prawdziwa wartość r wynosi zero. Otóż w naszym przykładzie to prawdopodobieństwo wynosi 3.850676e-26 (czyli jest ekstremalnie małe – r jest istotnie różne od zera).

Pomiar siły zależności: regresja liniowa

Regresja liniowa zakłada że istnieje związek przyczyna-skutek i ten związek można opisać linią prostą (stąd liniowa). Skutek jest jeden i nazywa się go zmienną zależną a przyczyn może być wiele i noszą nazwę zmiennych niezależnych. W przypadku gdy związek dotyczy dwóch zmiennych mówi się o regresji prostej a równanie prostej można zapisać jako:

\[Y = a \cdot X + b\]

Regresja wieloraka to zależność pomiędzy jedną zmienną objaśnianą a wieloma zmiennymi objaśniającymi:

\[Y = a_1 \cdot X_1 + a_2 \cdot X_2 + a_3 \cdot X_3 + ... + b\]

Omówimy tylko model regresji prostej.

Współczynnik \(a\) ma jasną interpretację: jeżeli wartość zmiennej \(X\) rośnie o jednostkę to wartość zmiennej \(Y\) zmienia się przeciętnie o \(a\). Wyraz wolny zwykle nie ma sensu (formalnie jest to wartość zmiennej \(Y\) dla \(X=0\))

Oznaczmy przez \(y_i\) wartości obserwowane (empiryczne) a przez \(\hat y_i\) wartości teoretyczne (leżące na prostej linii regresji). Wartości \(a\) oraz \(b\) wyznacza się minimalizując: \(\sum_{i=1}^N (\hat y_i - y_i)^2\) (suma kwadratów odchyleń wartości empirycznych od wartości teoretycznych ma być minimalna)

Ocena dopasowania: średni błąd szacunku

Średni błąd szacunku (MSE mean squared error). Oznaczając resztę jako: \(e_i = y_i - \hat y_i\), definiujemy średni błąd szacunku jako: \(S_e = \sqrt {\sum (e_i^2 / n -k)}\). (Przy okazji \(S_e^2\) nazywamy wariancją resztową)

Ocena dopasowania: współczynniki zbieżności i determinacji

Suma kwadratów reszt (albo odchyleń wartości teoretycznych od wartości empirycznych, albo suma kwadratów błędów vel resztowa suma kwadratów): \(RSK = \sum (y_i - \hat y_i)^2\). Suma kwadratów odchyleń wartości empirycznych od średniej (ogólna suma kwadratów): \(OSK = \sum (y_i - \bar y)^2\) Suma kwadratów odchyleń wartości teoretycznych od średniej (wyjaśniona suma kwadratów): \(WSK = \sum (\hat y_i - \bar y)^2\)

Można wykazać, że \(OSK = WSK + RSK\). Współczynnik zbieżności to \(R^2 = WSK/OSK\). Współczynnik determinacji to \(\Phi^2 = RSK/OSK\).

Współczynniki przyjmują wartość z przedziału \([0,1]\) lub \([0, 100]\)% Interpretacja współczynników: procent zmienność wyjaśniane/nie wyjaśniane przez linię regresji.

Ocena dopasowania: istotność parametru \(a\)

Zerowa wartość parametru \(a\)

Jeżeli: \(Y= 0 \cdot X + b\), to \(Y = b\) czyli nie ma zależności pomiędzy \(XY\). Wartości \(a\) bliskie zero wskazują na słabą zależność pomiędzy cechami. Przyjmijmy

\[\hat Y = a \cdot \hat X + b + e\]

w tym zapisie obserwowana linia regresji jest równa linii teoretycznej \(Y = a \cdot X + b\) plus błąd losowy. (przy czym losowość błędu przejawia się m.in w tym, że nie ma on charakteru systematycznego, tj. jest określona jako zero plus/minus cośtam)

Można teraz zadać pytanie jeżeli \(a=0\) to jakie jest prawdopodobieństwo, że współczynnik \(\hat a\) oszacowany na podstawie \(n\) obserwacji mierzonych z błędem \(e\) będzie (co do wartości bezwzględnej) większy niż \(a_e\). Albo inaczej: otrzymaliśmy \(a_e\), jakie jest prawdopodobieństwo otrzymania takiej wartości (lub mniejszej co do wartości bezwzględnej) przy założeniu że istotnie \(a=0\).

Jeżeli takie prawdopodobieństwo jest duże to zakładamy że być może \(a=0\) a jeżeli małe to że \(a \not0\). Duże/małe przyjmujemy arbitralnie, zwykle jest to \(0,1\), \(0,05\) lub \(0,01\). Tak zgadza się to prawdopodobieństwo to poziom istotności

Testowanie istotności współczynnika regresji jest ważnym kryterium oceny jakości dopasowania. Regresja z nieistotnym współczynnikiem nie może być podstawą do interpretowania zależności pomiędzy XY.

W każdym arkuszu kalkulacyjnym jest funkcja umożliwiająca oszacowanie regresji liniowej. Funkcja ta oblicza parametry prostej, średni błąd szacunku, współczynnik determinacji oraz testuje istotność współczynnika \(a\).

Przykład: GDP a spożycie mięsa

Podobno im większe GDP tym więcej ludzie jedzą mięsa:

\[MPC = 34.08 + 0,001 \cdot GDP\] Każdy USD /per capita więcej dochodu narodowego (GDP) oznacza przeciętny wzrost spożycia mięsa o 0,001 kg. Przeciętna różnica wartości teoretycznych od empirycznych wynosi 21,04 kg (średni błąd szacunku). Współczynnik zbieżności wynosi 41%. Współczynnik nachylenia prostej (mimo że jego wartość wynosi zaledwie 0,001) jest statystycznie istotny.

Dwie zmienne co najmniej porządkowe

Pomiar siły zależności: współczynnik korelacji rang

Współczynnik korelacji rang (Spearmana vel Spearman’s Rank-Order Correlation) może być stosowany w przypadku gdy cechy są mierzone w skali porządkowej (lub lepszej)

Obliczenie współczynnika Spearmana dla \(N\) obserwacji na zmiennych \(XY\) polega na zamianie wartości \(XY\) na rangi (numery porządkowe od \(1...N\)). Następnie stosowana jest formuła Pearsona, tj. (\(\tau_x\) oraz \(\tau_y\) oznaczają rangi):

\[\rho_{xy} = cov(\tau_x, \tau_y) / ( S_{\tau_x} \cdot S_{\tau_y} )\]

Współczynnik \(\rho_{xy}\) to także miara niemianowana, o wartościach ze zbioru [-1;1];

W podręcznikach można też spotkać formułę alternatywną:

\[r_s = 1 - \frac{6 \cdot \sum_\limits{i=1}^N d_i^2}{ N(N^2 -1)}\]

gdzie \(d_i\) oznacza różnicę między rangami dla cechy \(X\) oraz \(Y\), tj. \(d_i = \tau_{x_i} - \tau_{y_i}\). Wymagany jest aby \(d_i\) były różne od siebie. Można też spotkać (w podręcznikach) przykłady użycia formuły alternatywnej dla rang identycznych, wówczas (dla \(k\) identycznych rang, zamienia się je na wartości średnie) \(d_i = \frac{1}{k}\sum_\limits{i=1}^k d_i^k\). Raczej nie zalecane…

Przykład: spożycie mięsa

Współczynnik Pearsona i Spearmana dla zależności między spożyciem mięsa w 1980 a spożyciem mięsa w 2013 roku (zmienna objaśniana):

## [1] "współczynnik Pearsona: 0.68"
## [1] "współczynnik Spearmana: 0.68"

Nie ma sensu liczenia współczynnika korelacji rang w przypadku kiedy obie cechy są liczbami, bo wtedy należy użyć normalnego współczynnika Pearsona. Ale nie jest to też błędem więc w powyższym przykładzie go liczymy :-)

Współczynnik korelacji liniowej Spearmana wynosi 0.6845429 (umiarkowana korelacja).

Czy ta wartość jest istotnie różna od zera? Jest na to stosowny test statystyczny, który sprowadza się do określenia jakie jest prawdopodobieństwo otrzymania \(r_s\) = 0.6845429 przy założeniu że prawdziwa wartość \(r_s\) wynosi zero. Otóż w naszym przykładzie to prawdopodobieństwo wynosi 2.302116e-26 (czyli jest ekstremalnie małe – \(r_s\) jest istotnie różne od zera).

Zmienna liczbowa i zmienna nominalna

W tym przypadku najprostsza metoda ustalenia związku polega na porównaniu średnich w grupach określonych przez wartości zmiennej nominalnej.

Przykładowo chcemy zweryfikować czy istnieje zależność pomiędzy wielkością stresu (mierzoną za pomocą skali liczbowej) a rodzajem miejsca pracy (skala nominalna; dwie wartości: szpital, przychodnia) w pewnej grupie pielęgniarek.

Czy przeciętna wielkość stresu jest różna w szpitalach i przychodniach co by świadczyło o zależności pomiędzy stresem a rodzajem miejsca pracy.

W tego typu sytuacjach posługujemy się testami weryfikującymi czy zaobserwowane różnice pomiędzy średnimi są istotne czy nie. Jak są istotne to jest zależność a jak nie są to jej nie ma.

Jeżeli wartości zmiennej nominalnej są dwie co się często zdarza możemy zastosować test zwany U Manna-Whitneya (dostępny w każdym arkuszu kalkulacyjnym); jeżeli wartości zmiennej nominalnej są więcej niż dwie stosujemy test Kruskalla-Wallisa (także dostępny w każdym arkuszu kalkulacyjnym)

Przykład: depresja a miejsce pracy

Studenci pielęgniarstwa i ratownictwa PSW w 2023 roku wypełnili ankietę zawierającą test depresji Becka (wartość liczbowa) oraz pytanie o rodzaj miejsca pracy (skala nominalna):

m-pracy średnia n
Przychodnia 7.833333 12
Szpital 8.450549 91

Ponieważ grup jest dokładnie 2 stosujemy test U Manna-Whitneya.

Prawdopodobieństwo wystąpienia tak dużej różnicy przy założeniu, że średnie w obu grupach są identyczne wynosi 0.8528214 (różnica jest zatem nieistotna; obie średnie są identyczne–nie ma zależności)

W ankiecie było także pytanie o staż pracy:

staż (lata) średnia n
<06 8.512821 39
>19 8.533333 45
07-12 7.857143 7
13-18 7.666667 12

Ponieważ grup jest > 2 stosujemy test Kruskalla-Wallisa.

Prawdopodobieństwo tak dużych różnic w średnich przy założeniu, że średnie we wszystkich grupach są identyczne wynosi 0.678923 (różnice są zatem nieistotne; wszystkie średnie są identyczne–nie ma zależności)

Podsumowanie

Przedstawiono 6 następujących metod (każda jest dostępna w dowolnym arkuszu kalkulacyjnym):

  1. korelogram

  2. tablica korelacyjna/test chi-kwadrat

  3. współczynnik korelacji Pearsona

  4. współczynnik korelacji Spearmana

  5. regresję liniową

  6. test U Manna-Whitneya albo test Kruskalla-Wallisa