wtorek, 15 maja 2012

Statystyka opisowa: dane oraz przykładowe zadania

Zaliczenia przedmiotu Statystyka opisowa

Warunkiem zaliczenia przedmiotu Statystyka opisowa jest przygotowanie w grupie maksimum 3 osób pracy, w której:

  • Należy wszechstronnie porównać (średnie, rozproszenie, asymetria, wykresy) rozkład cechy dla dwóch (lub więcej) różnych zbiorowości. Można wykorzystać dane indywidualne (przykład 1 poniżej) lub pogrupowane (przykład 2), albo:

  • Zamiast badania danych przekrojowych można przeprowadzić analizę dynamiki szeregu czasowego (obejmującą: średnie, różne wskaźniki dynamiki, trend+sezonowość (jeżeli występuje), wykresy), por. przykład 3.

Uwaga ważne: należy wykonać niezbędne obliczenia za pomocą arkuszy OOCalc (zalecany) lub MS Excel a opis i interpretację uzyskanych wyników zawrzeć w dokumencie w formacie OOWriter/MS Word. Opis powinien bezwzględnie zawierać stronę tytułową z tytułem pracy oraz  nazwiskami i imionami autorów.

Oba pliki (arkusz+opis) należy wysłać elektronicznie na wiadomy adres. Dodatkowo proszę o wydrukowanie i dostarczenie opisu w formie papierowej. Sam plik arkusza bez niezbędnego opisu/interpretacji nie wystarczy do uzyskania zaliczenia.

Przykładowe zadania

Przykład 1: Liczba głosów oddanych na posłów wybranych do Sejmu 7 kadencji. Porównamy dwie partie określające się jako lewicowe (cokolwiek to znaczy), tj. RP (Ruch Poparcia Palikota) oraz SLD (Sojusz Lewicy Demokratycznej)

Komentarz do obliczeń wykonanych w programie oocalc

Pobieramy dane i zapisujemy w odpowiednim formacie. Arkusz zawiera dane dotyczące wszystkich posłów. Aby oddzielić posłów RP/SLD od reszty można po prostu posortować (Dane→Sortuj) względem kolumny ,,I''. Następnie metodą kopiuj-wklej przenosimy, to co jest potrzebne do oddzielnego arkusza. Analizę przeprowadzamy oddzielnie dla posłów RP i SLD.

Przykładowy arkusz w formacie OOCalc jest tutaj.

Aby wykreślić histogram należy dane pogrupować w szereg rozdzielczy. Do tego celu służy funkcja CZĘSTOŚĆ, której składnia jest następująca:

CZĘSTOŚĆ (obszar-danych ; obszar-końców-przedziałów)

W omawianym przykładzie obszar-danych, to kolumna zawierająca liczbę oddanych na posła głosów (log). Argument obszar-końców-przedziałów to obszar górnych końców przedziałów szeregu rozdzielczego. Funkcja CZĘSTOŚĆ zwraca obszar o jedną komórkę większy od obszaru górnych końców przedziałów -- liczebność tej komórki, to liczba elementów większych od wartości ostatniego górnego końca przedziału. Krótki film ilustruje jak należy używać funkcji CZĘSTOŚĆ.

Uwaga: funkcja CZĘSTOŚĆ jest specjalna (bo zwraca obszar): po jej wpisaniu do komórki należy nacisnąć Ctrl-Shift-Enter a nie zwyczajne Enter.

Kształt histogramu (a co za tym idzie wnioski dotyczące kształtu rozkładu) zależą w dużym stopniu od sposobu pogrupowania danych w szereg rozdzielczy. Kolejny arkusz przedstawia kilka histogramów dla różnych wariantów grupowania danych.

Do obliczenia średniej, mediany, kwartyli, wariancji, odchylenia standardowego itp. służą odpowiednie funkcje. Przykładowo poniższa funkcja wyznaczy wartość mediany:

MEDIANA (obszar-danych)

Interpretacja otrzymanych wyników:

Średni liczba głosów oddanych na posła SLD wyniosła 13075 a na posła RP 14894. Jeżeli chodzi o posłów SLD to 50% zdobyło nie więcej niż 9941 głosów, 25% z nich zdobyło nie więcej niż 8041 głosów a 75% -- zdobyło nie więcej niż 15597 głosów Jeżeli chodzi o posłów RP to 50% zdobyło nie więcej niż 12828 głosów, 25% z nich zdobyło nie więcej niż 9984 głosów a 75% -- zdobyło nie więcej niż 15837 głosów Wielkości wszystkich miar średnich (oraz analiza histogramu) wskazuje, iż przeciętnie posłowie RP wybrani do Sejmu 7 kadencji zdobyli więcej głosów.

Przeciętne odchylenie liczby zdobytych głosów od średniej arytmetycznej wyniosło dla posłów SLD 9313,7 głosów a dla posłów RP 12825,5 głosów. W wartościach bezwzględnych bardziej jednorodną grupą wydają się być zatem posłowie SLD. Także wartości względne, tj. wartości współczynników zmienności (klasycznych) wynoszące odpowiednio 71,2% (SLD) oraz 90,6% (RP) wskazują na większe zróżnicowanie grupy posłów RP.

Wartości współczynników zmienności opartych o medianę i kwartyle wskazują natomiast, że zróżnicowanie liczby głosów w grupie posłów SLD ($V_q = 76%$ ) jest większe niż w grupie posłów RP ($V_q = 45%$). Podobnie wygląda wskaźnik definiowany jako $V_Q=(Q_3-Q_1)/(Q_3+Q_1)$, którego wartość dla posłów SLD wynosi 32% a dla posłów RP 22,7%.

Skąd taka dziwna wartość współczynnika $V_s$?

Obliczmy średnie i miary rozproszenia dla obu grup posłów pomijając obserwacje skrajnie nietypowe (cf. jeszcze inny arkusz). W każdej z obu grup jest jedna nietypowa obserwacja: wynik lidera listy (J. Palikot i R. Kalisz). Po usunięciu tej jednej obserwacji wszystkie współczynniki zmienności, bez wyjątku wykazują na większą zmienność grupy posłów SLD. Wniosek: miary klasyczne nie są właściwe w przypadku rozkładów zawierających wielkości nietypowe i/lub rozkładów znacznie odbiegających kształtem od rozkładu jednomodalnego.

Uwaga dla studentów: proszę obliczyć wszystkie miary a ew. zdyskwalifikować miary klasyczne jako wniosek z analizy a nie z góry założyć że miary klasyczne są do kitu i ich nie liczyć. Za takie postępowanie zostanie obniżona ocena!

Kształt rozkładów (histogram) wskazuje, że rozkłady cechuje asymetria dodatnia, przy czym asymetria w przypadku posłów SLD jest zauważalnie większa a dla posłów RP wydaje się być nieduża. Obliczone wielkości klasycznego współczynnika asymetrii ($\mu$) wskazują na ogromną asymetrię obu rozkładów ($\mu> 3$ dla SLD oraz $\mu> 5$ dla RP). Współczynniki Pearsona są znacznie niższe ($\mu = 38$ dla SLD oraz $\mu> 0,03$ dla RP) i wskazują na -- zgodną ze stanem faktycznym (wykres) sytuacją.

Wniosek: także klasyczne miary asymetrii nie są właściwe w przypadku rozkładów zawierających wielkości nietypowe i/lub rozkładów znacznie odbiegających kształtem od rozkładu jednomodalnego. Wniosek ten potwierdza analiza skośności po usunięciu wartości nietypowych (J. Palikot i R. Kalisz).

Przykład 2: porównanie struktury wieku posłów dwóch największych partii wybranych do Reichstagu w wyborach z 1930 r. (źródło cf. Młodzi, wykształceni z wielkich miast).

Obliczenia są w miarę oczywiste i nie wymagają dodatkowego komentarza. Stosowny arkusz jest tutaj.

Uwaga: analizowany szereg ma nieokreślone: dolny koniec pierwszego oraz górny koniec ostatniego przedziału. Arbitralnie przyjęliśmy 20 lat jako dolny koniec pierwszego a 80 lat jako górny koniec ostatniego przedziału. Ponieważ liczebności tych przedziałów są niewielkie ewentualny popełniony błąd także nie będzie duży.

Interpretacja otrzymanych wyników:

Średni wiek posła NSDAP wynosił 38,8 lat a posła KPD był nieco niższy bo wynosił 37,35 lat. 50% posłów NSDAP miało 37,1 lat i mniej. 25% posłów NSDAP nie było starszych niż 32,5 lat (pierwszy kwartyl, $Q_1$) a 75% nie było starszych niż 44,1 lat (trzeci kwartyl, $Q_3$). 50% posłów KPD miało 36,7 lat i mniej. 25% posłów KPD nie było starszych niż 32,4 lat a 75% nie było starszych niż 42,1 lat.

Wielkości wszystkich miar średnich (oraz analiza histogramu) wskazuje, iż przeciętnie posłowie KPD wybrani do Reichstagu w wyborach z 1930 r. byli nieco młodsi od posłów NSDAP.

Przeciętne odchylenie wieku od średniej arytmetycznej wyniosło dla posłów NSDAP 9,1 lat a dla posłów KPD 7,12 lat co wskazuje że posłowie KPD są grupą bardziej jednorodną. Potwierdzają to wartości względne, tj. wartości współczynników zmienności (klasycznych) wynoszące odpowiednio 23,5% (NSDAP) oraz 19,1% (KPD). Wartości współczynników zmienności opartych o medianę są nieco wyższe i wynoszą odpowiednio: 31,1% oraz 26,3%.

Kształt rozkładów (histogram) oraz różnica pomiędzy wielkościami miar średnich wskazuje, że rozkłady cechuje niewielka asymetria dodatnia. Co potwierdzają obliczone wielkości współczynników zarówno klasycznego ($\mu$) jak i Pearsona wykorzystującego różnice pomiędzy średnimi ($W_s$)

Przykład 3: Kwartalny skup mleka w mln litrów w Polsce w latach 1996--2003 (32 obserwacje)

Komentarz do obliczeń wykonanych w programie oocalc

Obliczenia są zawarte w arkuszu dostępnym tutaj.

Ponieważ z wykresu wynika, że mamy do czynienia z sezonowością roczną obliczamy średnią ruchomą 4-okresową (kolumna ,,D'' arkusza).

Kolumna ,,M'' zawiera obliczenia potrzebne do wyznaczenia wartości średniego błędu kwadratowego

Kolumny ,,N'' do ,,T'' zawierają obliczenia niezbędne do wyznaczenia parametrów trendu liniowego (i oceny dopasowania tegoż trendu) metodą najmniejszych kwadratów.

Kolumny ,,V'' do ,,Y'' zawierają obliczenia niezbędne do wyznaczenia wskaźników sezonowości. Wskaźniki surowe zawiera wiersz 37 a oczyszczone wiersz 40 (w kolumnach ,,V''--,,Y'')

Kolumna ,,Z'' zawiera obliczone wskaźniki sezonowości.

Skorygowane o wskaźniki sezonowości wartości funkcji trendu liniowego zawiera kolumna ,,AA''. Kolumna ,,AB'' zaś niezbędne obliczenia dla wyznaczenia RMSE (pierwiastek błędu średniokwadratowego) dla szeregu skorygowanego przez przemnożenie przez odpowiednie wskaźniki sezonowości.

Interpretacja otrzymanych wyników:

Funkcja trendu skupu mleka w Polsce w latach 1993--2003 jest następująca $$ \hat y_t = 1591849,4 + 5582,4 t $$ Z funkcji trendu wynika, że w latach 1996--2003 skup mleka w Polsce wzrastał z kwartału na kwartał średnio o 5582,4 mln litrów. Wyraz wolny funkcji trendu (1591849,4) informuje o teoretycznej wielkości skupu mleka w pierwszym kwartale 1996 r.

Odchylenie standardowe składnika losowego wynosi 258448,5 mln litrów. Średnie błędy szacunku parametrów są zaś równe $\bar \beta = 4948,3$ oraz $\bar \alpha = 93559,9$. Wielkość standardowe składnika losowego zwłaszcza w przypadku parametru $\beta$ jest znaczna co świadczy o słabym dopasowaniu do danych empirycznych. Można oczekiwać, iż test istotności parametru strukturalnego $\beta \neq 0$ wykaże iż nie ma podstaw do odrzucenia $H_0$ (sprawdzić samodzielnie!)

