Centrální limitní věta (Cetral Limit Theorem – CLT)

“Rozdělení součtu mnoha vzájemně nezávislých veličin konverguje k normálnímu rozdělení.”

Pozn. Přičemž jednotlivé dílčí veličiny mohou mít libovolné rozdělení.

Galton’s quincunx: kuličkový stroj demonstrující CLT

Ukažme to na praktickém příkladu. Kuličky se náhodně odrážejí od mnoha překážek, důsledek každého odrazu je dán tolika příčinami, do kterých pořádně nevidíme, že je pro nás výsledek zkrátka náhodný (ať už má statistické rozdělení jakékoliv). A tyto dílčí výsledky se sčítají, celkové rozdělení končí téměř normálním rozdělením.

https://www.youtube.com/watch?v=1DTRzPRfu6s

https://www.youtube.com/watch?v=6YDHBFVIvIs

Sčítání proměnných s rovnoměrným rozdělením

Jedna proměnná

Vyrobíme 10 tisíc realizací proměnné s rovnoměrným rozdělením (v rozsahu 0 až 1). Histogram by měl odpovídat hustotě pravděpodobnosti rovnoměrného rozdělení.

x1 <- runif(10000)
hist(x1, 30)

Dvě proměnné

Co když sečteme výsledky dvou proměnných s rovnoměrným rozdělením x1 + x2? To bude opět jedno číslo, ale bude mít stále rovnoměrné rozdělení? Zkusíme vygenerovat opět 10 tisíc součtů takovýchto dvojic.

x1 <- runif(10000)
x2 <- runif(10000)
hist(x1+x2, 30)

Ó jé, co se nám to děje? Kdo myslel, že dostaneme opět rovnoměrné rozdělení, akorát možná s dvojnásobným rozsahem hodnot, je jistě překvapen. Jak je možné, že to nyní vypadá jako trojúhelník?

Zkusme úvahu selským rozumem. Kolika způsoby můžeme získat ze součtu dvou proměnných (s rozsahem 0 až 1) výsledný součet 0? Jedině jedinou kombinací 0 + 0. Musí padnout dvě nuly, což je velice nízká pravděpodobnost.

A kolika způsoby dostaneme druhý okraj rovný 2? Jedině 1 + 1, opět mizivá pravděpodobnost.

A kolika způsoby se můžeme dostat k prostřední hodnotě 1? Tak to máme buď 0 + 1 a nebo 1 + 0, ale také 0.5 + 0.5 nebo 0.4 + 0.6 či 0.6 + 0.4, samozřejmě též 0.1 + 0.9 a spoustou dalších variant. Přestože pravděpodobnost každé jednotlivé vyjmenované možnosti je také nepatrná, těch možností je evidentně mnohem více, než pro výše zmiňované součty 0 nebo 2. Takže histogram smysl dává.

Tři proměnné

Jak to dopadne pro součet tří proměnných x1 + x2 + x3? Opět zkusme 10 tisíc krát.

x1 <- runif(10000)
x2 <- runif(10000)
x3 <- runif(10000)
hist(x1+x2+x3, 30)

Deset proměnných

Toto už přeci nebudeme dělat ručně.

k <- 10     # počet sčítaných proměnných
n <- 10000  # počet generovaných součtů

cisla <- matrix(runif(n*k), nrow = k) # vygenerujeme n*k čísel a uložíme je postupně do matice
                                      # s 'k' řádky. Každý řádek tak bude odpovídat 'n' pokusům
                                      # jedné proměnné s rovnoměrným rozdělením.

soucty <- colSums(cisla)  # a sečteme vše po sloupcích - každý sloupec odpovídá součtu 'k' proměnných
                          # v rámci jednoho pokusu

hist(soucty, 30)

To už vypadá hodně jako normální rozdělení.

Interaktivní experiment

Pokud máte RStudio, můžete si zkusit zkopírovat následující příklad. Vyrobí posuvník, kterým můžete měnit počet sčítaných proměnných a sledovat výsledný histogram.

library(manipulate) # pokud není nainstalována, tak nejdříve spustit
                    # install.packages("manipulate") - vyžaduje internet

mujHistogram <- function(k) {  # k ... počet sčítaných proměnných
    n <- 10000  # počet generovaných součtů
    
    cisla <- matrix(runif(n*k), nrow = k)
    soucty <- colSums(cisla)
    hist(soucty, 30, col = "lightblue")
}

manipulate(mujHistogram(k), k = slider(1, 100, 1, step = 1))

Pozn. Už od velice malého počtu sčítaných proměnných to vypadá jako normální rozdělení. Spíše budete pozorovat, jak se s námi funkce pro histogram v R “pere” a nastavuje si počet sloupců podle sebe. Parametr 30 je podle nápovědy skutečně jen doporučený a ona si to stejně udělá po svém.

Sčítání hodů kostkou

Tak dobře, s výše uvedeným by se ještě nějak smířit dalo. Součty mnoha spojitých veličin vedou na spojité normální rozdělení, protože prostředek lze vyrobit více variantami hodnot než kraje.

Ale v úvodu bylo řečeno, že se to týká jakýchkoliv rozdělení. Tedy i diskrétních?

library(manipulate) # pokud není nainstalována, tak nejdříve spustit
                    # install.packages("manipulate") - vyžaduje internet

mujHistogram <- function(k) {  # k ... počet sčítaných proměnných
    n <- 10000  # počet generovaných součtů
    
    hodyKostkou <- matrix( sample(1:6, n*k, replace = TRUE) , nrow = k)
    soucty <- colSums(hodyKostkou)
    barplot(table(soucty), col = "lightblue") # jsou to diskrétní kategorie, lépe vykreslovat takto
}

manipulate(mujHistogram(k), k = slider(1, 100, 1, step = 1))

Že by ti matematici měli přeci jen pravdu?

Sčítání hodů mincí

Zkusme to s mincí. Rubu přiřadíme číslo 0, líci 1.

library(manipulate) # pokud není nainstalována, tak nejdříve spustit
                    # install.packages("manipulate") - vyžaduje internet

mujHistogram <- function(k) {  # k ... počet sčítaných proměnných
    n <- 10000  # počet generovaných součtů
    
    hodyMinci <- matrix( sample(c(0, 1), n*k, replace = TRUE) , nrow = k)
    soucty <- colSums(hodyMinci)
    barplot(table(soucty), col = "lightblue")
}

manipulate(mujHistogram(k), k = slider(1, 100, 1, step = 1))

Hm. A kolik že ta CLV stojí, že bychom ji na Ústav také zakoupili?

Úlohy

  1. Zkuste spustit uvedené interaktivní experimenty v RStudiu.

© 7. 3. 2015 Tomáš Bořil, borilt@gmail.com