Bootstrap replikáty průměru a SEM
V tomto cvičení vypočítáš bootstrap odhad funkce hustoty pravděpodobnosti průměrných ročních srážek na meteorologické stanici Sheffield. Odhadujeme průměrné roční srážky, které bychom dostali, kdyby stanice Sheffield mohla donekonečna opakovat všechna měření z let 1883–2015. Jde o pravděpodobnostní odhad průměru. Výsledné rozdělení zobrazíš jako histogram a uvidíš, že má normální tvar.
Dá se teoreticky dokázat, že za ne příliš omezujících podmínek má hodnota průměru vždy normální rozdělení. (To neplatí obecně, pouze pro průměr a několik dalších statistik.) Směrodatná odchylka tohoto rozdělení se nazývá střední chyba průměru (SEM) a rovná se směrodatné odchylce dat vydělené druhou odmocninou počtu datových bodů. Tedy pro daný dataset platí: sem = np.std(data) / np.sqrt(len(data)). Pomocí hacker statistics dostaneš stejný výsledek bez nutnosti ho odvozovat — a ověříš ho právě na svých bootstrap replikátech.
Dataset je pro tebe předem načten do pole rainfall.
Toto cvičení je součástí kurzu
Statistical Thinking in Python (Part 2)
Pokyny k cvičení
- Pomocí funkce
draw_bs_reps()a polerainfallvygeneruj10000bootstrap replikátů průměru ročních srážek. Nápověda: Pro výpočet průměru předej jakofunchodnotunp.mean.- Připomínka:
draw_bs_reps()přijímá 3 argumenty:data,funcasize.
- Připomínka:
- Vypočítej a vypiš střední chybu průměru hodnot
rainfall.- Vzorec pro výpočet:
np.std(data) / np.sqrt(len(data)).
- Vzorec pro výpočet:
- Vypočítej a vypiš směrodatnou odchylku svých bootstrap replikátů
bs_replicates. - Vytvoř histogram replikátů s argumentem
density=Truea50sloupci. - Klikni na tlačítko Odeslat a podívej se na výsledný graf!
Interaktivní cvičení na vyzkoušení si v praxi
Vyzkoušejte si toto cvičení dokončením tohoto ukázkového kódu.
# Take 10,000 bootstrap replicates of the mean: bs_replicates
bs_replicates = ____
# Compute and print SEM
sem = ____ / np.sqrt(____)
print(sem)
# Compute and print standard deviation of bootstrap replicates
bs_std = ____
print(bs_std)
# Make a histogram of the results
_ = plt.hist(____, ____=50, ____=True)
_ = plt.xlabel('mean annual rainfall (mm)')
_ = plt.ylabel('PDF')
# Show the plot
plt.show()