Rozdělení náhodných veličin

d Hustota pravděpodobnosti (probability density function – PDF)
p Distribuční funkce (cumulative distribution function – CDF)
q Kvantily – inverze distribuční funkce (quantiles – inverse CDF)
r Generátor náhodných čásel (random number generator)

Rozdělení	Hustota pravděpodobnosti	Distribuční funkce	Kvantily	Generátor
rovnoměrné spojité	dunif	punif	qunif	runif
normální	dnorm	pnorm	qnorm	rnorm
Studentovo	dt	pt	qt	rt
lognormální	dlnorm	plnorm	qlnorm	rlnorm
binomické	dbinom	pbinom	qbinom	rbinom
Poissonovo	dpois	ppois	qpois	rpois
chi-kvadrát	dchisq	pchisq	qchisq	rchisq
exponenciální	dexp	pexp	qexp	rexp

Další rozdělení: beta (beta), F-rozdělení (f), gamma (gamma), geometrické (geom), hypergeometrické (hyper), negativně binomické (nmbinom), Weibullovo (weibull).

Příklady vykreslení hustoty pravděpodobnosti a distribuční funkce

Normální rozdělení

par(mfrow = c(1, 2))  # Kreslení více obrázků: 1 řádek a 2 sloupce

p = seq(from = 40, to = 160, length.out = 1000)  # rozsah pravděpodobností (s jemností 1000 bodů)

plot(p, dnorm(x = p, mean = 100, sd = 15), main = "Hustota pravděpodobnosti",   # normální rozdělení hodnot IQ
     type = "l", frame = FALSE)                                      # volíme typ "spojitá čára" a bez rámečku
plot(p, pnorm(q = p, mean = 100, sd = 15), main = "Distribuční funkce",
     type = "l", frame = FALSE)

Rovnoměrné rozdělení

par(mfrow = c(1, 2))

p = seq(from = -10, to = 70, length.out = 1000)  # rozsah pravděpodobností (s jemností 1000 bodů)

plot(p, dunif(x = p, min = 0, max = 60), main = "Hustota pravděpodobnosti",   # rovnoměrné rozdělení od 0 do 60
     type = "l", frame = FALSE)
plot(p, punif(q = p, min = 0, max = 60), main = "Distribuční funkce",
     type = "l", frame = FALSE)

Kvantily standardního normálního rozdělení

q1 <- qnorm(0.025, mean = 0, sd = 1)  # 2.5% kvantil
q2 <- qnorm(0.975, mean = 0, sd = 1)  # 97.5% kvantil
print(q1)

## [1] -1.959964

print(q2)

## [1] 1.959964

p <- seq(from = -5, to = 5, length.out = 1000)
plot(p, dnorm(p), type = "l", frame = FALSE)

p95 <- seq(from = q1, to = q2, length.out = 200)  # rozsah hodnot mezi oběma kvantily
polygon(c(q2, q1, p95), c(0, 0, dnorm(p95)), col = "lightblue")
text(0, 0.15, "95 %")    # text na uvedených souřadnicích

Binomické rozdělení

Opakované pokusy (1 pokus má dvě možnosti výsledku).

Př. Pravděpodobnost, že se narodí holka, je 50 %. V rodině je 8 dětí.

Jaká je pravděpodobnost, že je méně nebo rovno 6 dívek?

pbinom(q=6, size=8, prob=.5) # pravděpodobnost "hodnoty méně nebo rovno" => z klasické distribuční funkce

## [1] 0.9648438

Pozn. je použit defaultní lower.tail=TRUE, kdy je otázka definována Pravděpodobnost[hodnota <= q].

Jaká je pravděpodobnost, že je více než 6 dívek (tedy 7 nebo 8)?

pbinom(q=6, size=8, prob=.5, lower.tail=FALSE) # opačná otázka => distribuční funkce pro opačnou otázku

## [1] 0.03515625

