Univerzita Karlova v Praze Matematicko-fyzikální fakulta. Testy normality. Katedra pravděpodobnosti a matematické statistiky

Univerzita Karlova v Praze Matematicko-fyzikální fakulta

BAKALÁŘSKÁ PRÁCE

Lukáš Kotlorz Testy normality Katedra pravděpodobnosti a matematické statistiky

Vedoucí bakalářské práce: prof. RNDr. Jiří Anděl, DrSc. Studijní program: Obecná matematika Studijní obor: Matematika

Praha 2012

Rád bych poděkoval vedoucímu své práce panu prof. RNDr. Jiřímu Andělovi, DrSc. za pomoc, kterou mi poskytl při vypracovávání bakalářské práce, a také za jeho čas, který mi věnoval.

Prohlašuji, že jsem tuto bakalářskou práci vypracoval samostatně a výhradně s použitím citovaných pramenů, literatury a dalších odborných zdrojů. Beru na vědomí, že se na moji práci vztahují práva a povinnosti vyplývající ze zákona č. 121/2000 Sb., autorského zákona v platném znění, zejména skutečnost, že Univerzita Karlova v Praze má právo na uzavření licenční smlouvy o užití této práce jako školního díla podle § 60 odst. 1 autorského zákona.

V . . . . . . . . . . . . . . dne . . . . . . . . . . . . . . .

Podpis autora

Název práce: Testy normality Autor: Lukáš Kotlorz Katedra: Katedra pravděpodobnosti a matematické statistiky Vedoucí bakalářské práce: prof. RNDr. Jiří Anděl, DrSc., Katedra pravděpodobnosti a matematické statistiky Abstrakt: Cílem práce zaměřené na testování normality je popsání nejen statistických testů, ale i grafických metod. V první kapitole práce jsou popsány grafické metody, které se používají při testování normality zejména histogram, boxplot, Q-Q plot. Ve druhé části jsou popsány testy např. Shapirův-Wilkův, Kolmogorovův-Smirnovův, Lillieforsův, Andersonův-Darlingův, chí-kvadrát používané na testování shody rozdělení náhodného výběru s normálním rozdělením. U každého testu je uvedena testová statistika, kritický obor, případně odkaz na literaturu, kde lze najít tabulku kritických hodnot. Ve třetí části je provedena simulační studie, zda náhodný výběr pochází z normálního rozdělení. Výběry z různých rozdělení byly generovány programem R. Klíčová slova: testování normality, normální rozdělení, knihovna nortest

Title: Normality tests Author: Lukáš Kotlorz Department: Department of Probability and Mathematical Statistics Supervisor: prof. RNDr. Jiří Anděl DrSc., Department of Probability and Mathematical Statistics Abstract: The aim of this thesis focused on testing normality is to describe both statistical tests and graphical methods. The first part is devoted to graphical methods used to testing normality (particularly Histogram, Boxplot and Q-Q Plot). The tests used for testing the conformity of random sample distribution with normal distribution, e.g., Shapiro-Wilk, Kolmogorov-Smirnov, Lilliefors, AndersonDarling, Chi-squared, are described in the second part. The test statistics, the critical region and alternatively the link for tabulated critical values are listed for each test. The simulations, whether the random sample comes from normal distribution, are described in the third part. The samples from different distributions were generated by Program R. Keywords: testing normality, normal distribution, library nortest

Obsah Úvod

2

1 Graf ické metody k ověřování normality 1.1 Histogram . . . . . . . . . . . . . . . . . 1.1.1 Konstrukce a základní vlastnosti 1.1.2 Optimální volba počtu tříd . . . . 1.1.3 Histogram v programu R . . . . . 1.2 Stem-and-leaf plot . . . . . . . . . . . . 1.3 Box-and-whiskers plot . . . . . . . . . . 1.4 Quantile-quantile (Q-Q) plot . . . . . . . 1.4.1 Kostrukce a základní vlastnosti . 1.4.2 Q-Q plot v programu R . . . . . 1.5 Percent-percent (P-P) plot . . . . . . . .

. . . . . . . . . .

. . . . . . . . . .

. . . . . . . . . .

. . . . . . . . . .

. . . . . . . . . .

. . . . . . . . . .

. . . . . . . . . .

. . . . . . . . . .

. . . . . . . . . .

. . . . . . . . . .

. . . . . . . . . .

. . . . . . . . . .

. . . . . . . . . .

. . . . . . . . . .

4 4 4 4 5 5 6 7 7 8 8

2 Regresní testy 2.1 Shapirův-Wilkův test . . . . . . . 2.2 Shapirův-Franciův test . . . . . . 2.3 Fillibenův test . . . . . . . . . . . 2.4 Percent-Percent Plot Correlation 2.5 LaBrecqueovy testy normality . .

. . . . .

. . . . .

. . . . .

. . . . .

. . . . .

. . . . .

. . . . .

. . . . .

. . . . .

. . . . .

. . . . .

. . . . .

. . . . .

. . . . .

12 12 12 13 14 14

. . . . . . . . .

16 16 17 17 18 19 19 19 20 20

. . . . . . . .

21 22 22 22 22 23 23 23 24

5 Chí-kvadrát test dobré shody 5.1 Chí-kvadrát test při známých parametrech . . . . . . . . . . . . . 5.2 Chí-kvadrát test při neznámých parametrech . . . . . . . . . . . . 5.3 Počet a šířka tříd . . . . . . . . . . . . . . . . . . . . . . . . . . .

25 25 25 26

. . . . .

. . . . .

. . . . .

. . . . .

3 Testy založené na momentech 3.1 Základní test šikmosti a špičatosti . . . . . . . . . . . 3.2 D’Agostinův test šikmosti . . . . . . . . . . . . . . . 3.3 Anscombeho-Glynnův test špičatosti . . . . . . . . . 3.4 Jarqueho-Berův test . . . . . . . . . . . . . . . . . . 3.5 Bowmanův-Shentonův test šikmosti a špičatosti . . . 3.6 Obdélníkový test šikmosti a špičatosti . . . . . . . . . 3.7 Testy založené na absolutních centrálních momentech 3.8 Gearyho test . . . . . . . . . . . . . . . . . . . . . . . 3.9 Bonettův-Seierové test . . . . . . . . . . . . . . . . .

. . . . . . . . .

4 Testy založené na empirických distribučních funkcích 4.1 Kolmogorovův-Smirnovův test a Lillieforsova varianta . 4.1.1 Kolmogorovův-Smirnovův test . . . . . . . . . . 4.1.2 Lillieforsova varianta . . . . . . . . . . . . . . . 4.2 Kuiperův test . . . . . . . . . . . . . . . . . . . . . . . 4.3 Cramérův-von Misesův test . . . . . . . . . . . . . . . 4.4 Watsonův test . . . . . . . . . . . . . . . . . . . . . . . 4.5 Andersonův-Darlingův test . . . . . . . . . . . . . . . . 4.6 Stephensovy modifikace . . . . . . . . . . . . . . . . .

1

. . . . . . . . .

. . . . . . . .

. . . . . . . . .

. . . . . . . .

. . . . . . . . .

. . . . . . . .

. . . . . . . . .

. . . . . . . .

. . . . . . . . .

. . . . . . . .

5.4

Yarnoldovo pravidlo . . . . . . . . . . . . . . . . . . . . . . . . . .

6 Další testy normality 6.1 Test založený na rozpětí 6.2 Uthoffův test . . . . . . 6.3 D’Agostinův test . . . . 6.4 Ojovy testy . . . . . . . 6.5 Vašíčkův test . . . . . .

. . . . .

. . . . .

. . . . .

. . . . .

. . . . .

. . . . .

. . . . .

. . . . .

. . . . .

. . . . .

. . . . .

. . . . .

. . . . .

. . . . .

. . . . .

. . . . .

. . . . .

. . . . .

. . . . .

. . . . .

. . . . .

. . . . .

. . . . .

27 28 28 28 28 29 29

7 Testy na odlehlá pozorování 7.1 Grubbsovy testy na odlehlá pozorování . . . . . . . . . . . . . . . 7.1.1 Test na jedno odlehlé pozorování . . . . . . . . . . . . . . 7.1.2 Test na více odlehlých pozorování na konkrétním chvostu . 7.1.3 Test na jedno odlehlé pozorování na obou chvostech . . . . 7.1.4 Grubbsovy testy v programu R . . . . . . . . . . . . . . . 7.2 Dixonův test na jedno odlehlé pozorování . . . . . . . . . . . . . . 7.3 Tietjenův-Mooreův test na více odlehlých pozorování na jednom nebo obou chvostech . . . . . . . . . . . . . . . . . . . . . . . . .

31 32 32 33 34 34 34

8 Porovnání některých testů normality na testovacích 8.1 Normální rozdělení . . . . . . . . . . . . . . . . . . . 8.2 Exponenciální rozdělení . . . . . . . . . . . . . . . . 8.3 Cauchyho rozdělení . . . . . . . . . . . . . . . . . . . 8.4 Laplaceovo rozdělení . . . . . . . . . . . . . . . . . . 8.5 Rovnoměrné rozdělení . . . . . . . . . . . . . . . . . 8.6 Kontaminované rozdělení . . . . . . . . . . . . . . . . 8.7 Shrnutí . . . . . . . . . . . . . . . . . . . . . . . . . . 8.8 Komentář k Softwarovému kódu . . . . . . . . . . . .

36 37 38 38 38 40 40 41 42

datech . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

. . . . . . . .

. . . . . . . .

35

Závěr

43

Seznam použité literatury

44

2

Úvod Hodně statistických testů je založeno na předpokladu, že náhodný výběr pochází z normálního rozdělení. Některé testy jsou na tento předpoklad velmi citlivé, v jiných případech mírné odchylky od normality ovlivní výsledek jen nepatrně. Bohužel žádný test, který by určil, že daný výběr pochází z normálního rozdělení neexistuje, ale již byla publikována spousta testů, které mohou prokázat, že výběr z normálního rozdělení nepochází. V bakalářské práci jsou představeny nejen testy, které se dnes běžně používají, ale i testy méně známé. V kapitole 1 jsou představeny grafické metody, které nám mohou před samotným testováním prozradit, co od jednotlivých testů očekávat. V kapitole 2 se věnujeme regresním testům, ve 3. kapitole jsou představeny testy založené na momentech. Kapitola 4 je věnována testům založeným na empirické distribuční funkci, v 5. kapitole popíšeme χ2 test dobré shody. V kapitole 6 jsou představeny další testy normality, jedná se o testy, které nebylo rozumné zařadit do některé z předchozích kapitol. V kapitole 7 se zabýváme testy na odlehlá pozorování, protože v 8. kapitole, ve které porovnáme některé testy na datech, zjistíme, že přítomnost odlehlého pozorování může mít velký vliv na výsledek. Testy byly rozděleny do jednotlivých kapitol stejně, jak to ve své knize učinil Thode [1]. U každého testu je uvedena testová statistika, informace, kdy zamítáme nulovou hypotézu, zdroj, kde najdeme tabulku kritických hodnot, pokud testovou statistku neporovnáváme s kvantily známých rozdělení, a zdroj, kde lze o daném testu najít další informace. U testů, které jsou naprogramovány ve statistickém programu R, je uveden i balíček, kde daný test najdeme, a příkaz, kterým daný test dostaneme, přičemž c představuje proměnnou, ve které máme uložen náhodný výběr, a tři tečky označují možnost přidání volitelných parametů, jejichž přehled lze u příslušného příkazu najít v dokumentaci k programu R [3]. V práci budeme předpokládat, že máme náhodný výběr x1 , x2 , ..., xn , kde n je rozsah výběru. Všechny testy normality jsou založené na nulové hypotéze, že náhodný výběr pochází z rozdělení N(µ, σ 2 ) oproti alternativě, že pochází z jiného rozdělení. Testovat budeme na hladině α, přičemž nejčastěji volíme α = 0,05. Tato testová hladina je použita i u simulační studie. V tabulce 1 je uveden přehled značení, která budu v práci používat a v tabulce 2 přehled hustot a značení rozdělení, která jsou v práci použita. Kritické hodnoty u(α) a χ2k (α) jsou definovány jako čísla, která náhodná veličina s rozdělením N(0, 1) resp. χ2k překročí s pravděpodobností α.

3

Tabulka 1: Značení. Značení x(i) x

Význam i-tý prvek po uspořádání náhodného výběru n X výběrový průměr, tj. n1 xi i=1

µ µ b σ2 a σ ≥ 0 σ b α3 α4 Φ Φ−1 (α) u(α) χ2k (α) dae bac

střední hodnota odhad střední hodnoty rozptyl a směrodatná odchylka odhad směrodatné odchylky E(X − EX)3 šikmost, tj. σ3 E(X − EX)4 špičatost, tj. σ4 distribuční funkce standardizovaného normálního rozdělení α-kvantil standardizovaného normálního rozdělení kritická hodnota standardizovaného normálního rozdělení kritická hodnota χ2 rozdělení o k stupních volnosti horní celá část a dolní celá část a

Tabulka 2: Hustoty a značení rozdělení. Rozdělení

Hustota

Značení

Normální

  (x − µ)2 √ exp − 2σ 2 2πσ 2

N(µ, σ 2 )

Exponenciální Laplaceovo

1

n xo 1 exp − I(x > 0) λ λ   1 |x − a| exp − 2b b

Exp(λ) DEx(a, b)

Rovnoměrné

1 I(a < x < b) b−a

Cauchyho

1 1 π 1 + x2

Chí kvadrát

1 xk/2−1 e−k/2 I(x > 0) χ2k k k/2 2 Γ( 2 )

4

R(a, b)