Słabe dopasowanie potwierdzają wielkości współczynników $\phi^2$ oraz $R^2$ wynoszące odpowiednio 96% oraz 4%, tj. 96% procent zmienności skupu mleka nie jest objaśniane przez liniową funkcję trendu. Albo: zaledwie 4% procent zmienności skupu mleka jest objaśniane przez liniową funkcję trendu.

Porównując wielkości współczynnika RMSE (pierwiastek błędu średniokwadratowego) dla liniowej funkcji trendu (RMSE_tl) oraz średniej ruchomej otrzymamy: $$ \begin{aligned} RMSE_{sr} &= \sqrt{53886700/30} =1796223.33 \\ RMSE_{tl} &= \sqrt{2003868574371/32} = 62620892949.09 \end{aligned} $$ co świadczy o dużo lepszym dopasowaniu średniej ruchomej do danych empirycznych (powodem jest oczywiście sezonowość--trend liniowy w takiej sytuacji dużo gorzej objaśnia kształtowanie się skupu mleka)

Sezonowość. Wielkość skupu mleka w pierwszym kwartale stanowiła 82,7% przeciętnej kwartalnej, zaś w drugim, trzecim i czwartym kwartale odpowiednio 110,6%, 116,9 oraz 89,8% przeciętnej kwartalnej.