Pozn. u lower.tail=FALSE je otázka definována Pravděpodobnost[hodnota > q].

Pozn. součet možností a) a b) musí být 100 %.

Jaká je pravděpodobnost, že je přesně 6 dívek?

dbinom(x=6, size=8, prob=.5)  # pravděpodobnost přesné hodnoty => z hustoty pravděpodobnosti

## [1] 0.109375

Pozn. součet dílčích pravděpodobností pro x = 0 až 8 musí být 100 %.

sum(dbinom(x=0:8, size=8, prob=.5))

## [1] 1

Největší pravděpodobnost má situace, že se narodí přesně polovina dívek (čtyři), protože pravděpodobnost holky je 0.5 (nastavili jsme prob=.5). Viz maximum následujícího grafu.

plot(0:8, dbinom(x=0:8, size=8, prob=.5))

Př. Pořádáme kolaudaci nového bytu. Každý člověk má pravděpodobnost, že polije koberec nápojem, 3 procenta. Celkem přijde 20 lidí. Jaká je pravděpodobnost, že nikdo nepolije koberec?

dbinom(x = 0, size = 20, prob = 0.03) # pravděpodobnost konkrétního počtu => z hustoty pravděpodobnosti

## [1] 0.5437943

Př. Pravděpodobnost, že student napoprvé dostane zápočet ze statistiky, je 80 procent. Jaká je pravděpodobnost, že napoprvé projde všech 16 studentů?

dbinom(x = 16, size = 16, prob = 0.80) # pravděpodobnost konkrétního počtu => z hustoty pravděpodobnosti

## [1] 0.0281475

Jaká je pravděpodobnost, že napoprvé projde více než 8 studentů (tedy 9 až 16)?

# pravděpodobnost "hodnoty více než" => distribuční funkce s "higher tail"
pbinom(q = 8, size = 16, prob = 0.80, lower.tail = FALSE)

## [1] 0.9929964

Bonus – kolik pokusů provést, aby se něco alespoň jednou povedlo?

Na tuto úlohu přímo funkci nemáme, ale připadá mi natolik důležitá, že ji tu musím zmínit. Odvození vyplývá také z Binomického rozdělení, ale řešíme opačnou úlohu – odpověď na otázku, kolik opakování provést, abychom s určitou jistotou alespoň jednou uspěli.

Př. Pravděpodobnost, že vyklíčí semínko, je 40 %. Kolik jich zasadit, aby s 99% jistotou alespoň jedna rostlinka vyrostla?

Uvedeme rovnou “magický vzorec” – zájemci o odvození se na mě obraťte, mám ho ve svých přípravách.

minimalniPocet <- log(1 - 0.99) / log(1 - 0.4)
print(minimalniPocet)

## [1] 9.015151

Chce to tedy alespoň 9 semínek, aby se teoreticky v 99 % případů experiment zdařil.

Poissonovo rozdělení

Aproximace binomického rozdělení pro situace, kdy je size obrovské a prob dílčích částí malinké. Velice dobře se hodí pro modelování pozorovaných počtů za nějaké období. Na rozdíl od binomického rozdělení je “výpočet uvnitř” výrazně jednodušší (binomické totiž počítá s faktoriály), a proto lépe proveditelný právě pro velká čísla.

Př. Na opuštěné železniční zastávce poblíž málo navštěvované zříceniny se v turistické sezóně objeví průměrně 3 lidé za hodinu. Tuto situaci můžeme modelovat Poissonovým rozdělením, protože je to výsledek obrovského množství lidí, kteří všichni s velice malou pravděpodobností vyrazí právě sem na výlet. Nicméně ono “hodně krát málo” způsobí, že se zde průměrně pár lidí zjeví.

Jaká je pravděpodobnost, že se během 2 hodin na zastávce objeví max. 8 lidí?

ppois(q = 8, lambda = 3*2)

## [1] 0.8472375

Jaká je pravděpodobnost, že se během 2 hodin na zastávce objeví více než 8 lidí?