1. Graf ické metody k ověřování normality Grafické metody patří mezi důležité způsoby prezentace statistických údajů. Využívají se nejen při testování normality, ale např. i v regresní analýze a testech dobré shody. Někdy jsou i vyžadovány jako součást statistické analýzy. Grafické metody sice nejsou tak přesné jako testové, ale mohou nám poskytnout určitou informaci o rozdělení dat. Bohužel neposkytují objektivní důvod k zamítnutí nebo nezamítnutí H0 , neboť z grafů nevyčteme, je-li již testová statistika signifikatní na námi zvolené hladině α. U velkých výběrů nám na oplátku mohou prozradit, zda zamítnutí H0 bylo způsobeno pouze drobnou odchylkou od normality, nebo byla-li H0 zamítnuta z důvodu, že náš výběr evidentně nepochází z normálního rozdělení, a mohou pomoci určit, z jakého rozdělení může výběr pocházet. Grafické metody na rozdíl od některých testů nekladou omezení na rozsah výběru. Dnes již je většina grafických metod implementována v každém statistickém softwaru, včetně programu R.

1.1 1.1.1

Histogram Konstrukce a základní vlastnosti

V minulosti byl histogram nejoblíbenější grafickou metodou při testování normality, neboť jej bylo snadné sestrojit ručně. Ve statistice nám umožňuje vizuálně porovnat tvar hustoty četností s tvarem hustoty námi zvoleného rozdělení. Při konstrukci nejprve určíme rozpětí naměřených hodnot R = x(n) − x(1) . Nechť k je počet tříd a h šířka tříd. Možnosti pro optimální volbu počtu tříd najdeme v sekci 1.1.2. Po zvolení hodnot k a h provedeme úpravu tak, aby hraniční body jednotlivých tříd byla zaokrouhlená čísla. Zastupitelné hodnoty zi spočteme jako aritmetický průměr horní a dolní hranice každé třídy. Poté zaznamenáme do tabulky četností počet výskytů pozorování v jednotlivých třídách. Při konstrukci vyneseme na vodorovnou osu hodnoty zi a na svislou osu hodnoty četností v jednotlivých třídách. Na obrázcích 1.5 - 1.8 postupně vidíme po šesti náhodných výběrech z rozdělení N(50, 200) postupně o rozsahu n = 50, 100, 150, 200. U výběrů o velikosti n = 50 by náš odhad, že na histogramu 1 se jedná o výběr z rozdělení s kladnou špičatostí a u histogramů 3 a 5 o rozdělení se zápornou špičatostí, nemusel být chybný. U výběrů o rozsahu n = 100 lze pozorovat větší odchylku jen u histogramu 10, kde máme málo pozorování v intervalu [50, 60]. U výběrů s větším rozsahem se odchylky od normality ještě více zmenšují, u histogramů 19 a 22 můžeme říci, že shoda je již skoro dokonalá.

1.1.2

Optimální volba počtu tříd

Pokud z nějakého důvodu nechceme sami rozhodnout o počtu tříd, doporučuje se použít některé z níže uvedených pravidel: • Sturgesovo pravidlo: Pravidlo k = d1 + log2 ne bylo publikováno v práci Sturges [6]. Ten uvažoval histogram o k třídách, kde počet pozorování v jed5


notlivých třídách je dán binomickým koeficientem k−1 , i = 0, 1, ..., k − 1. i S rostoucím k by se tvar tohoto ideálního histogramu blížil hustotě normálního rozdělení. Počet prvků n dostaneme jako n=

 k−1  X k−1 i=0

i

= (1 + 1)k−1 = 2k−1 .

(1.1)