Obliczenie RMSE dla danych skorygowanych dało w wyniku $232363367977,161/32 = 85213,6$ (dla danych nieskorygowanych było to $2003868574371,33/32 = 250241,7$ czyli znacznie więcej.

Dane

Pliki z przykładowymi danymi (w formacie CSV dla ułatwienia importu do różnych arkuszy i/lub innych narzędzi) oraz wyżej cytowane arkusze (w formacie OOCalc) są dostępne tutaj.

Literatura

Dokumentacja OOffice (niestety w języku).

poniedziałek, 14 maja 2012

Analiza dynamiki zjawisk

Średnia ruchoma (moving average)

Oblicznie $k$-okresowej średniej ruchomej w przypadku gdy $k$ jest liczbą nieparzystą: $$ y_{t + (k-1)/2} = \frac{1}{k} \sum_{i=t}^{t + k -1} y_i \quad t = 1,\dots,n - \frac{k-1}{2} \label{S-MovingAverage} $$

Przykładowo powyższy wzorek można rozpisać dla $k=3$ następująco: $$ \begin{aligned} \bar y_2 &= \frac{y_1 + y_2 + y_3}{3} \\ \bar y_3 &= \frac{y_2 + y_3 + y_4}{3} \\ \dots & \dots \\ \bar y_{n-1} &= \frac{y_{n-2} + y_{n-1} + y_n}{3} \end{aligned} $$

Oblicznie $k$-okresowej średniej ruchomej w przypadku gdy $k$ jest liczbą parzystą (średnia ruchoma scentrowana): $$ y_{t + k/2} = \frac{1}{k} (0,5 y_t + \sum_{i=t}^{t+k-1}y_i + 0,5 y_{t+k}) \quad t=1,\dots,n - \frac{k}{2} \label{S-MovingAverage-C} $$

Przykładowo powyższy wzorek można rozpisać dla $k=4$ następująco: $$ \begin{aligned} \bar y_3 &= \frac{0,5 y_1 + y_2 + y_3 + y_4 + 0,5 y_5}{4} \\ \bar y_4 &= \frac{0,5 y_2 + y_3 + y_4 + y_5 + 0,5 y_6}{4} \\ \dots & \dots \\ \bar y_{n-2} &= \frac{0,5 y_{n-4} + y_{n-3} + y_{n-2} + y_{n-1} + 0,5 y_n}{4} \end{aligned} $$

Ocena jakości dopasowania

Błędy prognoz ex-post

Średni błąd kwadratowy ($n^*=n -(k-1)/2$ lub $n^*=n -k/2$) w zależności od rodzaju średniej ruchomej: $$ MSE = \frac{1}{n^*} \sum_{i=1}^{n^*} (y_i -\bar y_i)^2 $$

Pierwiastek błędu średniokwadratowego: $$ RMSE = \sqrt{MSE} $$

Wykorzystanie modelu średniej ruchomej do prognozowania polega na wyznaczeniu prognozy jako średniej arytmetycznej zwykłej bądź ważonej z $k$ ostatnich wartości zmiennej, tj. $y^*_t = \frac{1}{k}\sum_{i=t-k}^{t-1}y_i$ przy czym $k$ należy wyznaczyć tak aby średni kwadratowy błąd ex post był minimalny.

Metoda analityczna

Liniowa funkcja trendu ma postać: $$y_t = \alpha + \beta t + \xi_t$$ gdzie: $y_t$ -- poziom zjawiska w okresie $t$; $\alpha$, $\beta$ -- parametry; $t$ -- zmienna czasowa (np. $t= 1,...n$) oraz $\xi_t$ -- składnik losowy.

Wartości (oceny) parametrów wyznaczone metodą najmniejszych kwadratów są następujące: $$ \begin{aligned} \hat \beta &= \frac{\sum_{t=1}^n (y_t -\bar y)(t - \bar t) }{ \sum_{i=1}^n (t - \bar t)^2 } \\ \hat \alpha &= \bar y - \bar \beta t \end{aligned} $$

Ocena jakości dopasowania

Standardowe odchylenie składnika resztowego: $$ s_{\xi} = \sqrt { \frac{1}{n-2} \sum_{t=1}^n (y_t - \hat y_t)^2 } $$

Odchylenie ocen parametrów (błąd standardowy oceny) $$ \begin{aligned} s(\beta) &= \sqrt { \frac{ s_{\xi}^2 }{\sum (t-\bar t)^2} } \\ s(\alpha) &= s_{\xi} \sqrt { \frac{\sum t^2}{n \sum (t-\bar t)^2} } \end{aligned} $$

Istotność parametrów strukturalnych ($H_0: \beta=0$). Statystyka: $$T_{n-2} = \frac{\hat \beta}{ s(\beta) } $$ ma rozkład $t$-Studenta z $n-2$ stopniami swobody.

Współczynnik zbieżności: $$ \phi^2 = \frac{\sum_{t=1}^n (y_t -\hat y_t)^2 }{ \sum_{t=1}^n (y_t -\bar y_t)^2 } \cdot 100% \quad 0\leq \phi^2 \leq 100 $$ Im $\phi^2$ jest bliższy 0, tym dopasowanie jest lepsze. Wartość $\phi^2$ interpretuje się jako ,,procent zmienności zmiennej $y$ nie objaśniony przez liniową funkcję trendu.''

Współczynnik determinacji $R^2=100-\phi^2$ interpretuje się jako ,,procent zmienności zmiennej $y$ objaśniony przez liniową funkcję trendu.''

Analiza wahań sezonowych

Zakładamy, że szereg dzieli się na $s$ powtórzeń a każde powtórzenie składa się z $k$ faz. Np. w szeregu 12 elementowym zawierającym obserwacje kwartalne $k=4$ a $s=3$.

Działania mające na celu wyodrębnienie wahań sezonowych są następujce:

Wygładzamy szereg czasowy analitycznie lub mechanicznie (średnia ruchoma).

Uwalniamy szereg czasowy od trendu. W tym celu wyliczamy wielkości $w_t = y_t/\hat y_t$

Pozbywamy się wahań przypadkowych w wielkościach $w_t$. W tym celu obliczamy średnie dla okresów jednoimiennych (surowe wskaźniki sezonowości): $$ c_i' = \frac{\sum_{j=1}^{s-1} w_{i+j\cdot k} }{s} \quad i=1, 2, \ldots ,k $$ Interpretacja: $(wskaźnik-surowy -1)\cdot 100%$ oznacza o ile procent poziom zjawiska w danej fazie jest wyższy/niższy od poziomu wyznaczonego przez trend.

Suma wskaźników surowych powinna być równa $k$, tj.  $\sum c_i'=k$. Jeżeli tak nie jest, to należy wskaźniki surowe pomnożyć przez współczynnik korygujący: $$ wk=\frac{k}{\sum c_i'} $$ Otrzymując w ten sposób czyste wskaźnik sezonowości.

Jeżeli pomnożymy w każdym okresie teoretyczny poziom zjawiska $\hat y_t$ przez odpowiedni dla danego okresu wskaźnik sezonowości, to otrzymamy teoretyczny poziom zjawiska uwzględniający wahania sezonowe.