ppois(q = 8, lambda = 3*2, lower.tail = FALSE)

## [1] 0.1527625

Pozn. součet pravděpodobností z obou otázek musí dát pochopitelně 100 %.

Normální rozdělení

Př. Předpokládejme, že diastolický tlak pro muže ve věku 30 – 40 let je přibližně normálně rozdělen se středem 78 (mm Hg) a směrodatnou odchylkou 12. Jaká je pravděpodobnost, že náhodně vybraný muž v tomto věku má tlak nižší nebo roven 70?

pnorm(q = 70, mean = 78, sd = 12) # podíváme se do klasické distribuční funkce

## [1] 0.2524925

Jaká je pravděpodobnost, že bude mít tlak vyšší než 86?

# podíváme se do distribuční funkce řešící opačně položenou otázku
pnorm(q = 86, mean = 78, sd = 12, lower.tail = FALSE)

## [1] 0.2524925

Pozn. normální rozdělení je vůči středu symetrické, takže obě předchozí otázky vedou na stejnou pravděpodobnost.

Jaká je pravděpodobnost, že bude mít tlak nižší nebo roven 86?

pnorm(q = 86, mean = 78, sd = 12) # klasická distribuční funkce

## [1] 0.7475075

Pozn. toto musí být samozřejmě doplněk k předchozímu, takže součet celkem 100 %.

Jaká je pravděpodobnost, že bude mít tlak v rozmezí 54 až 102?

pnorm(q = 102, mean = 78, sd = 12) - pnorm(q = 54, mean = 78, sd = 12)

## [1] 0.9544997

Jaká je pravděpodobnost, že bude mít tlak úplně přesně 78?

pnorm(q = 78, mean = 78, sd = 12) - pnorm(q = 78, mean = 78, sd = 12)

## [1] 0

Jaká je pravděpodobnost, že bude mít tlak po zaokrouhlení na celé roven 78?

pnorm(q = 78.49999, mean = 78, sd = 12) - pnorm(q = 77.5, mean = 78, sd = 12)

## [1] 0.03323524

Př. Objem mozku dospělých žen se zhruba řídí normálním rozdělením se střední hodnotou 1100 cm^3 a směrodatnou odchylkou 75 cm^3. Určitá žena se dozvěděla, že její mozek odpovídá 95. percentilu (kvantil 95 %), neboli 95 % žen má menší objem mozku a jen 5 % žen má větší objem mozku. Jaký má ona objem?

qnorm(p = 0.95, mean = 1100, sd = 75)  # Hledáme hodnotu kvantilu 0.95.

## [1] 1223.364

Žena má objem mozku 1050 cm^3, kolikaprocentnímu kvantilu to odpovídá? Neboli: kolik % žen má nižší objem?

pnorm(q = 1050, mean = 1100, sd = 75) # podíváme se do distribuční funkce

## [1] 0.2524925

Zdá se, že je přesně na úrovni prvního kvartilu (kvantil 0.25 – 25 %, neboli 25. percentil).

Exponenciální rozdělení

Životnost mnoha výrobků (které nejsou zmetky hned od výroby) se řídí přibližně exponenciálním rozdělením. Pro tento předmět to nepovažuji za stěžejní, ale pro zájemce (budoucí manažery) uvádím příklad, jak s ním počítat.

Př. Kolik % se rozbije po 2 letech, když střední doba životnosti je 3 roky?

# Otázka: kolik %? Použijeme tedy distribuční funkci, rate = průměrná "frekvence" rozbití. Jestliže průměrné
# trvání je 3 roky, tak frekvence je převrácená hodnota.
pexp(q = 2, rate = 1/3)

## [1] 0.4865829

Pozn. zajímavé je, že zbylých 51 % přežije ještě déle než jednou tolik (je zde nesymetrie), proto je střední doba 3 roky větší číslo než první 2 roky, kdy se rozbije zhruba půlka zařízení.