Vyjádříme-li ze vztahu (1.1) k, dostaneme Sturgesovo pravidlo. • Scottovo pravidlo: h = 3,5 σ bn−1/3 . Bylo publikováno v práci Scott [7]. Je vhodné pro výběry z normálního rozdělení, neboť minimalizuje střední čtvercovou chybu odhadu hustoty. −Q25 ) , kde Qi je i% kvan• Freedmanovo-Diaconisovo pravidlo: h = 2 (Q75√ 3n til. Na Wikipedii [52] lze najít odkaz na práci Freedman a Diaconis [8], kde bylo pravidlo publikováno. Je podobné Scottovu pravidlu, ale je méně citlivé na odlehlá pozorování. √ • Odmocninové kritérium: k = n. Jedná se o triviální odhad k, nicméně není vhodný, neboť pro výběr větší než 2500 dává více než 50 tříd. Lze jej nalézt na Wikipedii [52] a používá jej program MS Excel.

Počet nebo velikost tříd dále dopočteme dle vzorce R = k h.

1.1.3

Histogram v programu R

V programu R histogram dostaneme příkazem hist(c,...). počet tříd použito Sturgesovo pravidlo, ale volbou parametru žádat jiné. Použijeme breaks="Scott" pro Scottovo pravidlo, Friedmanovo-Diaconisovo, breaks=k pro námi zvolený počet program R upraví tak, aby h ∈ M, kde   m 10 , m ∈ N , j ∈ {1, 2, 5} , M= j

Implicitně je pro breaks si lze vybreaks="FD" pro tříd, který si ale

(1.2)

případně breaks=d, kde d je vektor dělících bodů s tím, že je nutné jej zadat tak, aby platily nerovnosti d1 ≤ x(1) a x(n) ≤ dm .

1.2

Stem-and-leaf plot

Stem-and-leaf plot obdobně jako histogram umožňuje vizuálně porovnat tvar hustoty četností s tvarem hustoty námi zvoleného rozdělení, přičemž dává nějaké informace navíc. Při konstrukci nejprve uspořádáme náš výběr, provedeme rozumné zaokrouhlení dat a určíme hodnoty, které vyneseme na stonek (stem). Šířka tříd je omezená na množinu M (viz (1.2)), přičemž za j volíme nejčastěji hodnotu 1. Hranice jednotlivých tříd jsou určeny šířkou tak, abychom mohli snadno vynést hodnoty na stonek. Na listy (leaves) vyneseme následující číslici z dekadického zápisu všech pozorování k příslušné hodnotě stonku. Postupujeme od prvku x(1) k prvku x(n) , tedy na každém listu budou nižší hodnoty blíže stonku. Stem-andleaf plot ukazuje konkrétní hodnoty jednotlivých pozorování, proto lze snadno 6

určit medián, extrémní hodnoty, kvartily, případně i jiné percentily. Prvek, který odpovídá i-tému percentilu, vypočítáme dle vzorce i(n + 1)/100. V programu R stem-and-leaf plot použijeme příkazem stem(c,...). Nevýhodou je, že výstupem není obrázek, ale text. Na obrázku 1.1 vidíme příklad stem-and-leaf plotu.

Obrázek 1.1: Stem-and-leaf plot náhodného výběru z N(50, 200), rozsah n = 100

1.3

Box-and-whiskers plot

Často je používán zkrácený název boxplot, případně krabičkový graf. Ukazuje jiné informace než histogram nebo stem-and-leaf plot. Boxplot poskytuje informace o symetrii a rozložení dat a umožňuje snadno odhalit odlehlá pozorování. Příklad vidíme na obrázku 1.2. Hodnota mediánu M je na obrázku znázorněna tlustou čarou uprostřed. Hodnoty D a H dolního a horního  závěsu, znázorněny n e , jako konce „krabiceÿ, se vypočtou ze vztahů D = med x ; i = 1, 2, ..., d (i) 2  n n H = med x(i) , i = b 2 c; b 2 c + 1, ..., n . Dolní vnitřní hradby jsou dány nejmenší hodnotou j splňující x(j) ≥ D − 1.5(H − D) pro dolní, resp. největší hodnotou l splňující x(l) ≤ H − 1.5(H − D) pro horní. Na obrázku jsou znázorněny jako konce čárkovaných úseček. Kolečkem jsou znázorněná odlehlá pozorování, jedná se o všechny hodnoty ležící za vnitřními hradbami. Pokud pochází náš výběr z normálního rozdělení, pak by měl být graf symetrický, protože hustota normálního rozdělení je symetrická. V programu R se boxplot vykreslí použitím příkazu boxplot(c,...). Do jednoho obrázku lze vykreslit histogram i boxplot příkazem simple.hist.and.boxplot(c), který najdeme v balíčku UsingR. ●

4 3 2 1 0

Obrázek 1.2: Boxplot náhodného výběru z Exp(1) o rozsahu n = 40

7

1.4

Quantile-quantile (Q-Q) plot

1.4.1

Kostrukce a základní vlastnosti

Tato grafická metoda porovnává kvantily uspořádaného výběru s teoretickými kvantily předpokládaného rozdělení, v našem případě normálního. Pro konstrukci budeme potřebovat uspořádaný výběr. Hodnota x(i) představuje qi % kvantil, kde qi = 100(i − 0,5)/n. Dále potřebujeme hodnoty pi , ve kterých vypočteme kvantily i−β , kde 0 ≤ β < 1 je normálního rozdělení. Ty se vypočtou ze vztahu pi = n+1−2β 1 1 korigující faktor. Nejčastěji za parametr β volíme hodnoty 2 , 3 nebo 38 . Graf, kde za β volíme hodnotu 0,317 5 bývá někdy označován jako Normal-Probability plot. Další volby parametru β lze nalézt v knize Thode [1], tabulka 2.1, str. 20. Poté do
grafu vyneseme dvojice Φ−1 (pi ), x(i) . Pokud výběr pochází z normálního rozdělení, měly by body ležet přibližně na přímce, jak vidíme na obrázku 1.3a. U Q-Q plotu se dá v porovnání s normálním rozdělením odhadnout šikmost a špičatost rozdělení našeho výběru. Při výběru z rozdělení s těžkými chvosty (α4 > 3) body vytvoří křivku odpovídající funkci tg x (obrázek 1.3b), pro rozdělení s lehkými chvosty (α4 < 3) bude křivka odpovídat funkci arctg x (obrázek 1.3c). Při výběru z rozdělení s kladnou šikmostí (α3 > 0) body vytvoří kovexní křivku (obrázek 1.3d), rozdělení se zápornou šikmostí (α3 < 0) odpovídá konkávní křivka (obrázek 1.3e). Na obrázku 1.3f vidíme Q-Q plot výběru z N(0, 1) se dvěma odlehlými pozorováními, kde je možnost záměny za výběr z rozdělení s těžkými chvosty, ale na rozdíl od obrázku 1.3b leží ostatní pozorování na přímce. b

c 4

●

4

4

a

●

2 0

2 0

●

−2

−2

0

●● ●● ●●●● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ●●●●● ● ●●

−2

2

● ● ● ●● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ●●● ● ● ● ●

●

● ●●●●●●● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ●●●●●●●

●

−2

−1

0

1

2

−4

−4

−4

●

−2

−1

0

2

−2

−1

0

1

2

4

f

4

e

4

d

1

●

0

1

2

●

−2

−1

0

1

2

●●● ●●●● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ●●●● ●●●

0

●

−2

−2 −4 −1

●●●●●●● ● ● ● ● ● ●● ● ●● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ●● ● ● ●

−4

−2

−2

2

2

●

−4

●

● ● ●●●●● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ●● ● ●●●●●●●

0

0

2

● ●

●

●

−2

−1

0

1

2

Obrázek 1.3: Q-Q plot náhodného výběru z (a) N(0, 1), (b) DEx(0, 1), (c) R(−2, 2), (d) Exp(1) − 1, (e) 1 − Exp(1) (f) N(0, 1) se dvěma odlehlými pozorováními

8

1.4.2

Q-Q plot v programu R

V programu R dostaneme Q-Q plot příkazem qqnorm(c,...). Zde je standardně pro výpočet hodnot pi použita pro výběr o velikosti n ≤ 10 hodnota β = 3 , při větším výběru β = 12 . Chceme-li jinou hodnotu, musíme použít příkaz 8 qqplot(d,c,...), kde d je vektor kvantilů. Ty lze vypočítat pomocí příkazu qnorm(ppoints(n,β)). Chceme-li zkoumat shodu s jiným rozdělením, použijeme z knihovny car příkaz qqPlot(x, distribution="name"), kde za name dosadíme námi požadované rozdělení.

1.5

Percent-percent (P-P) plot

Můžeme jej nalézt rovněž pod názvem Probability-probability plot. Stejně jako Q-Q plot se používá k porovnání kvantilů výběru a teoretického rozdělení. Opět budeme potřebovat pořádkové statistiky x(1) ≤ x(2) ≤ ... ≤ x(n) . Rozdíl je v tom, 

že u P-P plotu se do grafu vynášejí body (pi , Zi ), kde Zi = Φ

x(i) −b µ σ b

, proto

0.0

0.2

0.4

0.6

0.8

1.0

i všechny body leží ve čtverci [0 ; 1] × [0 ; 1]. Zde se používají hodnoty pi = n+1 . Pokud náš výběr pochází z normálního rozdělení, měly by body vytvořit úsečku s krajními body [0 ; 0] a [1 ; 1], takže není nutné hledat nejlepší přímku, která dané body aproximuje. Šikmost a špičatost rozdělení našeho výběru poznáme podle stejných pravidel jako u Q-Q plotu. Je-li náš výběr ze symetrického rozdělení a rozdělení, se kterým jej porovnáváme, je rovněž symetrické, pak křivka vytvořená těmito body bude procházet bodem [0,5 ; 0,5]. V programu R žádná funkce na vykreslení P-P plotu není implementována, ale lze jej vytvořit pomocí příkazu plot(ppoints(n,0),pnorm((sort(c) - mean(c))/sd(c))). Příklad P-P plotu vidíme na obrázku 1.4. Chceme-li porovnávat shodu s jiným rozdělením, nahradíme příkaz pnorm, který počítá hodnoty distribuční funkce z normálního rozdělení, příkazem počítajícím hodnoty distribuční funkce příslušného rozdělení.

● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ●

0.0

0.2

0.4

0.6

0.8

1.0

Obrázek 1.4: P-P plot náhodného výběru z N(0, 1) o rozsahu n = 125

9

20

40

60

15 0

5

10

15 10 5 0

0

5

10

15

20

Histogram 3

20

Histogram 2

20

Histogram 1

80

20

60

80

60

60

80

20 0

5

10

15

20 15 5 0 40

40

Histogram 6

10

15 10 5 0

20

20

Histogram 5

20

Histogram 4

40

80

20

40

60

80

20

40

60

80

Obrázek 1.5: Histogramy šesti náhodných výběrů z N(50, 200) o rozsahu n = 50

0

20

40

60

80 100

10 15 20 25 30 5 0

20

40

60

80 100

40

60

80 100

20

40

60

80 100

10 15 20 25 30

Histogram 12

0

5

10 15 20 25 30 5 20

0

Histogram 11

0

0

5

10 15 20 25 30

Histogram 10

0

Histogram 9

0

5

10 15 20 25 30

Histogram 8

0

0

5

10 15 20 25 30

Histogram 7

0

20

40

60

80 100

0

20

40

60

80 100

Obrázek 1.6: Histogramy šesti náhodných výběrů z N(50, 200) o rozsahu n = 100

10

0

20

40

60

80 100

50 40 30 20 10 0

20

40

60

80 100

40

60

80 100

20

40

60

80 100

50

Histogram 18

0

10

20

30

40

50 40 30 20 10 20

0

Histogram 17

0

0

10

20

30

40

50

Histogram 16

0

Histogram 15

0

10

20

30

40

50

Histogram 14

0

0

10

20

30

40

50

Histogram 13

0

20

40

60

80 100

0

20

40

60

80 100

Obrázek 1.7: Histogramy šesti náhodných výběrů z N(50, 200) o rozsahu n = 150

0

20

40

60

80 100

10 20 30 40 50 60 0

20

40

60

80 100

40

60

80 100

20

40

60

80 100

10 20 30 40 50 60

Histogram 24

0

10 20 30 40 50 60 20

0

Histogram 23

0

0

10 20 30 40 50 60

Histogram 22

0

Histogram 21

0

10 20 30 40 50 60

Histogram 20

0

0

10 20 30 40 50 60

Histogram 19

0

20

40

60

80 100

0

20

40

60

80 100

Obrázek 1.8: Histogramy šesti náhodných výběrů z N(50, 200) o rozsahu n = 200

11

2. Regresní testy I když nám grafické metody mohou velmi pomoci při určení charakteristiky rozdělení našeho výběru, přesto potřebujeme metody, které mohou poskytnout
infor−1 maci, můžeme-li zamítnout H0 . Předpokládaný vztah mezi Φ (pi ), x(i) v Q-Q plotu (viz sekce 1.4) nás svádí k použití lineární regrese. Testy by šlo založit na vzdálenosti jednotlivých bodů od regresní přímky.

2.1

Shapirův-Wilkův test

Tento test byl publikován v práci Shapiro a Wilk [9]. Vzhledem k dobré vypovídající hodnotě se jedná o nejpoužívanější test pro malé výběry, dokonce pro výběry o velikosti n ≤ 50 je předepisován i normou ČSN 01 0225 [5]. Nechť y(1) ≤ y(2) ≤ ... ≤ y(n) je uspořádaný náhodný výběr z rozdělení N(0, 1), m> = (m1 , m2 , ..., mn ) = E(y(i) ) je vektor středních hodnot i-té pořádkové statistiky standardizovaného normálního rozdělení a V = (vij ) = cov(y(i) , y(j) ) je kovarianční matice pořádkových statistik. Nechť pozorované hodnoty x(i) pocházejí z rozdělení N(µ, σ 2 ), potom x(i) = µ + σy(i) . Nejlepší nestranný odhad σ získaný metodou nejmenších čtverců je roven b = a> x, kde a= √

m> V−1 m> V−1 V−1 m

, a> a = 1, x = (x(1) , x(2) , ..., x(n) )> .

Pro vektor a platí an = −an+1−i , a navíc při lichém n je a n+1 = 0. Proto lze 2 b spočítat ze vzorce bn c 2 X b= an−i+1 (x(n−i+1) − x(i) ). i=1

Pro zjednodušení výpočtu byly hodnoty ai pro n ≤ 50 tabelovány a lze je najít v [9], tabulka 5, str. 603–604. Testová statistika je definovaná vztahem W =

b2 n X

.

(xi − x)2

i=1

Hypotézu H0 zamítáme, pokud je hodnota W menší než tabelovaná kritická hodnota. Tabulku kritických hodnot najdeme v [9], tabulka 6, str. 605. Problémem tohoto testu je, že prvky matice V jsou přesně známy pouze pro výběry o rozsahu ne větším než 20. Pro výběry od 21 do 50 prvků byly uvedeny odhady ai v [9], str. 596–600. V práci Royston [10] byla pro výběry o rozsahu 7 ≤ n ≤ 2000 prezentována transformace testové statistiky W . Po transformaci bude mít W normální rozdělení. V programu R dostaneme Shapirův-Wilkův test pro výběry o rozsahu 3–5 000 příkazem shapiro.test(c).

2.2

Shapirův-Franciův test

Tento test byl publikován v práci Shapiro a Francia [11] jako modifikace ShapirovaWilkova testu pro výběry o rozsahu n > 50. Při takto velkém výběru mohou být 12

jednotlivá pozorování považována za nezávislá, proto v tomto případě můžeme použít substituci V−1 = I jak je navrženo v práci Gupta [12]. Po substituci dostaneme testovou statistiku m (d> x)2 , d= √ , W0 = n X m> m 2 (xi − x) i=1

kde m je opět vektor středních hodnot i-tých pořádkových statistik standardizovaného normálního rozdělení. Hypotézu H0 zamítáme, pokud je hodnota W 0 menší než tabelovaná kritická hodnota (viz Shapiro, Francia [11]). Tento test lze v programu R najít v balíčku nortest pod příkazem sf.test(c) a lze jej použít pro výběry o rozsahu 5–5 000 prvků.

2.3

Fillibenův test

Tento test byl publikován v práci Filliben [13]. Nechť E(x(i) ) jsou střední hodnoty pořádkových statistik a Mi = med(x(i) ) mediány pořádkových statistik. Zde uvažujme mediány ze standardizovaného normálního rozdělení. Filliben navrhl použít místo E(x(i) ) hodnoty Mi , neboť spočítat Mi je snazší než spočítat hodnoty E(x(i) ). Navíc hodnoty Mi lze spočítat i z rozdělení, které nemají střední hodnotu definovanou, např. Cauchyho rozdělení. Máme-li hodnoty Mi spočteny,
jsme schopni sestrojit P-P plot (viz sekce 1.5), kde vyneseme hodnoty x(i) , Mi . Testová statistika r je definována jako korelační koeficient mezi Mi a x(i) : n X (x(i) − x)(Mi − M )

. r = cor(Mi , x(i) ) = s ni=1 n X X (x(i) − x)2 (Mi − M )2 i=1

i=1

Vzhledem k tomu, že u standardizovaného normálního rozdělení platí M = 0 a Mi = −Mn−i+1 , lze testovou statistiku r zjednodušit do tvaru n X

r=s n X

x(i) Mi

i=1

(x(i) − x)2

i=1

n X

. Mi2

i=1

Nyní nám zbývá určit hodnoty Mi z N(0, 1). Nechť mi jsou mediány pořákových statistik z R[0, 1]. Aplikováním transformace Mi = Φ−1 (mi ) dostaneme hodnoty Mi . Filliben doporučuje mi spočítat ze vztahu   1 − mn pro i = 1, i−0,317 5 pro 1 < i < n, mi =  n+0,365 1/n 0,5 pro i = n, které se od skutečných hodnot liší nejvýše v řádu desetitisícín (viz Filliben [13], str. 115–116). Hypotézu H0 zamítáme, pokud je hodnota r menší než tabelovaná kritická hodnota. Tabulku kritických hodnot lze nalézt v práci Filliben [13], tabulka 1, str. 113. 13

2.4

Percent-Percent Plot Correlation

Tyto dva testy byly publikovány v práci Gan a Koehler [14]. Jsou založeny na rozdíl od předchozích testů na vzdálenosti bodů od regresní přímky v P-P plotu (viz sekce 1.5). Hodnoty pi vypočteme, obdobně jako u P-P plotu, ze vztahu i . Testové statistiky mají tvar pi = n+1 " n #2 X (Z(i) − Z)(pi − p)   2 , k 2 = cor(Zi , p(i) ) = n i=1 n X X 2 2 (pi − p) (Z(i) − Z) i=1

i=1

" k02 =

n X

#2 (Z(i) − 0,5)(pi − 0,5)

i=1 n X

n X (Z(i) − 0,5) (pi − 0,5)2

i=1

kde Zi = Φ



x(i) −b µ σ b



,

2

i=1

. Pokud se hodnota Z blíží k hodnotě 0.5, což by mělo nastat

pro symetrická rozdělení, mají testové statistiky k 2 a k02 podobné vlastnosti. Statistika 1 − k 2 dává podíl chyb při použití metody nejmenších čtverců mezi body (pi , Zi ) a regresní přímky oproti chybám při použití metody nejmenších čtverců mezi body (pi , Zi ) a Z. Statistika 1 − k02 dává odpovídající podíl v případě, že spojnice bodů v P-P plotu musí procházet bodem [0,5 ; 0,5]. Hypotézu H0 zamítneme, pokud je hodnota k 2 nebo k02 menší než kritická hodnota, která se vypočte pro příslušnou testovou hladinu α ze vztahu kα2 = 1 −

1 , aα + nbα

kde koeficienty aα a bα pro normální, exponenciální a Gumbelovo rozdělení lze nalézt v práci Gan a Kohler [14], tabulka 2, str. 294.

2.5

LaBrecqueovy testy normality

Testové statistiky F1 , F2 a F3 byly publikovány v práci LaBrecque [15]. Nechť x1 , x2 , ..., xn je uspořádaný náhodný výběr. Je-li z normálního rozdělení, potom E(xi ) = µ + σξi ,

(2.1)

kde ξi = E(yi ) a yi = xiσ−µ je výběr ze standardnizovaného normálního rozdělení. Pokud by data nebyla z normálního rozdelění potom E(xi ) = µ + σξ + ζφ2 (ξi ) + ηφ3 (ξi ),

(2.2)

kde φj (ξi ) je polynom j-tého stupně proměnné ξi . Pokud je model (2.1) správný, potom ζ = η = 0. Pokud není, potom rozšířený model (2.2) je lepší než model (2.1). Dále mějme ještě modely E(xi ) = µ + σξ + ζφ2 (ξi ), E(xi ) = µ + σξ + ηφ3 (ξi ). 14

(2.3) (2.4)

Model (2.3) je vhodný pro nesymetrická rozdělení, kde je Q-Q plot ve tvaru jako na obrázku 1.3d a 1.3e, model (2.4) je vhodný pro rozdělení s těžkými nebo lehkými chvosty, kde je Q-Q plot ve tvaru jako na obrázku 1.3b a 1.3c. Jejich kombinací vznikne model (2.2). Testové statistiky mají tvar ( cζζ )2 + ( cηbη )2 b

F1 =

2s2 ( cζζ )2

,

b

F2 = F3 =

s2 ( cηbη )2 s2

, .

Obdobně jako u Shapirova-Wilkova testu jsou tyto testové statistiky lineární kombinací pořádkových statistik. Navíc ai = an−i+1 , bi = bn−i+1 a pro liché n platí bb n+1 c = 0, takže 2

      

n 2 X

ai (xi − xn−i+1 ) pro n sudé, ζb i=1 = bn c 2  cζ X    ai (xi − xn−i+1 ) + a n+1 x n+1 pro n liché,   2 2 i=1

bnc

2 X ηb = bi (xi − xn−i+1 ). cη i=1

Hodnoty ai a bi pro 4 ≤ n ≤ 64 jsou tabelovány v LaBrecque [15], tabulka 1, str. 294–299, jejich odvození najdeme na str. 295-297. Hypotézu H0 zamítneme, pokud je hodnota odpovídající testové statistiky větší než kritická hodnota, kterou najdeme pro 4 ≤ n ≤ 12 √ v tabulce 2, str. 301 a pro hodnoty n > 12 je lze vypočítat ze vztahu a0 + a1 n, kde hodnoty a0 a a1 jsou v tabulce 3, str. 302.

15

3. Testy založené na momentech Jednou z možností, jak lze od sebe odlišit dvě různá rozdělení, je porovnání jejich centrálních momentů. Centrální moment k-tého řádu je definován vztahem µk = E(X − EX)k , existuje-li konečná střední hodnota EX. Dále definujme šikmost α3 = µ3 /σ 3 a špičatost α4 = µ4 /σ 4 . Pro normální rozdělení platí α3 = 0 a α4 = 3. Bohužel normální rozdělení není jediné, které jednu z rovností splňuje, protože např. liché momenty všech symetrických rozdělení jsou rovny 0. Navíc existují rozdělení, která nejsou symetrická, ale všechny liché momenty jsou rovny 0 (viz Ord [16]). Pro většinu symetrických nenormálních1 rozdělení platí α4 6= 3, ale existují i symetrická nenormální rozdělení s α4 = 3 (viz Johnson, Tietjen a Beckman [17]; Balanda a MacGillivrayová [18]).

3.1

Základní test šikmosti a špičatosti

Pro k ≥ 2 je výběrový centrální moment k-tého řádu z náhodného výběru x1 , x2 , ..., xn definován vztahem n

1X (xi − x)k . mk = n i=1 Pro test šikmosti použijeme testovou statistiku p m3 b1 = p 3 , m2 √ přičemž H0 zamítneme, pokud | b1 | je větší než tabelovaná kritická hodnota.2 Tabulku kritických hodnot najdeme v knize Thode √ √ [1], tabulka B5, str. 381. Je-li b1 > 0, pak jsou data zešikmena doprava, při b1 < 0 jsou zešikmena doleva. Je možné použít i jednostrannou alternativu testu. Test špičatosti je určen testovou statistikou b2 =

m4 , m22

kde H0 zamítneme, pokud b2 neleží v intervalu, který je dán horní a dolní kritickou hodnotou.3 Tabulku dolních a horních kritických hodnot najdeme v knize Thode [1], tabulka B6, str. 382. Tyto testy byly jedny z nejstarších a zároveň nejsilnějších testů normality. Pokud je ovšem náš výběr z rozdělení, pro které rovněž platí α3 = 0 nebo α4 = 3 (příklady viz 3), potom příslušný test bude pro tento výběr velmi slabý. Řešením by mohlo být použít oba testy, ale dohromady nebudou testy provedeny na hladině α, 1

Nenormálním rozdělením je myšleno rozdělení, které je jiné než normální. √ I když by se dle značení dalo usoudit, že hodnota testové statistiky b1 je vždy nezáporná, není tomu tak. Toto značení je pro test šikmosti běžně používané, proto jsem se rozhodl, že jej také použiji. 3 Interval daný horní a dolní kritickou hodnotou není symetrický kolem hodnoty 3, i když by se to v porovnání s testem šikmosti mohlo chybně očekávat. 2

16

ale na hladině větší. Pokud chceme zároveň provést test na šikmost a špičatost na hladině α, použijeme buď Jarqueho-Berův test (viz 3.4), Bowmanův-Shentonův test (viz 3.5) nebo obdélníkový test (viz 3.6). Dnes se při samotném testování √ šikmosti nebo špičatosti používají transformace, při kterých se b1 a b2 převedou do standardizovaného normálního rozdělení. Tyto transformace, které používá i program R, jsou popsány u následujících testů.

3.2

D’Agostinův test šikmosti

√ Za předpokladu nulové hypotézy má b1 asymptoticky N(0, n6 ). V pracích Pearson √ [19] a [20] byly pro konečné výběry vypočteny s přesností n−3 momenty b1 rozdělení. Ty jsou dány vztahy: 6(n − 2) , (3.1) (n + 1)(n + 3) p 36(n − 7)(n2 + 2n − 5) α4 ( b1 ) = 3 + . (3.2) (n − 2)(n + 5)(n + 7)(n + 9) √ √ √ Ze symetrie rozdělení b1 ihned plyne µ( b1 ) = 0 a α3 ( b1 ) = 0. V práci √ D’Agostino [21] byla za platnosti H0 pro n ≥ 8 představena transformace b1 , po které dostaneme testovou statistiku, která má asymptoticky N(0, 1). Nechť √ b1 Y = p √ , σ 2 ( b1 ) q p W2 = 2(α4 ( b1 ) − 1) − 1, p σ 2 ( b1 ) =

√ √ √ kde σ 2 ( b1 ) a α4 ( b1 ) jsou po řadě odhady rozptylu a šikmosti b1 určené vztahy (3.1) a (3.2). Potom náhodná veličina "r # r p 1 Y 2 (W 2 − 1) Y 2 (W 2 − 1) X( b1 ) = √ log + +1 (3.3) 2 2 log W √ má asymptoticky rozdělení N(0, 1). Platí-li |X( b1 )| > u(α/2), zamítneme H0 . Lze použít i jednostrannou alternativu. Tento test lze v programu R najít v balíčku moments pod příkazem agostino.test(c, ...). Chceme-li použít jednostrannou alternativu, použijeme volitelný parametr alternative="less" případně alternative="greater".

3.3

Anscombeho-Glynnův test špičatosti

Obdobně jako u D’Agostinova testu šikmosti má za předpokladu nulové hypotézy ). Konvergence je zde pomalá, neboť test je konstruován b2 asymptoticky N(3, 24 n tak, aby byl kladný, ale shora není nijak omezen. V pracích Pearson [19] a [20] byly pro konečné výběry vypočteny s přesností n−3 i momenty b2 rozdělení. Ty jsou dány vztahy 17

3(n − 1) , (3.4) n+1 24n(n − 2)(n − 3) σ 2 (b2 ) = , (3.5) (n + 1)2 (n + 3)(n + 5) 216(n − 3)(n − 5)(n2 − 5n + 2)2 α3 (b2 ) = , (3.6) (n − 3)(n − 2)n(n + 5)2 (n + 7)2 36(15n6 − 36n5 − 628n4 + 982n3 + 5 777n2 − 6 402n + 900) . α4 (b2 ) = 3 + (n − 3)(n − 2)n(n + 7)(n + 9)(n + 11)(n + 13) µ(b2 ) =

V práci Anscombe a Glynn [22] byla za platnosti H0 pro n ≥ 20 představena transformace b2 , po které dostaneme testovou statistiku, která má asymptoticky N(0, 1). Nechť ζ =

b2 − µ(b2 ) , σ(b1 ) 8

A = 6+ p α3 (b2 )

2

p + α3 (b2 )

s

4 1+ α3 (b2 )

! ,

kde jsme opět využili odhady momentů b2 určené vztahy (3.4), (3.5) a (3.6). Potom náhodná veličina    31  r 2 1− A 9A  2   q (3.7) Z(b2 ) = −  1 − 2 2 9A 1+ζ A−4

má asymptoticky rozdělení N(0, 1). Platí-li |Z| > u(α/2), zamítneme hypotézu H0 . Velké kladné hodnoty Z předpovídají rozdělení s těžkými chvosty, velké záporné hodnoty Z předpovídají rozdělení s lehkými chvosty. Opět můžeme použít jednostrannou alternativu. Tento test lze v programu R najít v balíčku moments pod příkazem anscombe.test(c, ...). Chceme-li použít jednostrannou alternativu, použijeme volitelný parametr alternative="less" případně alternative="greater".

3.4

Jarqueho-Berův test

Tento test, který testuje zároveň šikmost i√špičatost, byl publikován v práci Jarque a Bera [23]. Pokud by testové statistiky b1 a b2 byly nezávislé a normálně rozdělené, potom by náhodná veličina √ 2  √ 2 b1 − µ( b1 ) b2 − µ(b2 ) √ J= + σ(b2 ) σ( b1 ) měla√χ2 rozdělení o dvou stupních volnosti. Jak již bylo v sekcích 3.2 a 3.3 zmíněno, b1 a b2 mají asymptoticky normální rozdělení a jsou i asymptoticky nezávislé, proto můžeme použít testovou statistiku   n m23 (m24 − 3)2 JB = + . 6 m32 4 18

Neplatí-li nerovnost χ22 (1 − α/2) ≤ JB ≤ χ22 (α/2), potom zamítneme H0 . Jelikož konvergence b2 k N(3, 24 ) je pomalá, doporučuje se test používat až pro výběry n o rozsahu n ≥ 200. Tento test lze v programu R najít v balíčku moments pod příkazem jarque.test(c).

3.5

Bowmanův-Shentonův test šikmosti a špičatosti

Tento test, který je kombinací testů 3.2 a 3.3, byl publikován v práci Bowman a Shenton [24]. Testová statistika vyplývá z (3.3) a (3.7) a je dána vztahem p Ke2 = X 2 ( b1 ) + Z 2 (b2 ). √ Protože X( b1 ) a Z(b2 ) mají asymptoticky rozdělení N(0, 1), testová statistika Ke2 má χ2 rozdělení o 2 stupních volnosti. H0 zamítneme, pokud neplatí nerovnost χ22 (1 − α/2) ≤ Ke2 ≤ χ22 (α/2). V práci Bowman a Shenton [24], obrázek 2, str. 248 √ najdeme graficky znázorněný kritický obor pro X( b1 ) a Z(b2 ) při použití testové statistiky Ke2 pro α = 0,1; 0,05; 0,01.

3.6

Obdélníkový test šikmosti a špičatosti

√ Tento test je založen na použití testových statistik b1 a b2 resp. jejich transformací (3.3) a (3.7). V práci Pearson, D’Agostino a Bowman [25] byla vypoč√ 0 tena testová hladina α = 1 − 1 − α, při které nezávisle použijeme základní test šikmosti a špičatosti (viz 3.1), případně D’Agostinův test šikmosti (viz 3.2) a Anscombeho-Glynnův test špičatosti (3.3), přičemž dohromady dostaneme celkovou testovou hladinu α.

3.7

Testy založené na absolutních centrálních momentech

Obdobně jako centrální momenty můžeme definovat absolutní centrální moment k-tého řádu z rozdělení s hustotou f (x) vztahem Z ∞ 0 |x − µ|k f (x)dx, k ∈ R r {0} . vk = −∞

Výběrový absolutní centrální moment k-tého řádu je n

vk =

1X |xi − x|k , k ∈ R r {0} . n i=1

Testové statistiky pro testování absolutních momentů jsou definovány analogicky jako u testu šikmosti a špičatosti a(k) =

vk k/2 v2

19

=

vk k/2

m2

,

kde k může být libovolné reálné číslo kromě 0 a 2, obvykle volíme přirozené číslo. Hypotézu H0 zamítneme, pokud vk neleží v intervalu daném dolní a horní kritickou hodnotou. Tabulku kritických hodnot najdeme v knize Thode [1], tabulka B7, str. 383–384.

3.8

Gearyho test

Tento test byl publikován v práci Geary [26] jako alternativa k b2 , protože při malém výběru byl přesnější. Gearyho testová statistika (původně označovaná wn0 ) má tvar n X |xi − x| . a = a(1) = s i=1n X n (xi − x)2 i=1

pπ

V optimálním případě by mělo platit a = 2 ≈ 0,7979. Hypotézu H0 zamítneme, pokud a neleží v intervalu daném horní a dolní kritickou hodnotou. Tabulku kritických hodnot lze nalézt v práci Geary [27], tabulka 3, str. 303. Protože korelační koeficient mezi a a b2 je záporný (viz Geary [28]), tak zde naopak při rozdělení s těžkými chvosty bude a nabývat malých hodnot, při rozdělení s lehkými chvosty, bude a nabývat velkých hodnot. V práci D’Agostino [29] byla představena transformace testové statistiky a, která √ . Náhodná využívá odhady střední hodnoty (0,7979) a směrodatné odchylky 0,2123 n veličina √ n(a − 0,7979) z= 0,2123 má za předpokladu platnosti H0 standardizované normální rozdělení. Hypotézu H0 zamítneme, pokud |z| > u(α/2). V programu R lze v balíčku moments najít příkaz geary(c), který vypočte hodnotu testové statistiky a. Porovnání s kritickou hodnotou, případně použití D’Agostinovy transformace, nebylo prozatím implementováno.

3.9

Bonettův-Seierové test

Tento test byl publikován v práci Bonett a Seierová [30]. V práci Geary [27] byla definovaná míra p šikmosti jako τ /σ, kde τ = E|X − µ|. Za předpokladu normality platí τ /σ = 2/π ≈ 0,7979 a tedy ω = 13,29(log τ − log σ) ≈ 3. Při použití √ odhadu ω b = 13,29(log v1 − log m2 ) dostaneme (po úpravách uvedených v práci Bonett a Seierová [30], str. 437) testovou statistiku √ n + 2(b ω − 3) zω = , 3,54 která má standardizované normální rozdělení. Hypotézu H0 zamítneme, pokud |zω | > u(α/2). Tento test lze v programu R najít v balíčku moments pod příkazem bonett.test(c).

20

4. Testy založené na empirických distribučních funkcích Jedná se o testy dobré shody, které jsou založeny na porovnání empirické a hypotetické distribuční funkce. Existují dva hlavní typy testů, a to lineární a kvadratické. Nechť x1 , x2 , ..., xn je náhodný výběr. Potom empirická distribuční funkce z tohoto výběru je schodovitá funkce definovaná vztahem   0, x < x(1) , i , x(i) ≤ x < x(i+1) pro i = 1, 2, ..., n − 1, Fn (x) =  n 1, x(n) ≤ x. Empirická distribuční funkce Fn (x) udává pravděpodobnost počtu pozorování, které jsou menší nebo rovny x, se skoky o velikosti n1 v pozorovaných hodnotách. Na obrázku 4.1 vidíme příklad empirické distribuční funkce výběru z R(0, 1) o rozsahu n = 20. Tečkovaně je znázorněna distribuční funkce R(0, 1). Protože distribuční funkce normálního rozdělení je spojitá, při testování stačí znát hodnoty   x(i) − µ b p(i) = Φ , (4.1) σ b neboť vzdálenost empirické distribuční funkce od distribuční funkce spojitého rozdělení je v bodech, kde nenastává skok, menší než v bodech skoku. Hodnoty p(i) dány vztahem (4.1) budou používány v celé kapitole.

1.0 0.8 0.6 0.4 0.2 0.0 0.0

0.2

0.4

0.6

0.8

1.0

Obrázek 4.1: Empirická distribuční funkce náhodného výběru z R(0, 1) o rozsahu n = 20

21

4.1 4.1.1

Kolmogorovův-Smirnovův test a Lillieforsova varianta Kolmogorovův-Smirnovův test

Tento test je založen na maximální lineární vzdálenosti empirické a hypotetické distribuční funkce. Testová statistika je dána vztahem   i + − p(i) , D = max i=1,...n n   i−1 − , D = max p(i) − i=1,...n n  (4.2) D = max D+ , D− . Velkou nevýhodou tohoto testu je, že pro výpočet p(i) musíme znát přesné hodnoty µ baσ b. Použitím odhadů střední hodnoty a směrodatné odchylky není zaručena přesná testová hladina α. V programu R dostaneme tento test příkazem ks.test(c, pnorm, µ b, σ b). Chceme-li testovat shodu s jiným rozdělením, použijeme místo pnorm, µ b, σ b příslušnou distribuční funkci požadovaného rozdělení s odpovídajícími parametry. Dnes má větší uplatnění dvouvýběrový Kolmogorovův-Smirnovův test, který testuje, mají-li dva různé nezávislé výběry stejnou distribuční funkci, více viz Anděl [2], kapitola 11.2.4, str. 240–243.

4.1.2

Lillieforsova varianta

V práci Lilliefors [31] byla představena varianta Kolmogorova-Smirnova testu, která nevyžaduje přesnou znalost µ baσ b, ale stačí jejich odhady. Test byl pro výběry o rozsahu n ≥ 100 pro dobrou vypovídající schopnost předepisován normou ČSN 01 0225 [5]. Pro výpočet hodnot p(i) použijeme známé odhady n

1 X (xi − x)2 . µ b = x, σ b = n − 1 i=1 2

Testová statistika D se dostane opět ze vztahu (4.2). Hypotézu H0 zamítneme, pokud D je větší než tabelovaná kritická hodnota. Tabulku kritických hodnot najdeme v práci Lilliefors [31], tabulka 1, str. 400 s opravou uvedenou v práci Lilliefors [32]. Tento test lze v Programu R najít v balíčku nortest pod příkazem lillie.test(c).

4.2

Kuiperův test

Tento test byl stejně jako Kolmogorovův-Smirnovův test založen na kombinaci veličin D+ a D− , a tedy na vzdálenosti empirické a hypotetické distribuční funkce. Testová statistika je dána nikoliv maximem, ale součtem těchto veličin, a tedy má tvar V = D+ + D− .

22

Jako test normality se dnes již nepoužívá, ale uplatnění má v případě, že chceme testovat, zdali jsou jednotlivá pozorvání rovnoměrně rozdělená na kružnici. Pro tento případ jej v programu R najdeme v balíčku circular pod příkazem kuiper.test(c).

4.3

Cramérův-von Misesův test

V práci Anderson a Darling [33] byla publikována třída kvadratických testů založených na empirické distribuční funkci definovaná vztahem Z ∞ [Fn (x) − F (x)]2 ψ(F (x))dF (x), (4.3) n −∞

kde ψ(F (x)) je váhová funkce. Cramérova-von Misesova testová statistika je dána váhovou funkcí ψ(F (x)) = 1 a je tvaru 2 n  X 1 2i − 1 W = + p(i) − . 12n i=1 2n 2

Tento test lze v programu R nalézt v balíčku nortest pod příkazem cvm.test(c).

4.4

Watsonův test

Tento test je alternativou ke Cramérovu-von Misesovu testu, když chceme testovat, zdali jsou jednotlivá pozorování rovnoměrně rozdělená na kružnici. Test byl publikován v práci Watson [34]. Testová statistika je dána vztahem 2 n  X 2i − 1 1 + p(i) − p + 0,5 − . U = 12n i=1 2n 2

Vztah mezi U 2 a W 2 lze vyjádřit ve tvaru U 2 = W 2 − n(p − 0,5)2 . Tento test najdeme v programu R v balíčku circular pod příkazem watson.test(c)

4.5

Andersonův-Darlingův test

Další z třídy kvadratických testů založených na empirické distribuční funkci definovaných vztahem (4.3) byl publikován v práci Anderson a Darling [35]. Anderson a Darling použili váhovou funkci ψ(F (x)) = [p(1 − p)]−1 a dostali testovou statistiku n

  1X A = −n − [2i − 1] log p(i) + log(1 − p(n−i+1) ) . n i=1 2

Pro dobrou vypovídající schopnost je tento test doporučován pro výběry o rozsahu n ≥ 50 normou ČSN 01 0225 [5]. V programu R tento test najdeme v balíčku nortest pod příkazem ad.test(c).

23

4.6

Stephensovy modif ikace

V práci Stephens [36] byly představeny modifikace testových statistik D, V , W 2 , U 2 a A2 tak, aby kritické hodnoty nezávisely na velikosti výběru n. Upravíme-li testové statistiky dle vztahů   √ 0,85 ∗ , D = D n − 0,01 + √ n   √ 0,82 ∗ V = V , n + 0,05 + √ n   0,5 2∗ 2 W = W 1+ , n   0,5 2∗ 2 U = U 1+ , n   0,75 2,25 2∗ 2 A = A 1+ + 2 , n n potom stačí použít kritické hodnoty uvedené v práci Stephens [36], tabulka 1A, část 1.3, str. 732. Tyto modifikace používá ve všech testech i program R. Hypotézu H0 zamítneme, pokud je hodnota testové statistiky větší než tabelovaná kritická hodnota.

24

5. Chí-kvadrát test dobré shody Jedná se o nejstarší test shody výběru s hypotetickým rozdělením, tedy i test normality. Byl prezentován již roku 1900 v práci Pearson [37]. Ze začátku byl tento test velmi oblíbený a používaný. Později, když se objevilo více testů normality, přestával být tento test používaný, vzhledem k tomu, že nebyl moc silný (viz Shapiro a Wilk [9]). Dnes se při testování normality dává přednost testům založeným na momentech (viz kapitola 3). Chí-kvadrát test dobré shody je dobře použitelný pro testování shody s rozdělením, u kterého nemusíme odhadovat žádné parametry.

5.1

Chí-kvadrát test při známých parametrech

Tento test použijeme, chceme-li testovat shodu našeho výběru s rozdělením, jehož všechny parametry známe. V případě normality známe střední hodnotu a rozptyl. Nechť x1 , x2 , ..., xn je náhodný výběr. Nejprve vytvoříme třídy (−∞, b1 ], (b1 , b2 ], (b2 , b3 ], ..., (bk−2 , bk−1 ], (bk−1 , ∞), tak, aby pravděpodobnost pi , i = 1, 2, ..., k, že náhodná veličina padne do i-té třídy, byla nenulová. Pro přehlednost označme i-tou třídu symbolem Ji . Optimální volby počtu a šířky tříd jsou popsány v sekci 5.3, přičemž nezapomeneme na Yarnoldovo pravidlo (viz sekce 5.4). Dále spočteme náhodné veličiny Yi , které budou představovat počet pozorování v jednotlivých třídách, tedy Yi =

n X

I(xj ∈ Ji ), i = 1, 2, ..., k.

j=1

Potom Pearsonova statistika χ2 =

k X (Yi − npi )2

npi

i=1

(5.1)

má asymptoticky χ2 rozdělení o k − 1 stupních volnosti (viz Anděl [2], věta 12.5). Hypotézu H0 zamítneme, pokud χ2 ≤ χ2k−1 (α). Tento test v programu R dostaneme příkazem chisq.test(Y, p = r), kde Y = (Y1 , Y2 , ..., Yk )> a r = (p1 , p2 , ..., pk )> .

5.2

Chí-kvadrát test při neznámých parametrech

Většinou pravděpodobnosti pi , i = 1, 2, ..., k závisí na nějakém neznámém parametru a = (a1 , ..., am )> , tedy můžeme psát pi = pi (a), i = 1, 2, ..., k. Odhad b a parametru a získáme vyřešením soustavy k X Yi ∂pi (a) = 0. p (a) ∂a i j i=1

25

(5.2)

Postup, jak jsme dospěli k této soustavě, je uveden v knize Anděl [2], kapitola 12.3, str. 272–273. Máme-li odhad b a vypočtený, aplikujeme postup uvedený v sekci 5.1 s předpokladem, že náš výběr je z rozdělení se známým parametrem b a, ale v tomto případě bude mít testová statistika asymptoticky χ2 rozdělení o k −m−1 stupních volnosti. Při testování normality dostaneme z (5.2) po úpravě soustavu k

1 X Yi µ= n i=1 pi

k

Z

1 X Yi xf (x)dx, σ = n i=1 pi Ji 2

Z

(x − µ)2 f (x)dx.

Ji

Ta se zpravidla řeší iteračně. Tento test pro případ normality najdeme v programu R v balíčku nortest pod příkazem pearson.test(c). Program R používá pro odhad parametrů µ a σ 2 známé odhady n

b2 = µ b = x, σ

1 X (xi − x)2 , n − 1 i=1

jejichž výpočet je jednodušší, ale zvyšuje pravděpodobnost chyby prvního druhu.

5.3

Počet a šířka tříd

Nevýhodou χ2 testu je, že pro stejný výběr může být výsledek testů silně závislý na počtu a šířce vybraných tříd. Jedna z možností, jak zvolit šířku tříd, byla prezentována v práci Mann a Wald [38]. Ti navrhovali, aby každá z k tříd měla stejnou pravděpodobnost, tedy aby pi = k. V tomto případě aproximace Pearsonovy statistiky (5.1) je přesnější a lze místo (5.1) použít ekvivalentní vztah k

χ2 = −n +

kX 2 Y . n i=1 i

V knize Thode [1] je představena Kendallova a Stuartova volba k maximalizující sílu testu, za předpokladu, že síla testu je β, přičemž √

k=b

2n2 [Φ−1 (1 − α)]2 + T (β)

! 52 ,

(5.3)

kde T (β) je β-kvantil hypotetického rozdělení a b volitelný kladný parametr. Dále bylo doporučeno používat (5.3) pouze pro výběry o rozsahu n ≥ 200. Mann a Wald [38] uvažovali případ, kdy β = 0,5 a b = 4. V práci Koehler a Larntz [39] bylo ukázano, že Pearsonův χ2 test je adekvátní na hladině α = 0,05 a α = 0,01, pokud npi ≤ 0,25 pro všechna i = 1, 2, ..., k, když k ≥ 3, n ≥ 10 a n2 /k ≥ 10. Program R používá pro výpočet k vztah k = 2n2/5 . V případě testu, zdali náš výběr pochází z diskrétního rozdělení, které s pravděpodobností 1 nabývá konečně mnoha hodnot, případně pokud pozorovatel zapsal četnosti do konečně mnoha skupin, zařadíme, pokud je to možné, každou hodnotu x splňující P(X = x) > 0, případně každou skupinu, do zvláštní třídy. Mezi případy, kdy není možné zařadit každou hodnotu do zvláštní třídy řadíme případy, kdy počet hodnot x splňující P(X = x) > 0 je větší nebo blízký počtu 26

pozorování, případně máme-li data pro nějaká pozorování sloučena (např. pozorovatel již nerozlišoval, kolikrát nastal daný jev, nastal-li již více než t-krát). V případě, že by nebylo dodrženo Yarnoldovo pravidlo (viz sekce 5.4), slučujeme třídy s nejmenšími hodnotami npi se sousedními třídami do té doby, než bude dodrženo.

5.4

Yarnoldovo pravidlo

Protože χ2 test je asymptotický, je vhodné jej používat pouze při velkém rozsahu výběru n. Obvykle se uvádí, že musí platit npi ≥ 5 pro každé i = 1, 2, ..., n. Dnes se, v případě že neodhadujeme žádné parametry, spíše užívá Yarnoldovo kritérium (viz Yarnold [40] a Eaton [41]). V tomto případě stačí, aby platilo npi ≥ 5s/k pro všechna i = 1, 2, ..., n, kde s je počet tříd, pro něž platí npi < 5. Odhadujeme-li nějaké neznámé parametry před použitím χ2 testu, je vhodné upustit od Yarnoldova pravidla a použít přísnější původní variantu.

27

6. Další testy normality 6.1

Test založený na rozpětí

Tento test byl publikován v práci David, Hartley a Pearson [42]. Testová statistika má tvar
√ n − 1 x(n) − x(1) u= . n X (xi − x)2 i=1

Tento test je vhodný, očekáváme-li, že náš výběr je z rovnoměrného nebo jiného rozdělení s lehkými chvosty. Hypotézu H0 zamítneme, pokud je u větší než tabelovaná dolní kritická hodnota. Test založený na rozpětí lze použít i jako test na jedno odlehlé pozorování (viz kapitola 7), přičemž zde H0 zamítneme, pokud je u menší než tabelovaná horní kritická hodnota. Tabulku dolních a horních kritických hodnot najdeme v práci David a kol. [42], tabulka 6, str. 491.

6.2

Uthoffův test

Tento test byl pubilkován v práci Uthoff [43]. Test je podobný Gearyho testu (viz sekce 3.8), pouze výběrový průměr byl nahrazen mediánem. Testová statistika má tvar n X |xi − x˜| , U = s i=1n X 2 n (xi − x) i=1

kde x˜ = med {x1 , x2 , ..., xn }. V práci Uthoff [43] bylo dokázáno, že Uthoffův U test je asymptoticky ekvivalentní s Gearyho a(1) testem. Tento test je vhodný proti alternativě Laplaceova rozdělení. Hypotézu H0 zamítáme, pokud U neleží mezi dolní a horní kritickou hodnotou. Tabulku kritických hodnot najdeme v knize Thode [1], tabulka B9, str. 393.

6.3

D’Agostinův test

Tento test byl publikován v práci D’Agostino [44] jako rozšíření Shapirova-Wilkova testu (viz sekce 2.1) pro středně velké a velké výběry. Testová statistika D’Agostinova testu má tvar  n  X n+1 x(i) i− 2 i=1 . D= s n X 2 n3 (xi − x) i=1

28

Za předpokladu platnosti H0 jsou střední hodnota a směrodatná odchylka statistiky D rovny 1 √ ≈ 0,282 094 79, 2 π 0,029 985 98 √ σ(D) ≈ . n E(D) =

Kritické hodnoty pro D jsou prezentovány pro standardizovanou statistiku √ (D − E(D)) n(D − 0,282 094 79) Y = ≈ . σ(D) 0,029 985 98 Hypotézu H0 zamítneme, pokud Y neleží mezi horní a dolní kritickou hodnotou. Tabulku kritických hodnot najdeme pro n > 100 v práci D’Agostino [44], tabulka 1, str. 343, pro n ≤ 100 v práci D’Agostino [45], tabulka 1, str. 220.

6.4

Ojovy testy

Tato sada testů byla publikována v práci Oja [46]. Testy jsou silné v případě, že alternativní rozdělení má jinou šikmost nebo špičatost než teoretické rozdělení, v našem případě normální. K detekci šikmosti, je použita testová statistika T1 , k detekci špičatosti testová statistika T2 , které mají tvar  −1 X x(k) − x(j) n T1 = 3 x(k) − x(i) 1≤i<j u(α/2) při testu na různou špičatost, |T2 | > u(α/2) při testu na různou šikmost, nebo je-li S ≤ χ22 (α) při kombinovaném testu.

6.5

Vašíčkův test

Test publikovaný v práci Vašíček [47] je založen na odhahu výběrové entropie. Nechť f (x) je hustota rozdělení, potom entropie je dána vztahem Z ∞ H(f ) = f (x) log(f (x))dx. −∞

29

Odhad H(f ) lze vypočíst ze vztahu n

Hmn =

  n 1X (x(i+m) − x(i−m) ) , log n i=1 2m

kde m je přirozené číslo menší než n2 a x(j) = x(1) pro j < 1 a x(j) = x(n) pro j > n. Pro všechna rozdělení s rozptylem σ 2 √ je H(f ) maximalizována normálním √ rozdělením. V tomto případě je H(f ) = log( 2πeσ), takže exp[H(f )]/σ ≤ 2πe pro každou hustotu f (x). Testová statistika má tvar

Kmn

" n #1 √ Y
n n x(i+m) − x(i−m) . = n X 2m (xi − x)2 i=1 i=1

Hypotézu H0 zamítáme, pokud je Kmn menší než tabelovaná kritická hodnota. Tabulku kritických hodnot najdeme v práci Vašíček [47], tabulka 1, str. 58.

30

7. Testy na odlehlá pozorování V této kapitole budeme uvažovat nulovou hypotézu H0 , že ve výběru není žádné odlehlé pozorování, oproti alternativě H1 , že ve výběru je právě k odlehlých pozorování. Velmi často bývá důležité zjistit, jsou-li extrémní hodnoty v našem výběru odlehlé. Pokud existují v našem výběru odlehlé hodnoty, potom může nastat problém, neboť u některých testů může nastat významný rozdíl při zachování a odstranění extrémní hodnoty. Pro samotné testování normality se testy na odlehlá pozorování používají jen zřídka. Jejich hlavním úkolem je zjistit, zda-li nějaká pozorování mají velký vliv při odhadech nebo testech. Odlehlá pozorování vznikají z různých příčin, mezi nejčastější patří: 1. Náhoda. Je-li výběr z normálního rozdělení, dá se očekávat, že přibližně jedno z 20 pozorování bude od střední hodnoty vzdáleno o více než dvojnásobek směrodatné odchylky, jedno ze 100 může být vzdáleno i o více než 2,5 násobek střední hodnoty. Přítomnost odlehlého pozorování tedy nemusí nutně znamenat problém. 2. Chyba při generování dat. V tomto případě může být identifikace odlehlých pozorování velmi důležitá při zjišťování chyby v experimentu. 3. Měřený objekt není homogenní s ostatními. V tomto případě je rozumné předpokládat, že výběr pochází z kontaminovaného rozdělení. 4. Chyba měření, způsobena technickou vadou nebo selháním experimentátora. 5. Chyba v zaznamenání dat, způsobena např. chybným zápisem naměřených hodnot, přepisem dat do počítače apod. Pokud již v našem výběru máme nějaká odlehlá pozorování, zkusíme zjistit, z jaké příčiny se objevila. Bylo-li to způsobeno chybou v zaznamenání dat, stačí příslušné hodnoty opravit. Jsme-li si jisti, že daná pozorování jsou chybná, ale nemáme možnost získat správné hodnoty, potom jediná možnost je tyto hodnoty zahodit. V tomto případě ale musíme také zahodit všechna další pozorování, která byla danou chybou ovlivněna. Nebyla-li zjištěna žádna chyba, ale přesto se nám některá odlehlá pozorování zdají natolik extrémní, že bychom s nimi raději nepracovali, potom máme několik možností: 1. použít metody robustní analýzy, 2. odstranit odlehlá pozorování a pokračovat v analýze, 3. odstranit odlehlá pozorování a předpokládat, že redukovaný výběr je cenzorovaný, 4. nahradit odlehlá pozorování nejbližšími neodlehlými, 5. nahradit odlehlá pozorování novými 31

6. použít standardní metody analýzy na původní data a data bez odlehlých hodnot, poté publikovat oba výsledky. Zde může nastat problém, pokud se oba výsledky od sebe velmi významně liší. Pro zjišťování odlehlých pozorování používáme často i grafické metody. Histogram (viz 1.1), stem-and-leaf plot (viz 1.2) nebo boxplot (viz 1.3) mohou být použity pro určení odlehlých pozorování. Testy na odlehlá pozorování nám mohou pomoci určit, zda-li není náš výběr z kontaminovaného rozdělení s nízkým faktorem výskytu, případně mohou indikovat u výběru o malém rozsahu, že se jedná o rozdělení s těžkými chvosty, případně o zešikmené rozdělení. Testy na odlehlá pozorování uvedené v této kapitole předpokládají znalost přesného počtu možných odlehlých pozorování. Máme-li ve výběru 2 nebo více odlehlých pozorování, potom test na jedno odlehlé pozorování jej nemusí odhalit, neboť ostatní mohou jeho existenci „maskovatÿ. V případě, že testujeme přítomnost k > 1 odlehlých pozorování a v našem výběru jich je méně než k, nastane jedna z následujících možností: 1. Test zamítne H0 protože vliv skutečných odlehlých hodnot byl tak velký, že způsobil identifikaci více odlehlých pozorování, než jich skutečně je. 2. Skutečná odlehlá pozorování nemají takový vliv, aby byla testová statistika signifikantní, takže test nevyjde průkazně a nelze identifikovat žádná odlehlá pozorování.

7.1

Grubbsovy testy na odlehlá pozorování

Tyto testy byly představeny v práci Grubbs [48].

7.1.1

Test na jedno odlehlé pozorování

Existují tři varianty testu na jedno odlehlé pozorování, které se liší podle toho, zdali je odlehlé pozorování minimem nebo maximem ve výběru. Testovou statistiku T1 použijeme za předpokladu, že x(1) je odlehlé pozorování. Za předpokladu, že x(n) je odlehlé pozorování, použijeme Tn . Pokud nejsme schopni učinit předpoklad, která z hodnot x(1) a x(n) je případným odhlelým pozorováním, použijeme T . Testové statistiky jsou dány vztahem √ n − 1 (x − x(1) ) , T1 = s n X (xi − x)2 i=1

Tn

√ n − 1 (x(n) − x) s n , = X (xi − x)2 i=1

T = max {T1 , Tn } Hypotézu H0 zamítneme, pokud je příslušná testová statistika větší než tabelovaná kritická hodnota. Tabulku kritických hodnot najdeme v knize Thode [1], 32

tabulka B8, str. 386–392. Grubbs [48] dále ukázal, že Tn a T1 se dají vyjádřit pomocí n−1 X (xi − xn )2

Sn2 = S2

i=1 n X

= 1−

Tn2 , n−1

= 1−

T12 , n−1

(xi − x)2

i=1

n X (xi − x1 )2

S12 = S2

i=2 n X

(xi − x)2

i=1

kde xn resp. x1 je výběrový průměr s vynecháním x(n) resp. x(1) . Tyto testové statistiky se nepoužívají, protože výpočet samotných hodnot T1 a Tn je jednodušší, ale mají souvislost s dalšími Grubbsovými testy.

7.1.2

Test na více odlehlých pozorování na konkrétním chvostu

Testujeme-li existenci dvou odlehlých pozorování na jednom chvostu, využijeme testové statistiky n−2 X 2 Sn−1,n 2

S

=

(xi − xn−1 )2

i=1 n X

,

(xi − x)2

i=1 n X 2 S1,2 2

S

(xi − x2 )2

=

i=3 n X

, 2

(xi − x)

i=1

kde xn−1 resp. x2 je výběrový průměr s vynecháním x(n) a x(n−1) resp. x(1) a x(2) . Testové statistiky se dají zobecnit pro k odlehlých pozorování na jednom chvostu a mají tvar n−k X (xi − xn−k )2

Lk =

i=1 n X i=1 n X

L∗k =

, (xi − x)2

(xi − x∗k )2

i=k+1 n X

, 2

(xi − x)

i=1

33

kde xn−k resp. x∗k je výběrový průměr s vynecháním k posledních resp. k prvních pozorování. Protože Lk i L∗k jsou menší než 1, hypotézu H0 zamítáme, pokud je hodnota jednotlivé testové statistiky menší než tabelovaná kritická hodnota. Tabulku kritických hodnot najdeme v knize Thode [1], tabulka B18, str. 405–410.

7.1.3

Test na jedno odlehlé pozorování na obou chvostech

Obdobně jako další Grubssovy testy, má i testová statistika pro jedno odlehlé pozorování na obou chvostech tvar n−1 X (xi − x1,n )2

2 S1,n i=2 n 2 = X S

, 2

(xi − x)

i=1

kde x1,n je výběrový průměr s vynecháním prvního a posledního pozorování. 2 /S 2 menší než tabelovaná kritická hodnota. Hypotézu H0 zamítáme, pokud je S1,n Tabulku kritických hodnot najdeme v knize Thode [1], tabulka B19, str. 411.

7.1.4

Grubbsovy testy v programu R

V programu R lze najít všechny varianty Grubbsových testů v balíčku outliers pod příkazem grubbs.test(c, ...). Volbou prvního parametru type zvolíme typ, který chceme testovat. Implicitně je nastaven typ 10 pro jedno odlehlé pozorování, dále můžeme použít typ 11 pro jedno odlehlé pozorování na obou chvostech a typ 20 pro dvě odlehlá pozorování na jednom chvostu (pro tento případ je bohužel omezen rozsah výběru na 3 ≤ n ≤ 30). Test pro více odlehlých pozorování není prozatím v programu R implementován. Při použití typu 10 nebo 20 program R předpokládá, že odlehlé pozorování je to, které je nejvíce vzdáleno od x, v případě typu 20 k němu přidá jemu nejbližší pozorování. Chceme-li otestovat odlehlost pozorování na opačné straně, použijeme volitelný parametr opposite = TRUE. Chceme-li testovat pomocí testové statistiky T použijeme volitelný parametr two.sided = TRUE.

7.2

Dixonův test na jedno odlehlé pozorování

Tento test na jedno odlehlé pozorování byl publikován v práci Dixon [49]. Základní testová statistika pro testování odlehlosti x(n) je dána jako poměr mezi vzdáleností dvou největších pozorování a rozsahem, tedy ve tvaru r10 =

x(n) − x(n−1) . x(n) − x(1)

Obdobně testová statistika pro testování odlehlosti x(1) je tvaru 0 r10 =

x(2) − x(1) . x(n) − x(1) 34

V případě, že nejsme schopni učinit předpoklad, která z hodnot x(1) nebo x(n) 0 je odlehlým pozorováním, můžeme použít testovou statistiku r = max {r10 , r10 }. Vzhledem k tomu, že si Dixon byl vědom toho, že přítomnost dalšího pozorování blízko zkoumaného může „maskovatÿ odlehlost, navrhl i alternativní testové statistiky x(n) − x(n−j) , x(n) − x(k+1) x(j+1) − x(1) = , x(n−k) − x(1)

rjk = 0 rjk

které eliminují vliv j − 1 dalších extrémních pozorování. V praxi se používají pouze hodnoty j = 1, 2 a k = 0, 1, 2. Hypotézu H0 zamítneme, pokud je příslušná testová hodnota větší než tabelovaná kritická hodnota. Tabulky kritických hodnot najdeme v práci Dixon, [50], tabulky 1–6, str. 73–78. V programu R najdeme tyto testy v balíčku outliers pod příkazem dixon.test(c,...). Volbou prvního parametru type zvolíme typ, který chceme testovat. Číslo odpovídá oběma indexům u testové statistiky. Implicitně je nastaven typ 10 pro výběry o rozsahu 3–7, typ 11 pro výběry o rozsahu 8–10, typ 21 pro výběry o rozsahu 11–13 a typ 22 pro výběry o rozsahu 14–30. Pro větší výběry prozatím nejsou v programu R k dispozici kritické hodnoty. Volba, která pozorování budou zkoumána jako odlehlá, je stejná jako u Grubbsova testu. Pro test odlehlosti pozorování na opačné straně, použijeme volitelný parametr opposite = TRUE. Test založený na testové statistice r dostaneme volitelným parametrem two.sided = TRUE.

7.3

Tietjenův-Mooreův test na více odlehlých pozorování na jednom nebo obou chvostech

Tento test, který je obdobou Grubbsových statistik Lk a L∗k , vyšetřuje odlehlost pro k nejvzdálenejších pozorování od střední hodnoty bez ohledu na to, ze kterého chvostu pocházejí. Byl publikován v práci Tietjen a Moore [51]. Výběr x1 , x2 , ..., xn uspořádáme podle vzdálenosti od výběrového průměru a přeznačíme na z(1) , z(2) , ..., z(n) . Testová statistika má tvar n−k X (z(i) − z n−k )2

Ek =

i=1 n X

,

(x(i) − x)2

i=1

kde z n−k je výběrový průměr ze z(1) , z(2) , ..., z(n−k) . Hypotézu H0 zamítneme, pokud je hodnota testové statistiky menší než tabelovaná kritická hodnota. Tabulku kritických hodnot najdeme v knize Thode [1], tabulka 20, str. 412–417.

35

8. Porovnání některých testů normality na testovacích datech V této kapitole porovnáme sílu jednotlivých testů na základě simulací v programu R. Budeme porovnávat pouze testy, které jsou naprogramovány v programu R, tedy Shapirův-Franciův test (SF), D’Agostinův test šikmosti (DA), AnscombehoGlynnův test špičatosti (AG), Jarqueho-Berův test (JB), Bonettův-Seierové test (BS), Lillieforsovu variantu Kolmogorova-Smirnova testu (LI), Cramérův-von Misesův test (CM), Andersonův-Darlingův test (AD) a χ2 test (χ2 ). Možností, jak jednotlivé testy porovnat, je spousta, já jsem se rozhodl porovnat sílu jednotlivých testů na náhodných výběrech o rozsahu n = 100 a n = 200 ze standardizovaného normálního rozdělení, exponenciálního rozdělení se střední hodnotou 1, Laplaceova rozdělení s parametry a = 0 a b = 1, Cauchyova rozdělení, rovnoměrného rozdělení na intervalu (−1; 1) a ze standardizovaného normálního rozdělení, které bude kontaminováno Cauchyovým rozdělením s pravděpodobností p = 0,05. Z každého rozdělení vygenerujeme 1000 náhodných výběrů a zjistíme, jestli jednotlivý test zamítne hypotézu H0 na hladině α = 0,05. Dále jsem se rozhodl zjistit, jak jsou dané testy citlivé na přítomnost odlehlých pozorování, tedy jak se změní výsledky, rozhodneme-li se extrémní hodnotu zahodit. Pro zjištění přítomnosti odlehlého pozorování použijeme Grubbsův test. V případě nalezení odlehlého pozorování znovu provedeme stejný test na tentýž výběr, přičemž odlehlé pozorování z něj vyřadíme. V tabulce jsou uvedeny vždy pouze počty výběrů, které splnily příslušnou vlastnost. V přiloženém softwarovém kódu najdeme i pořadí testů, u kterých došlo ke splnění dané vlastnosti. U tabulek, kde jsou sloupce nadepsány pouze čísly, platí pro jednotlivé sloupce: 1. byla zamítnuta H0 , 2. platí 1. a existuje odlehlé pozorování, 3. platí 1., ale neexistuje odlehlé pozorování, 4. platí 2. a zamítáme H0 i bez odlehlého pozorování, 5. platí 2., ale nezamítáme H0 bez odlehlého pozorování, 6. nebyla zamítnuta H0 , 7. platí 6. a existuje odlehlé pozorování, 8. platí 6., ale neexistuje odlehlé pozorování, 9. platí 7., ale H0 zamítáme bez odlehlého pozorování, 10. platí 7. a H0 nezamítáme i bez odlehlého pozorování.

36

8.1

Normální rozdělení

Nejprve zjistíme, jak jednotlivé testy reagují na výběr z normálního rozdělení. V tomto případě by počet testů, u kterých zamítneme H0 , měl odpovídat přibližně αn = 50 testům. V tabulce 8.1 vidíme výsledky pro 1 000 testů o rozsahu n = 100, v tabulce 8.2 jsou výsledky pro 1 000 testů o rozsahu n = 200. Pokud si odmyslíme DA, jehož alternativou je šikmost vzdálená od 0, které jde u náhodného výběru z normálního rozdělení velmi těžko dosáhnout, tak až na drobné odchylky výsledky odpovídají našemu předpokladu o 50 zamítnutých výběrech. Dále vidíme, že SF, AG a JB mají větší tendenci zamítat výběry obsahující odlehlá pozorování, navíc u nich je i velké procento testů, kde po odstranění odlehlého pozorování již H0 zamítnuta nebyla. Počet testů, které s původním výběrem vyšly nevýznamně, ale po odstranění odlehlého pozorování již vyšly významně je u všech testů velmi nízký. Závěr: na normální rozdělení reagují testy tak, jak by měly. SF, AG a JB jsou citlivější na odlehlá pozorování. Tabulka 8.1: Výsledky 1 000 testů náhodného výběru z N(0, 1) o rozsahu n = 100 Test 1 SF 57 DA 8 AG 54 JB 53 BS 57 LI 52 CM 66 AD 57 2 χ 45

2 37 8 26 43 17 6 13 13 7

3 4 5 6 20 4 33 943 0 0 8 992 28 1 25 946 10 7 36 947 40 5 12 943 48 3 3 948 53 4 9 934 44 4 9 943 38 3 4 955

7 65 94 76 59 85 96 89 89 95

8 9 10 878 0 65 898 0 94 870 4 72 888 0 59 858 4 81 852 0 96 845 1 88 854 1 88 860 2 93

Tabulka 8.2: Výsledky 1 000 testů náhodného výběru z N(0, 1) o rozsahu n = 200 Test 1 2 3 4 5 6 7 SF 51 23 28 2 21 949 70 DA 1 0 1 0 0 999 93 AG 51 17 34 1 16 949 76 JB 40 24 16 1 23 960 69 BS 45 6 39 1 5 955 87 LI 59 6 53 4 2 941 87 CM 57 4 53 2 2 943 89 AD 54 4 50 2 2 946 89 2 χ 45 2 43 2 0 955 91

37

8 9 10 879 0 70 906 0 93 873 4 72 891 0 69 868 2 85 854 0 87 854 1 88 857 0 89 864 3 88

8.2

Exponenciální rozdělení

Nyní zjistíme, jak testy reagují na výběr z rozdělení, které je zešikmené. Očekáváme, že počet testů, které zamítnou H0 , by se ideálně měl blížit 1 000. Výsledky vidíme v tabulce 8.3. Je vidět, že většina testů zamítla H0 pro úplně všechny výběry a nezaváhaly i při odstranění odlehlého pozorování. U AG a BS rostla „úspěšnostÿ s rozsahem výběru. Při odstranění odlehlého pozorování docházelo často hlavně u BS k nezamítnutí H0 . Závěr: na exponenciální rozdělení reaguje většina testů výborně, AG test špičatosti bychom v tomto případě většinou nezkoušeli použít, BS je při alternativě exponenciálního rozdělení slabý. Tabulka 8.3: Výsledky 1 000 testů náhodného výběru z Exp(1) Test 1 DA, n = 100 999 AG, n = 100 870 AG, n = 200 988 BS, n = 200 547 BS, n = 100 791 ostatní, n = 100 1 000 ostatní, n = 200 1 000

8.3

2 902 851 972 541 787 902 976

3 97 19 16 6 4 98 24

4 5 6 7 8 9 10 897 5 1 0 1 0 0 665 186 130 51 79 0 51 928 44 12 4 8 0 4 289 252 453 361 92 8 353 557 230 209 189 20 4 185 902 0 0 0 0 0 0 976 0 0 0 0 0 0

Cauchyho rozdělení

Cauchyho rozdělení sice nemá žádné momenty, ale protože se jedná o rozdělení s hodně těžkými chvosty, tak očekáváme, že počet výběrů, kdy bude zamítnuta H0 se bude blížit 1 000. Všechny testy s výjimkou DA nakonec u všech 1 000 výběrů zamítly H0 . Po odstranění odlehlého pozorování, jehož existence byla zjištěna u všech výběrů, nedošlo ani v jednom případě k nezamítnutí H0 . DA zamítl H0 „ jenÿ 937-krát, s tím, že 161-krát nezamítl H0 po odstranění odlehlého pozorování. Ze 63 případů, kdy H0 zamítnuta nebyla, se 51-krát zamítla H0 po odstranění odlehlého pozorování. Nejpravděpodobněji se tak stalo z důvodu existence podobně odlehlého pozorování na druhém chvostu. Testy na výběry o rozsahu n = 200 nebyly vzhledem k jednoznačným výsledkům u menších výběrů provedeny. Závěr: až na AD test, který se ale na symetrické výběry nehodí, ale i přesto dopadl slušně, testy prokázaly schopnost odhalit, že výběr nepochází z normálního rozdělení.

8.4

Laplaceovo rozdělení

Nyní zjistíme, jak testy reagují na výběr z Laplaceova rozdělení. Toto rozdělení je, stejně jako normální, symetrické a exponenciálního typu. Liší se tím, že má těžší chvosty. Počet testů, které zamítnou H0 , by se měl blížit 1 000, ale protože Laplaceovo rozdělení není tak odlišné jako Cauchyho nebo exponenciální, očekáváme, že bude existovat nemalý počet výběrů, kde H0 nezamítneme. V tabulce 38

8.4 vidíme výsledky pro 1 000 testů o rozsahu n = 100, v tabulce 8.5 jsou výsledky pro 1 000 testů o rozsahu n = 200. Odmyslíme-li si DA, který nemá u Laplaceova rozdělení smysl, neboť jeho šikmost je rovna 0, je u ostatních testů patrné, že s rostoucím počtem pozorování ve výběru roste počet testů, kde zamítáme H0 . Při porovnání jednotlivých testů je vidět, že BS je výrazně nejsilnější, následován SF, CM a AD. Pro malé rozsahy výběru nejsou AG a JB založené na špičatosti u Laplaceova rozdělení silné. Výsledky LI dopadly obdobně jako u AG a JB, neboť distribuční funkce normálního a Laplaceova rozdělení jsou velmi podobné, takže u výběru s malým rozsahem nemůžeme očekávat velkou sílu LI. Nejslabší je u Laplaceova rozdělení χ2 test, dalo by se říci, že i nepoužitelný. Dále je vidět, že po odstranění odlehlého pozorování zejména u výběrů o rozsahu n = 100 již H0 ve velkém množství případů zamítnuta nebyla. Stejně jako u normálního rozdělení vidíme, že nejcitlivěji na odlehlé pozorování reagují AG a JB. V případě nezamítnutí H0 a odstranění odlehlého pozorování byla H0 výrazně zamítána pouze u χ2 testu. Závěr: až na AD, který u symetrických rozdělení nemá smysl používat, dopadly testy s výjimkou χ2 testu velmi dobře. Je vidět, že v tomto případě má velký vliv na výsledek rozsah výběru a ponechání nebo odstranění odlehlého pozorování. Tabulka 8.4: Výsledky 1 000 testů náhodného výběru z DEx(0, 1), rozsah n = 100 Test 1 SF 859 DA 227 AG 768 JB 796 BS 910 LI 709 CM 826 AD 832 2 χ 466

2 679 224 667 676 686 539 622 629 368

3 180 3 101 120 224 170 204 203 98

4 505 27 421 440 577 377 496 496 218

5 174 197 246 236 109 162 126 133 150

6 7 8 9 10 141 31 110 0 31 773 486 287 26 460 232 43 189 0 43 204 34 170 2 32 90 24 66 0 24 291 171 120 4 167 174 88 86 0 88 168 81 87 0 81 534 342 192 21 321

Tabulka 8.5: Výsledky 1 000 testů náhodného výběru z DEx(0, 1), rozsah n = 200 Test 1 SF 984 DA 258 AG 959 JB 959 BS 996 LI 952 CM 986 AD 984 2 χ 778

2 835 254 833 833 837 806 830 828 666

3 149 4 126 126 159 146 156 156 112

4 5 6 7 8 9 10 804 31 16 3 13 0 3 53 21 742 584 158 48 536 765 68 41 5 36 0 5 762 71 41 5 36 0 5 831 6 4 1 3 0 1 749 57 48 32 16 2 30 804 26 14 8 6 0 8 805 23 16 10 6 0 10 561 105 222 172 50 23 149

39

8.5

Rovnoměrné rozdělení

Nyní zjistíme, jak testy reagují na výběr z rozdělení s lehkými chvosty. Počet testů, které zamítnou H0 , by se měl blížit 1000, ale vzhledem k symetrii se dá očekávat, že budou existovat i výběry, kde H0 nazamítneme. V tabulce 8.6 vidíme výsledky pro 1000 testů o rozsahu n = 100 a n = 200. U rovnoměrného rozdělení není překvapením, že odlehlé pozorování neexistovalo ani v jednom případě. Z výsledků vidíme, že s rostoucím počtem pozorování ve výběru roste počet výběrů, kdy byla zamítnuta H0 , u JB je nárůst velmi velký. Z obdobných důvodů jako u Laplaceova rozdělení nedopadl LI nejlépe, χ2 test byl znovu nejslabší. Nejsilnější byly AG, SF a BS. Výsledky DA nepřekvapí, neboť šikmost je u rovnoměrného rozdělení rovněž 0. Závěr: až na AD, který u symetrických rozdělení nemá smysl používat, dopadly testy s výjimkou LI a χ2 testu velmi dobře. Je vidět, že i v tomto případě má velký vliv na výsledek rozsah výběru zvláště u JB. Tabulka 8.6: Výsledky 1 000 testů náhodného výběru z R(−1, 1) Test SF DA AG JB BS LI CM AD χ2

8.6

zamítnuta H0 n = 100 957 0 1 000 581 932 577 839 943 431

nezamít. H0 n = 100 43 1000 0 419 68 423 161 57 569

zamítnuta H0 n = 200 1 000 0 1 000 1 000 1 000 936 997 1 000 904

nezamít. H0 n = 200 0 1 000 0 0 0 64 3 0 96

Kontaminované rozdělení

Nakonec zjistíme, jak jednotlivé testy reagují na výběr z normálního rozdělení, které je kontaminováno Cauchovým s pravděpodobností p = 0,05. Testy na kontaminované rozdělení s touto pravděpodobností používali i Venables a Ripley [4], str. 110–111. V tomto případě je těžké předpovídat výsledky, neboť vliv Cauchyho rozdělení u jednotlivých výběrů může být velmi rozdílný, protože nevíme přesně, kolik pozorování z něj bude. Navíc hustota Cauchyho rozdělení je, stejně jako u normálního, nejvíce soustředěna kolem nuly, takže vliv těchto pozorování může být i zanedbatelný. V tabulce 8.7 vidíme výsledky pro 1 000 testů o rozsahu n = 100, v tabulce 8.8 jsou výsledky pro 1 000 testů o rozsahu n = 200. U výběrů o rozsahu n = 100 vidíme, že každý test zamítl přibližně polovinu náhodných výběrů, nejvíce jich zamítly SF, AG a JB, které mají zárovň i největší podíl zamítnutých výběrů s odlehlým pozorováním oproti výběrům bez odlehlého pozorování. Dále vidíme, že v tomto případě má velký vliv na výsledek odlehlé pozorování, neboť u každého testu nebyla u více než poloviny náhodných výběrů po odstranění odlehlého pozorování zamítnuta H0 . U výběrů o rozsahu n = 200 40

je rozpětí zamítnutých výběrů širší, od 577 u χ2 testu až po 816 u AG. Opět je vidět, že v tomto případě jsou nejsilnější AG, JB a SF, které jsou citlivé na odlehlá pozorování. Počet testů, kdy nebyla zamítnuta H0 , se po odstranění odlehlého pozorování snižuje, což může být způsobeno existencí druhého odlehlého pozorování. Závěr: nejsilnější se u normálního rozdělení kontaminovaného Cauchyovým jeví AG, JB a SF, které nejvíce reagují na odlehlé pozorování. Nejslabší je χ2 test následován LI. V tomto případě je vidět velmi velký vliv existence odlehlého pozorování na výsledek. Tabulka 8.7: Výsledky 1 000 testů náhodného výběru z N(0, 1) kontaminovaného Cauchyovým rozdělením s pravděpodobností p = 0,05 Test 1 SF 563 DA 448 AG 551 JB 556 BS 499 LI 424 CM 450 AD 475 χ2 386

2 547 448 544 552 481 400 427 453 360

3 16 0 7 4 18 24 23 22 26

4 219 122 206 217 173 117 152 164 92

5 328 326 338 335 308 283 275 285 268

6 437 552 449 444 501 576 500 525 614

7 56 155 59 51 122 203 176 150 243

8 9 381 0 397 13 390 4 393 1 379 1 373 6 374 3 375 1 371 5

10 56 142 55 50 121 197 173 149 238

Tabulka 8.8: Výsledky 1 000 testů náhodného výběru z N(0, 1) kontaminovaného Cauchyovým rozdělením s pravděpodobností p = 0,05 Test 1 SF 807 DA 681 AG 816 JB 814 BS 745 LI 647 CM 680 AD 710 2 χ 577

8.7

2 3 800 7 681 0 808 8 810 4 736 9 638 9 669 11 699 11 566 11

4 502 311 475 500 401 276 327 362 211

5 298 370 333 310 335 362 342 337 355

6 193 319 184 186 255 353 320 290 423

7 46 165 38 36 110 208 177 147 280

8 9 10 147 0 46 154 28 37 146 1 37 150 0 36 145 3 107 145 7 201 143 2 175 143 2 145 143 7 273

Shrnutí

Po provedení série testů bylo zjištěno, že SF, JB, BS, CM a AD patří mezi nejsilnější testy z těch, na kterých byla provedena studie. Výběry z normálního rozdělení zamítaly podle předpokladů, u výběrů z exponenciálního a Cauchyova rozdělení byly neomylné, u dalších tří rozdělení patřily mezi nejsilnější. Zhodnotit, který z těchto testů je nejsilnější, je velmi těžké, neboť počet různých rozdělení, u kterých byla provedena studie, byl malý a navíc každý test byl silnější u jiného 41

rozdělení a někde byl slabší, např. BS u exponenciálního rozdělení, JB u rovnoměrného. LI byl u výběrů s menším počtem pozorování, zejména u testů na rozdělení, které mají obdobný tvar distribuční funkce, slabší než výše uvedená pětice testů, ale jeho výsledky se dají označit za uspokojivé. Jako nejslabší dopadl χ2 test, ale je těžké zjistit, jak velký vliv na to měly odhady střední hodnoty a rozptylu, které zvyšují pravděpodobnost chyby prvního druhu. DA a AG testy byly uvedeny spíše pro úplnost. Porovnávat jejich sílu s jinými testy není možné, neboť testují šikmost nebo špičatost. Dají se použít pouze v případě, když výběr pochází z nesymetrického rozdělení nebo z rozdělení s výrazně lehčími nebo těžšími chvosty. V tomto případě je lepší místo nich použít JB, který testuje šikmost i špičatost zároveň. Testy na výběry z kontaminovaného rozdělení prokázaly, že na zamítnutí nebo nezamítnutí H0 může mít velký vliv odstranění odlehlého pozorování. Z tohoto důvodu je vhodné zjistit před případným odstraněním příčinu výskytu (viz kapitola 7). Ze všech testů nejsilněji na odlehlé pozorování reagují SF a JB, u kterých je zamítáno více výběrů s odlehlým pozorováním než výběrů bez odlehlého pozorování. JB test má největší počet testů, u kterých po odstranění odlehlého pozorování H0 zamítnuta nebyla, ačkoli před odstraněním byla zamítnuta.

8.8

Komentář k Softwarovému kódu

V příloze najdeme softwarový kód, který byl použit pro simulační studii. Byla použita hodnota set.seed(7777) pro simulaci výběrů o rozsahu n = 100 a hodnota set.seed(8888) pro simulaci výběrů o rozsahu n = 200. V kódu nahradíme text name.test příkazem pro příslušný test,12 který chceme použít (5-krát). Dále za YY dosadíme rozsah výběru. Pro generování ostatních náhodných výběrů nahradíme příkaz rbinom(YY, 1, 0.05) příkazem rnorm(YY), rexp(YY), rcauchy(YY), runif(YY, -1, 1) nebo rlaplace(YY), který najdeme v balíčku VGAM, a odstraníme následující dva řádky. Význam jednotlivých proměnných je popsán v následující tabulce, kde OOP znamená odstranění odlehlého pozorování. Tabulka 8.9: Význam jednotlivých proměnných Proměnná test p test i test o test io test n test in test op test oi test np test ni

Význam Počet testů, kdy zamítáme H0 Indexy testů, kdy zamítáme H0 Počet testů, kdy zamítáme H0 a existuje odlehlé pozorování Indexy testů, kdy zamítáme H0 a existuje odlehlé pozorování Počet testů, kdy nezamítáme H0 a existuje odlehlé pozorování Indexy testů, kdy nezamítáme H0 a existuje odlehlé pozorování Počet testů z test o, kdy zamítáme H0 po OOP Indexy testů z test io, kdy zamítáme H0 po OOP Počet testů z test n, kdy zamítáme H0 po OOP Indexy testů z test in, kdy zamítáme H0 po OOP

1

Jednotlivé příkazy byly popsány u každého testu. U Cauchyho, exponenciálního a rovnoměrného rozdělení byla u CM a AD srovnávána testová statistika s kritickou hodnotou. 2

42

Závěr Práce je zaměřena zejména na popisnou stránku testů, kterými můžeme ověřovat normalitu. V poslední části byla provedena simulační studie na vygenerovaných datech, která měla za cíl ukázat, že síla jednotlivých testů se může lišit v závislosti na tom, jaké je rozdělení náhodného výběru. Na základě této simulace nelze objektivně zjistit, který test je nejsilnější, neboť počet výběrů, na kterých byla studie provedena, by musel být větší. Obdobně by se musela studie provést i na výběrech o jiném rozsahu a výběry by měly pocházet z více rozdělení. Čtenář by se měl v práci seznámit s testovými i grafickými možnostmi, jak lze normalitu otestovat. Provedené simulace měly za cíl ukázat, že nelze používat výhradně jeden test. Před rozhodnutím, který z testů použijeme, je rozumné zvážit, jakou přednost by měl test splňovat. Je vhodné použít i grafické metody, abychom věděli, jaký výsledek očekávat. Výpočty probíhaly ve statistickém programu R, který mi připadal nejvhodnější, protože je veřejně dostupný a zároveň je k dispozici velké množství balíčků, které neustále přibývají.

43

Seznam použité literatury [1] Thode, H. C. (2002). Testing for normality. Marcel Dekker, New York. ISBN 0-8247-9613-6. [2] Anděl, J. (2005). Základy matematické statistiky. Matfyzpress, Praha. ISBN 80-86732-40-1. [3] R Development Core Team (2011). R: A language and environment for statistical computing. R Foundation for Statistical Computing. Vienna, Austria. ISBN 3-900051-07-0, URL http://www.R-project.org/. [4] Venables, V. N., Ripley, B. D. (2002). Modern applied statistics with S. Fourth edition. Springer, New York. ISBN 0-387-95457-0. [5] Jaroš, F., Rosa, Z. (1980). ČSN 01 0225 - Aplikovaná statistika. Testy shody empirického rozdělení s teoretickým. Český normalizační institut, Praha. [6] Sturges, H. A. (1926). The choice of a class interval. Journal of the American Statistical Association, 21(153), 65–66. [7] Scott, D. W. (1979). On optimal and data-based histograms. Biometrika, 66(3), 605–610. [8] Freedman, D., Diaconis P. (1981). On the histogram as a density estimator: L2 theory. Zeitschrift für Wahrscheinlichkeitstheorie und verwandte Gebiete, 57(4), 453–476. [9] Shapiro, S. S., Wilk M. B. (1965). An analysis of variance test for normality (complete samples). Biometrika, 52(3/4), 591–611. [10] Royston, J. P. (1982). An extension of Shapiro and Wilk’s W test for normality to large samples. Journal of the Royal Statistical Society, Series C (Applied Statistics), 31(2), 115–124. [11] Shapiro, S. S., Francia, R. S. (1972). An approximate analysis of variance test for normality. Journal of the American Statistical Association, 67(337), 215–216. [12] Gupta, A. K. (1952). Estimation of the mean and standard deviation of a normal population from a censored sample. Biometrika, 39(3/4), 260–273. [13] Filliben, J. J. (1975). The probability plot correlation coefficient test for normality. Technometrics, 17(1), 111–117. [14] Gan, F. F., Koehler, K. J. (1990). Goodness-of-fit tests based on P-P probability plots. Technometrics, 32(3), 289–303. [15] LaBrecque, J. (1977). Goodness-of-fit tests based on nonlinearity in probability plots. Technometrics, 19(3), 293–306.

44

[16] Ord, J. K. (1968). The discrete Student’s t distribution. The Annals of Mathematical Statistic, 39(5), 1513–1516. [17] Johnson, M. E., Tietjen, G. L., Beckman, R. J. (1980). A new family of probability distributions with applications to Monte Carlo studies. Journal of the American Statistical Association, 75(370), 276–279. [18] Balanda, K. P., MacGillivrayová, H. L. (1988). Kurtosis: a critical review. The American Statistician, 42(2), 111–119. [19] Pearson, E. S. (1930). A further development of tests for normality. Biometrika, 22(1/2), 239–249. [20] Pearson, E. S. (1931). Note on test for normality. Biometrika, 22(3/4), 423–424. [21] D’Agostino, R. B. (1970). Transformation to normality of the null distribution of g1 . Biometrika, 57(3), 679–681. [22] Anscombe, F. J., Glynn, W. J. (1983). Distribution of the kurtosis statistic b2 for normal samples. Biometrika, 70(1), 227–234. [23] Jarque, C. M., Bera, A. K. (1980). Efficient tests for normality, homoscedasticity and serial independence of regression residuals. Economics Letters, 6(3), 255–259. [24] Bowman, K. O., Shenton, L. R. √ (1975). Omnibus test contours for departures from normality based on b1 and b2 . Biometrika, 62(2), 243–250. [25] Pearson, E. S., D’Agostino, R. B., Bowman, K. O. (1977). Tests for departure from normality: Comparison of power. Biometrika, 64(2), 231–246. [26] Geary, R. C. (1935). The ratio of the mean deviation to the standard deviation as a test of normality. Biometrika, 27(3/4), 310–332. [27] Geary, R. C. (1936). Moments of the ratio of the mean deviation to the standard deviation for normal samples. Biometrika, 28(3/4), 295–307. [28] Geary, R. C. (1935). Note on the correlation between β2 and w0 . Biometrika, 27(3/4), 353–355. [29] D’Agostino, R. B. (1970). Simple compact portable test of normalty: Geary’s test revisited. Psychological Bulletin, 74(2), 138–140. [30] Bonett, D. G., Seierová E. (2002). A test of normality with high uniform power. Computational Statistics and Data Analysis, 40(3), 435–445. [31] Lilliefors, H. W. (1967). On the Kolmogorov-Smirnov test for normality with mean and variance unknown. Journal of the American Statistical Association, 62(318), 399–402. [32] Lilliefors, H. W. (1969). Corrigenda: On the Kolmogorov-Smirnov test for normality with mean and variance unknown. Journal of the American Statistical Association, 64(328), 1702. 45

[33] Anderson, T. W., Darling, D. A. (1952). Asymptotic theory of certain ”goodness of fit”criteria based on stochastic processes. The Annals of Mathematical Statistics, 23(2), 193–212. [34] Watson, G. S. (1961). Goodness-of-fit tests on a circle. Biometrika, 48(1/2), 109–114. [35] Anderson, T. W., Darling, D. A. (1954). A test of goodness of fit. Journal of the American Statistical Association, 49(268), 765–769. [36] Stephens, M. A. (1974). EDF statistics for goodness of fit and some comparisons. Journal of the American Statistical Association, 69(347), 730–737. [37] Pearson, K. (1900). On the criterion that a given system of deviations from the probable in the case of a correlated system of variables is such that it can be reasonably supposed to have arisen from random sampling. Philosophical Magazine, 50, 157–175. [38] Mann, H. B., Wald, A. (1942). On the choice of the number of class intervals in the application of the chi square test. The Annals of Mathematical Statistics, 13(3), 306–317. [39] Koehler, K. J., Larntz, K. (1980). An empirical investigation of goodness-of-fit statistics for sparse multinomials. Journal of the American Statistical Association, 75(370), 336–344. [40] Yarnold, J. K. (1970). The minimum expectation in χ2 goodness of fit tests and the accuracy of approximations for the null distribution. Journal of the American Statistical Association, 65(330), 864–886. [41] Eaton, P. W. (1978). Yarnold’s criterion and minimum sample size. The American Statistician, 32(3), 102–103. [42] David, H. A., Hartley, H. O., Pearson, E. S. (1954). The distribution of the ratio, in a single normal sample, of range to standard deviation. Biometrika, 41(3/4), 482–493. [43] Uthoff, V. A. (1973). The most powerful scale and location invariant test of the Normal versus the double exponential. The Annals of Statistics, 1(1), 170–174. [44] D’Agostino, R. B. (1971). An omnibus test of normality for moderate and large size samples. Biometrika, 58(2), 341–348. [45] D’Agostino, R. B. (1972). Small sample probability points for the D test of normality. Biometrika, 59(1), 219–221. [46] Oja, H. (1981). Two location and scale-free goodness-of-fit tests. Biometrika, 68(3), 637–640. [47] Vašíček, O. (1976). A test for normality based on sample entropy. Journal of the Royal Statistical Society. Series B, 38(1), 54–59.

46

[48] Grubbs, F. E. (1950). Sample criteria for testing outlying observations. The Annals of Mathematical Statistics, 21(1), 27–58. [49] Dixon, W. J. (1950). Analysis of extreme values. The Annals of Mathematical Statistics, 21(4), 488–506. [50] Dixon, W. J. (1951). Ratios involving extreme values. The Annals of Mathematical Statistics, 22(1), 68–78. [51] Tietjen, G. L., Moore, R. H. (1972). Some Grubbs-type statistics for the detection of several outliers. Technometrics, 14(3), 583–597. [52] http://en.wikipedia.org/.

47