Nyní se však nabízí otázka. Jsem výrobce a nechci čekat dlouhá léta, abych si udělal kompletní statistiku. Potřebuji nastavit cenu výrobku, abych měl peníze na případné reklamace. Jak odhadnout střední dobu životnosti, když teď je čas \(t\) a už se jich rozbilo \(r\) procent?

Vyloučím ty, které se rozbily hned (zmetky při výrobě).
Střední doba životnosti = \(-t/log(1-r)\)
Pak můžu počítat jiná \(r\) pro jiná \(t\) :-)

Vztahy dávám k dispozici bez odvození (opět: zájemci, ozvěte se). Dosazení je naštěstí jednoduché.

Např. Po 1 roce se rozbilo 20 % výrobků.

stredniDoba <-  -1/log(1 - 0.2)

Kolik % se rozbije do dvou let (od počátku)?

pexp(q = 2, rate = 1/stredniDoba)

## [1] 0.36

Po dvou letech se rozbije asi 36 %.

Podle tohoto procenta navýším cenu výrobku a mám vyhráno :-) Pozn. samozřejmě by to také chtělo uvážit daně, fixní náklady a další vychytávky.

Pozn. v R znamená log() přirozený logaritmus, což je pro tyto výpočty správně.

Úlohy

Jako učitel máte na starost třídenní výlet s 25 dětmi. Pravděpodobnost, že si jedno dítě zapomene zabalit s sebou zubní kartáček, je 5 %. Jaká je pravděpodobnost, že všechny děti budou mít zubní kartáčky? [Odpověď: zhruba 28 %]
Jaká je pravděpodobnost, že při 7 hodech mincí padne minimálně pětkrát rub? [Odpověď: cca 23 %]
Ve volební komisi jsou zvyklí, že se při sčítání průměrně objevují 3 neplatné lístky. Jaká je pravděpodobnost, že se jich objeví 5 a více? Pozn. U každé osoby je pravděpodobnost omylu velice nízká, nicméně v tom velkém množství se zkrátka skoro vždycky někdo najde. [Odpověď: asi 18 %]
V kadeřnictví se z hlediska dlouhodobého průměru objevují za hodinu 4 lidé. Pozorujeme situaci 3 hodiny. Jaká je pravděpodobnost, že se v tomto časovém úseku objeví 8 nebo méně lidí? [Odpověď: přibližně 16 %]
Počet dosažených bodů v testu se u studentů gymnázií řídí přibližně normálním rozdělením se střední hodnotou 60 a směrodatnou odchylkou 20. Určitý student se dozvěděl, že jeho výsledek odpovídá 92. percentilu, je tedy lepší než 92 % ostatních studentů, ale 8 % jiných studentů dosáhlo ještě lepší výsledky. Kolik bodů obdržel? [Odpověď: 88 bodů]
Mějme rozdělení z předchozí otázky. Jistý student obdržel 53 bodů. Jakému percentilu to odpovídá? [Odpověď: 36. percentil]
Denní návštěvnost webových stránek se řídí zhruba normálním rozdělením se střední hodnotou 3000 a směrodatnou odchylkou 1000. Jaká je pravděpodobnost, že během jednoho dne navštíví stránku více než 4000 lidí? [Odpověď: okolo 16 %]
Vygenerujte 10000 náhodných čísel a zobrazte histogram (s cca 20 hranicemi intervalů)
- pro normální rozdělení modelující hodnoty IQ (střední hodnota 100, směrodatná odchylka 15),
- pro rovnoměrné rozdělení v rozsahu hodnot 0 až 60.
Jestliže je IQ modelováno normálním rozdělením dle parametrů z předchozí otázky, kolik procent populace má IQ vyšší než 125? [Odpověď: zhruba 5 %]
Doplnění k předchozí otázce: je 5 % málo? [Odpověď: každý dvacátý není zrovna málo, v naší republice to odpovídá zhruba půl milionu lidí.]