Cvičenie 7 (praktická časť 1)

(Základy práce so štatistickým programom R)

Zdrojový soubor pro Knit: Rmd soubor

I. Štatistický program R (stručný úvod)

Štatistický program R (R Core Team, 2020) je výpočetný nástroj určený k práci s datami rôznych typov – od prípravy dat, cez štatistické spracovaníe (analýzu), až po následné reportovanie získaných výsledkov (vo forme nejakej výskumnej správy). Jedná sa o GNU projekt, podobný programu S a primárne je navrhnutý pre štatistickú analýzu dat – praktickú aplikáciu rôznych stochastických (pravdepodobnostných a štatistických) metód.
Program R (dostupný pod GNU GPL licenciou) je k dispozícii k stiahnutiu (free of charge) na adrese

https://www.r-project.org

K dispozícii sú distribúcie s priamou podporou pre OS Windows, Linux aj Macintosh.

Základnú inštaláciu programu R je možné následne rozšíriť pomocou dodatočných knižníc (balíčkov), ktoré sú k dispozícii na rôznych online repozitároch (zoznam hlavných repozitárov je na adrese https://cran.r-project.org/mirrors.html). Jednotlivé R knižnice sú ale tvorené samotnými užívateľmi softwaru R a ich správne fungovanie nie je garantované – je preto namieste určitá opatrnosť a hlavne aktívne a kritické premýšľanie pri ich používaní.

Pre užívateľov programu R sú k dispozícii aj rôzne grafické rozhrania, ktoré je možne dodatočne nainštalovať a umožňujú (v určitých smeroch) jednoduchšiu a prehľadnejšiu prácu. Najznámejší a pravdepodobne aj jeden z najčastejšie používaných R interfacov je RStudio (k dizpozícii pre OS Windows, Linux, aj Macintosh).

Užitočné doplňujúce materiály pre prácu so štatistickým programom R

Bína, V., Komárek, A. a Komárková, L.: Jak na jazyk R. (PDF súbor)
Komárek, A.: Základy práce s R. (PDF súbor)
Kulich, M.: Velmi stručný úvod do R. (PDF súbor)
De Vries, A. a Meys, J.: R for Dummies. ISBN-13: 978-1119055808 (link na PDF súbor)

Na následujúcich riadkoch je uvedených niekoľko jednoduchých príkladov a úloh, ktoré majú primárne slúžiť ako názorné ilustrácie fungovania programu R (implementácia jednotlivých príkazov a práca s vlastnými datovými súbormi). Nejedná sa ale o kompletný súhrn možnosti, ktoré sú v programe R k dispozícii. Zdrojový kod v šedých boxoch je možné priamo kopírovať do príkazovej konzoly programu R.

Je možné využiť napr. podrobnejšieho “sprievodcu” na cestu základmi Rka: Hrátky s R (autor: doc. Arnošt Komárek);
(k Hrátkam s R sú potrebné datove súbory: auta2004.dat, auta2004.csv a auta2004.xls)

1. Základné matematické operácie v R

Najjednoduchšie využitie programu R ako jednoduchá, ale sofistikovaná kalkulačka…

2 + (4 * 5)^2 ## za '#' nasleduje komentar, ktory R ignoruje
sqrt(256) ## druha odmocnina
1:5 * 2 ## skalarny sucin jednoducheho vektoru a skalaru
1 + 2 * (2 * (2 - 1) + 2)

Výsledky je možné uložiť do vhodne zvoleného objektu, na ktorý sa následne dá odvolávať/odkazovať.

v1 <- 2 + (4 * 5)^2 ## za '#' nasleduje komentar, ktory R ignoruje
v2 <- 1:5 * 2
print(v1)

Uložené objekty (výsledky, alebo aj komplexnejšie datové štruktúry ako vektor, matica, pole, list) je možné vhodné kombinovať a ďalej s nimi pracovať.

v1 + v2
v3 <- log(v1 + v2)
v4 <- v3 - 2
v5 <- sin(pi * v2)

Užitočné

Program R je objektový nástroj a zvláda aj matematické operacie, ktoré štandardne nie sú definované. Napr. súčet dvoch vektorov, ktoré majú rôznu dĺžku, operácie s maticami, ktoré majú rôzne typy/dimenzie. Je podstatné správne pochopiť, ako program R v takýchto prípadoch funguje. Zabráni sa možným problémom v budúcnosti a v mnohých prípadoch to zase na druhu stranu može zjednodušiť výpočet a ušetriť čas (viď niekoľko explicitných príkladov nižšie).

Štandardné operácie v programe R získame pomocou štandardných znamienok +, -, *, \, a symbolu ^ pro mocniny/odmocniny;
Komplexnejšie matematické výrazy môžeme formulovať pomocou zátvoriek. V Rku slúžia k tomuto účelu pouze okrúhle zátvorky (). Hranaté zátvorky [] a zátvorky {} maju v programe R vlastnú, špecifickú úlohu (viď nižšie);
Ku každému príkazu v R je k dispozícii návod, ako príkaz správne funguje. Príslušný help sa zavolá zadaním samotného príkazu, ktorému predchádza otazník: napr. ?c zobrazí help pre funkciu c(), analogicky ?seq zobrazí help k funkcii seq(). Help obsahuje popis samotnej funkcie, vysvetlenie implementácie a tiež niekoľko ilustratívnych príkladov a odkazov na dalšie podobné, alebo súvisejúce funkcie.
Podrobný návod je k dispozícii aj na internete a v prípade akéhokoľvek problému stačí použiť Google/ChatGPT a správne sa pýtať.

2. Práca s vektormi a maticami

Komplexnejšie datové štruktúry zahrňajú napr. matice, polia, listy…. V prípade, že máme v objektoch uložené informácie typu vektor/matica/pole/list, je možné sa v Rku odkazovať aj na jednotlivé elementy v daných objektoch s používať ich v dalších výpočtoch, alebo s nimi manipulovať podľa potreby. K tomuto účelu práve slúžia hranaté zátvorky [].

v3[1] ### zobrazi prvu zlozku vektora v3
v3[-c(2,3)] ### rovnaky vystup, vektor v3 bez druhej a tretej zlozky

A to isté platí aj pre matice (a analogicky aj pre pole).

m1 <- rbind(v3, v3 * 2, v3 * 3)
m1[1:2,1:2] ## elementy z prveho a druheho riadku a stlpca matice m1

Príkazy je možné vhodneˇ kombinovať a vytvárať nové, zložitejšie objekty/štruktúry:

m2 <- rbind(cbind(m1[1:2,1:2], c(0,0)), seq(10,11, length = 3))
l1 <- list(m1, m2) ## list, ktory obsahuje dve rozne matice, m1 a m2
l1[[2]][1,] # vypise na obrazovku prvy riadok matice m2

Objekty rôzneho typu môžu byť pohromade uložené napr. v objekte typu list()

li <- list(v3, m1) ## list, ktory obsahuje vektor v3 a maticu m1

Na jednotlivé elementy v liste sa odkazuje následovne:

li[[1]] ## zobrazi prvy element listu, vektor v3
li[[2]] ## zobrazi druhy element listu, maticu m1

Vždy je vhodné voliť kompaktnú a čo najjednoduchšiu formu zápisu. Pre porovnanie, všetky nasledujúce príkazy vytvoria vo výsledku jednotkovú maticu typu 4x4, ale najkompaktnejši a najstručnejší je pouze posledný zápis. Vyžaduje ale znalosť príkazu diag() a jeho spávnu implementáciu.

m1 <- rbind(c(1,0,0,0), c(0,1,0,0), c(0,0,1,0), c(0,0,0,1))

m2 <- cbind(c(1,0,0,0), c(0,1,0,0), c(0,0,1,0), c(0,0,0,1))

m3 <- matrix(rep(0,16), nrow = 4)
diag(m3) <- 1

m4 <- NULL
for (i in 1:4){
  m4 <- rbind(m4, c(rep(0, i - 1), 1, rep(0, 4 - i)))
}

m5 <- diag(1, 4, 4)

Užitočné

V programe R existuje niekoľko prikazov, ktoré vytvoria vektor: napr. c(), seq(), rep(), replicate() a pod.;
Podobne existuje niekoľko príkazov, ktoré vytvoria maticu: napr. cbind(), rbind(), matrix(), diag(), upper.tri(), lower.tri(), data.frame() a pod.;
Užitočné príkazy pre prácu s vektormi a maticami sú napr. as.numeric(), as.matrix(), array().
Pre prácu s objektom typu list() sú užitočné príkazy unlist(), as.list(), pairlist() a pod.
Pre maticové násobenie je nutné použiť operátor %*%, ktorý výnasobi dve matice (ak majú správne rozmery). Vhodnosť rozmerov sa jednoducho skontroluje pomocou príkazu dim() (užitočne je poznať tiež príkazy nrow(), ncol(), alebo dimnames()).
Inverznú maticu lze získať príkazom solve(), vlastne čísla/vektory matice pomocou príkazu eigen(), rozklady pomocou príkazov eiv(), svd(), qr(), alebo chol().
Štandardne binárne operátory +, -, *, /, ^ fungujú ‘po zložkách’ – toto ‘fungovanie’ je ale v programe R špecifické – porovnajte napr. nasledujúce výstupy (a vysvetlite spôsob fungovania jednotlivých maticových operácii):``
```
1:4 + 1 ### sucet po zlozkach
1:4 * 1:2 ### nasobenie po zlozkach v paroch (dlzka jedneho vektora je celociselny nasobok dlzky druheho)
1:4 / 1:4 ### nasobenie v paroch (rovnake dlzky oboch vektorov)
1:5 + 1:2 ## Error - preco?
```

3. Náhodné data – generátory (pseudo) náhodných čísel

Program R je primárne štatistický nástroj – umožňuje preto aj prácu s náhodnymi hodnotami. Implementované sú rôzne generátory (pseudo) náhodných čísel. Pre rôzne pravdepodobnostné rozdelenia sú k dispozícii príslušne generátory, ktoré simulujú náhodné hodnoty z daného pravdepodobnostného rozdelenia.

Diskrétne rozdelenia: rbinom(), rpois(), rgeom(), a ďalšie;
Spojité rozdelenia: runif(), rexp(), rnorm(), a ďalšie;

Pre správne použitie geneátorov je potrebné preštudovať prislušný návod – napr. ?rnorm zobrazí help pre generátor náhodných čísel s normálnym rozdelením. Je vhodné si tiež všimnúť, že existujú podobné príkazy, označené ako pnorm(), dnorm() a rnorm(), ktoré postupne značia distribučnú funkciu daného rozdelenia, hustotu a kvantilovú funkciu.

Užitočné

Je nutné dôsledne si preštudovať implementáciu jednotlivých generátorov a funkcií, ktoré súvisia s pravdepodobnostnými rozdeleniami – je totiž pomerne časté, že parametre, ktoré špecifikujú nejaké rozdelenie, majú iný význam. Napríklad:
```
x <- rnorm(10, 0, 2) ## simuluje 10 pozorovani z N(0, 4) - sigma^2 = 4
y <- rexp(10, 4) ## simuluje 10 pozorovani z Exp(lambda = 4), kde EX = 1/4
```
V prípade použítia nejakého generátora náhodných čísel je vhodné predom nastaviť tzv. seed. Umožní to následnu replikáciu generovania s garanciou rovnakých výsledkov.
```
set.seed(1234)
```
Výsledok použitia generátora náhodných čisel môže byť napríklad vektor, alebo matica (každý stĺpec ma stejný počet hodnôt), alebo list (umožňuje rôzny počet hodnôt v jednotlivých elemntoch).

set.seed(1234)
g1 <- rnorm(20, 160, 20) ## napr. vyska 20 nahodnych lidi s prumernou vyskou 160 cm a rozptylom 400

g2 <- cbind(g1, rexp(20, 1), rbinom(20, size = 1, prob = 1/2)) ## matica

g3 <- list(g1, rexp(10, 1), rbinom(50, size = 1, prob = 1/2)) ## list

4. Reálne data a vlastné datové súbory

V programe R je k dispozícii niekoľko vzorových ilustračných datových súborov. Ich zoznam je možne zobraziť zavolaním príkazu data(). Konkrétny datovy súbor je následne zobrazený zavolaním mena prislušného súboru z príkazoveho riadku v R – napr. Orange, alebo mtcars. Vo väčšíne prípadov sa jedná o objekt typu data.frame – čo ješpeciálny typ matice (pozri help ?data.frame). Na rozdiel od klasickej matice umožňuje type data.frame()použiť rôzne datové typy (integer, boolean, character, etc.) v jednotlivých stĺpcoch/premenných. Je tak možné vytvoriť datový súbor, ktorý bude obsahovať premenné rôznych typov - spojité veličiny, disrétne hodnoty, ale aj faktory, alebo štandardný popis vo forme textu.

stlpec1 <- seq(1:5)
stlpec2 <- c(rep("male", 3), rep("female", 2))
stlpec3 <- c("a", "b", "c", "d", "e")

### nie je mozne vytvorit maticu
cbind(stlpec1, stlpec2, stlpec3)
### ale je mozne vytvorit data.frame
data.frame(stlpec1, stlpec2, stlpec3)

S ľubovolným datovým súborom sa následne pracuje ako so štandardnou maticou - môžeme sa odkazovať ne jednotlivé elementy, celé stĺpce (väčšinou jednotlivé premenné), alebo konkrétne řádky (väčšinou jednotlivé pozorovania). Nasledujúce príkazy dávajú stejné výstupy, spôsob implementácie príkazu je ale zakaždým iný. ´

head(Orange) ### zaujima nas pouze promena circumference (treti sloupec)

attach(Orange)
circumference

Orange$circumference

Orange[,3]

K načítaniu vlastných datových súborov slúži niekoľko príkazov. V programe R je technicky možné načítať vpodstate akýkoľvek datovy súbor (txt, csv, xls, ale aj jpg, png, mp3 a mnoho ďalších). Najčastejšie používané príkazy sú read.table() a read.csv(), ktoré vyžaduju súbor uložený vo formáte txt resp. csv. Dodatočné parametre, ktoré upresňujú format podkladového spboru súboru a spôsob jeho načítania do R, je možné využiť pri volaní každého príkazu. Viac podrobnosti v príslušnom helpe ?read.table prípadne ?read.csv.

Analogicky je možné stejne súbory pomocou programu R uložiť a zapísať konkrétny subor na disk.

Užitočné

Program R umožňuje priamo načítať aj online zdroje. Potrebné je iba špecifikovať link a použiť vhodny príkaz. Napríklad:
```
test.data <- read.table("http://www.karlin.mff.cuni.cz/~maciak/NMSA407/NMSA407-1617-HW1.txt", header=T)
```

Pomocou príkazov dim(), str(), head(), alebo summary() preštudujte štruktúru datového súboru a pomocou vhodných indexov, názvu premenných, alebo pomocou príkazov subset() a select() sa pokúste definovať rôzne podsúbory (v zmysle pozorovaní – riadkov, aj v zmysle premenných – stĺpcov). Napríklad:

### pouze vzdalenosti do 50 km
data50_1 <- test.data[test.data$travelDistance == "< 50", ] 
data50_2 <- test.data[test.data[,2] == "< 50", ] 
data50_3 <- subset(test.data, travelDistance == "< 50")

### jedna sa o tri ekvivalentne datove subory
all.equal(data50_1, data50_2)

## [1] TRUE

all.equal(data50_1, data50_3)

## [1] TRUE

### nebo pouze data, kde stopovali vzdy dvaja a trip byl v CR
dataCR2_1 <- test.data[test.data$hitchhikers == "two guys" & test.data$country == "CZ", ]
dataCR2_2 <- test.data[test.data[,1] == "two guys" & test.data[,4] == "CZ", ]
dataCR2_3 <- subset(test.data, hitchhikers == "two guys" & country == "CZ")

### opat sa jedna o tri ekvivalentne datove podsubory
all.equal(dataCR2_1, dataCR2_2)

## [1] TRUE

all.equal(dataCR2_1, dataCR2_3)

## [1] TRUE

### vyberieme data, kde krajina je CZ nebo "other", stopuju pouze dva a zaujima nas  waiting time a vzdalenost
dataCZother_1 <- test.data[(test.data$country == "CZ" | test.data$country == "other") & test.data$hitchhikers == "two guys", c(2,5)]
dataCZother_2 <- test.data[(test.data$country == "CZ" | test.data$country == "other") & test.data$hitchhikers == "two guys", -c(1,3,4,6,7)]
dataCZother_3 <- test.data[(test.data$country == "CZ" | test.data$country == "other") & test.data$hitchhikers == "two guys", names(test.data) %in% c("travelDistance", "waitingTime")]
dataCZother_4 <- test.data[(test.data[,4] == "CZ" | test.data[,4] == "other") & test.data[,1] == "two guys", names(test.data) %in% c("travelDistance", "waitingTime")]
dataCZother_5 <- subset(test.data, (country == "CZ" | country == "other") & hitchhikers == "two guys", select = c("travelDistance", "waitingTime"))

### opat sa jedna o ekvivalentne datove subory
all.equal(dataCZother_1, dataCZother_2)

## [1] TRUE

all.equal(dataCZother_1, dataCZother_3)

## [1] TRUE

all.equal(dataCZother_1, dataCZother_4)

## [1] TRUE

all.equal(dataCZother_1, dataCZother_5)

## [1] TRUE

Použijte Google/ChatGPT a Rkový help ku k príslušným knižniciam a pokúste sa načítať do programu R napr. súbor typu xls/xlsx, alebo obrázok (jpg súbor). Podívajte sa, čo ste do Rka načítali (resp. aký objekt je v Rku po načítaní uložený).
Datové subory vo formáte xls/xlsx sú pomerne časté a v praxi je program Microsoft Excel hodne obľúbený a široko rozšírený. V programe R existuje niekoľko rôznych knižníc a funkcii určených k načítaniu datových súborov z Excelu (napr. knižnica readxl, knižnica xlsx, alebo knižnica XLConnect).

5. Príprava a programovanie vlastných R skriptov

Štatistický program R je možné plnohodnotne využiť aj ako štandardný programovací/developerský jazyk, Ppomocou Rka je možné naprogramovať takmer akúkoľvek úlohu. Je pri tom užitočné poznať rôzne nástroje a konkrétne funkcie, ktoré k tomuto účelu môžu poslúžiť (najdôležitejšie sú asi implementácia cyklov for a while, alebo overovanie podmienok pomocou if).

for cyklus:
```
for (i in 1:10){
  print(1:i)
}
```

while cyklus:

i <- 1
while (i < 10){
  print(1:i)
  i <- i + 1
}

overovanie podmienky pomocou if

N <- rpois(1, lambda = 10)
if (N %% 2 == 0){
  print("Sude cislo")
} else {
  print("Liche cislo")
}

Pri používaní cyklov je dôležité dôsledne skontrolovať skript - môže sa totiž veľmi jednoducho stáť, že program R neúmyselne zacyklíme. V takom prípade je dobré poznať klávesovú skratku CTR + C (v konzole programu R), prípadne ESC v interface RStudio.

Užitočné

Jednotlivé príkazy je možné v programe R vzájomne kombinovať a vytvárať z nich komplexné Rkové skripty. K dispozícii je samozrejme mnoho ďalších nástrojov určených k pokročilému programovaniu (viď napr. help a rôzne tutoriály).

II. Odhady různych charakteristik rozdělení pravděpodobnosti

Primárne je program R určený pro statistickou analýzu dat. Každá analýza libovolného datového souboru by měla pozostávat z dvou časti: popisná (resp. exploratívna) 0část a analyticka (resp. konfirmačná) čast. V popisnej časti se pomoci popisnych charakteristik (průměr, výberový rozptyl, vyběrový median, a pod.) a obrázků/grafů snážíme získat dostatečne kvalifikovaný pohled na data – t.j., jakou maji struktúru, co s čím vzájemně souvisí a pod. Na základne popisné části pak volíme vhodný pravděpodobnostní/statistický model, pomocou kterého data formálne analyzujeme (napr. dvouvýběrový t-test, analýza rozptylu, nebo model lineární egrese, a pod.).

1. Výpočet základních popisných charakteristik náhodného výběrů

Budeme uvažovať náhodný výběr z normálního rozdělení, který jsme vygenerovali v časti 3 a který máme v programe R uložený jako objekt x.

## nagenerovane hodnoty nehodneho vyberu
print(x)
## výběrový průměr
mean(x)
## výběrový medián
median(x)
## výběrový rozptyl
var(x)

## minimum
min(x)
## maximum
max(x)

## summary
summary(x)

## usporadany nahodny vyber
sort(x)
## 22. pořádková statistika
sort(x)[22]
## vektor pořadí
rank(x)

Užitočné

Pomocou helpu a niekoľkých príkladov sa pozreite ma rozdiel vo fungovaní následujúcich funkcií: rank(), sort() a order().

2. Tabulace náhodného výběru, empirické relativní četnosti

## najdi, v kterých intervalech pozorování jsou
tbl=cut(x,breaks=c(-Inf,-3:3,Inf))

## kolik pozorování je v kterém intervalu
table(tbl)
## empirické relativní četnosti
table(tbl)/length(x)

3. Grafické zobrazování náhodného výběru a výběrových charakteristik

Boxplot (krabicový diagram) zobrazuje krabici ohraničenou horním a dolním výběrovým kvartilem (v krabici sa teda nachádza 50 % pozorování), výběrový medián tlustou vodorovnou čárou uvnitř krabice, vousy, které jdou k najvzdialenejšiemu pozorovaniu v rámci jedenapůlnásobku délky boxu (mezikvartilového rozmezí) nahoru i dolů. Pozorování, která se nevejdou do vousů, jsou nakreslena zvlášť (odlehlá pozorování/outliers).

boxplot(x)

Stejný boxplot, ale s väčším dôrazom na estetickú formu (prezentovaný grafický výstup by nikdy neměl byt vytvořen pouze defaultním volaním nějaké funkce). Dobrý obrázek musí obsahovat všechny podstatné informace: popisky na na oboch osích (ve stejném jazyku jako je text celého dokumentu) a případně informaci o použitých jednotkách.

boxplot(x, col = "lightblue", xlab = "Popis pro os x", ylab = "Popis pro os y", main = "Boxplot pro X")

Histogram ukazuje absolutní nebo relativní četnosti pozorování v určitých intervalech (lze měnit jejich počet i umístění). Histogram s relativními četnosti odhaduje hustotu pomocí po částech konstantní funkce. Histogram ale není konsistentním odhadem hustoty.

## histogram s absolutními četnostmi
hist(x, col = "lightblue")

Stanovíme si rozsah osy x, namaluje e reativní četnosti místo absolutních a proložíme teoretickou hustotu rozdělení N(0,1).

## Stanovíme rozsah osy x na (-3.5,3.5), 
## malujeme relativní četnosti místo absolutních
hist(x,freq=FALSE,xlim=c(-3.5,3.5), col = "lightblue")

## Spočítej a domaluj hustotu N(0,1) na intervalu (-3.5,3.5)
xpts=seq(-3.5,3.5,length=500)
lines(xpts,dnorm(xpts,mean=0,sd=1), col = "red", lwd = 2)

Nyní zvětšíme rozsah výběru na n=5000.

x=rnorm(n=5000,mean=0,sd=1)

hist(x,freq=FALSE,xlim=c(-3.5,3.5), ylim = c(0, 0.6), col = "lightblue")
xpts=seq(-3.5,3.5,length=500)
lines(xpts,dnorm(xpts,mean=0,sd=1), col = "red")

Podobá se tvar histogramu hustotě rozdělení, z kterého data pocházejí?

Ešte raz zväčšíme rozsah vyberu na n=50000 a zjemníme historam.

x=rnorm(n=50000,mean=0,sd=1)

hist(x,freq=FALSE,xlim=c(-3.5,3.5), ylim = c(0, 0.6), breaks = 40, col = "lightblue")
xpts=seq(-3.5,3.5,length=500)
lines(xpts,dnorm(xpts,mean=0,sd=1), col = "red")

Empirická distribuční funkce je nestranným a konsistentním odhadem skutečné distribuční funkce (viz větník a poznámky z přednášek).

\[ \widehat{F}_n(u)=\frac{1}{n}\sum_{i=1}^n\mathbb{I}_{(-\infty,u\rangle}(X_i) \]

Ukažme si empirickou distribuční funkci pro výběr z N(0,1) o rozsahu n=30.

## Rozsah výběru 30
x=rnorm(n=30,mean=0,sd=1)

## Výpočet a malování empirické distribuční funkce
plot(ecdf(x))

## Stanovíme rozsah osy x na (-3.5,3.5), zmenšíme puntíky, přepíšeme titulek
plot(ecdf(x),xlim=c(-3.5,3.5),cex.points=0.5, 
     main="Empiricka distribucni funkce N(0,1), n=30")

## Spočítej a namaluj skutečnou d.f. N(0,1) na intervalu (-3.5,3.5)
xpts=seq(-3.5,3.5,length=500)
lines(xpts,pnorm(xpts,mean=0,sd=1), col = "red")

Nyní zvětšíme rozsah výběru na n=200.

## Rozsah výběru 200
x=rnorm(n=200,mean=0,sd=1)

## Empirická distribuční funkce
plot(ecdf(x),xlim=c(-3.5,3.5),cex.points=0.1, 
     main="Empiricka distribucni funkce N(0,1), n=200")

## Spočítej a namaluj d.f. N(0,1) na intervalu (-3.5,3.5)
xpts=seq(-3.5,3.5,length=500)
lines(xpts,pnorm(xpts,mean=0,sd=1), col = "red")

Podobá se empirická distribuční funkce distribuční funkci, z které data pocházejí?

4. Nestrannost a konsistence bodových odhadů (pomoci simulace)

Pre nagenerovaný náhodný výber (uložený v objekte x) si spočítame jednu konkrétnu výberovú charakteristiku (napr. výberový rozptyl) a pomocou simulácii a obrázkov si ilustrujme pojmy nestrannosť a konsistence. Pre výberový rozptyl platí nasledujúci vzťah:

\[ S^2_n=\frac{1}{n-1}\sum_{i=1}^n (X_i-\overline{X}_{n})^2 \]

a) Nestrannost výběrového rozptylu

## Spočítej 1000x výběrový rozptyl výběru o rozsahu 100
nsim=1000
nobs=100
xmat=rnorm(n=nsim*nobs,mean=0,sd=1)
xmat=matrix(xmat,nrow=nobs,ncol=nsim)

dim(xmat)

## [1]  100 1000

## V i-tém řádku je 1000 opakování pozorování X_i.
## V j-tém sloupci je náhodný výběr o rozsahu 100.
## Z každého výběru spočítáme výběrový rozptyl, bude jich 1000.
rozptyly=apply(xmat,2,var) 

## Podívejte se na prvních 50 odhadů rozptylu.
## Podobají se skutečnému rozptylu?
rozptyly[1:50]

##  [1] 1.0676066 0.9815202 1.0174550 1.0086134 1.0603656 0.9888448 0.9298011
##  [8] 0.8698954 0.9076893 0.9081903 0.8595972 0.8219792 0.9054385 1.3799720
## [15] 0.9000663 1.0610294 1.0370524 0.9762047 1.0340270 1.0140795 0.9056098
## [22] 1.0929593 1.0718781 0.8551343 1.1811898 0.8473784 0.7400396 1.2659727
## [29] 1.0051653 0.7469222 1.1531585 1.0664515 0.9683733 0.9200197 1.0125006
## [36] 1.1579836 1.1982492 1.1159301 0.8902754 1.4508194 1.3865924 0.8998533
## [43] 1.0361340 1.0644491 1.1896048 0.9891542 1.0022059 1.0414074 0.9181797
## [50] 0.7876316

## Jaký je průměr z 1000 odhadů rozptylu?
mean(rozptyly)

## [1] 1.006246

Odpovídá tento výsledek tomu, co bychom očekávali od nestranného odhadu? Jak by ste prezentovali tento výsledek graficky?

b) Konsistence výběrového rozptylu

## Spočítej výběrový rozptyl z výběru o rozsahu 10, 100, 200, 1000
rozsahy=c(10,100,1000,10000, 100000)

## zopakuj nezavisle 5 krat
set.seed(123)
opakovat <- 5

## Sem dáme výsledky
vybrozp <- matrix(rep(0, opakovat * length(rozsahy)), nrow = opakovat)
colnames(vybrozp)=paste("n",rozsahy,sep="=")

for (o in 1:opakovat){
  x <- rnorm(max(rozsahy),mean=0,sd=1)
  for(n in rozsahy)
  {
    vybrozp[o, paste("n",n,sep="=")]=var(x[1:n])
  }
}

Tohle lze zobrazit graficky napr. pomocou nasledujúceho obrázku.

plot(0,0, xlim = c(1,length(rozsahy)), ylim = c(0.5,1.6), pch = "", xlab = "Logaritmus rozsahu výběru", ylab = "Vyběrový rozptyl")
for (o in 1:opakovat){
  lines(vybrozp[o,] ~ seq(1,length(rozsahy)), col = "red")
}
abline(1,0, col = "black", lty = 2)

Užitočné

Pomocou helpu naštudujte, k čomu slúžia funcie apply(), tapply(), aggregate() a by();
Použijte tieto funkcie v súvislosti s datovým súborom ‘test.data’ a pomocou nich spočítajte niektoré empirické charakteristiky;
V pripade, že pracujeme v programe R len s jednym datovým súborom, je možné prácu (resp. implementáciu príkazov) mierne zjednodušiť - pomocou príkazu attach(test.data) program R zprístupni jednotlivé premenné, ktoré je potom možne volať priamo pomocou ich názvu:
```
attach(test.data)
mean(waitingTime[driverGender == "male"])
mean(waitingTime[driverGender == "female"])
```
V prípade, že takto zprístupnene premenné už neplánujeme využívať, je dobre využiť komplementárny príkaz detach().
Je dôležité správne pochopiť fungovanie Rka a odkazovanie na premenné po použití príkazu attach(). Program R totiž zprístupni premenne v stave, v akom sa v momente volania príkazu attach() aktuálne nachádzaju v data.frame. Akýkoľvek následný zásah do data.framu a zmena hodnôt v ňom, nebude reflektovaná ak sa na premennu odkážeme pomocou jej názvu a nie priamo cez odkaz na data.frame;

Samostatná úloha

Zopakujte Časť 4 (t.j. nestrannosť a konzistenci bodových odhadov) pre náhodný výber generovaný z exponenciálneho rozdelenia s parametrom \(\lambda = 2\) a strednou hodnotiou \(1/2\) a z Cauchyho rozdelenia C(0,1). Zamerajte sa na nestrannosť a konzistenci výberového rozptylu (alebo výberovej strednej hodnoty, čo si vyžaduje určité modifikácie vyššie uvedeného R-kového kódu)-

## Generování veličin z Exp(2): 
rexp(n,2)
## Hustota Exp(2): 
dexp(x,2)
## Distribuční funkce Exp(2): 
pexp(x,2)

## Generování veličin z C(0,1): 
rcauchy(n,0,1)
## Hustota C(0,1): 
dcauchy(x,0,1)
## Distribuční funkce C(0,1): 
pcauchy(x,0,1)

NMFM301, ZS 2024/25

Cvičenie 7 (praktická časť 1)

(Základy práce so štatistickým programom R)

I. Štatistický program R (stručný úvod)

1. Základné matematické operácie v R

Užitočné

2. Práca s vektormi a maticami

Užitočné

3. Náhodné data – generátory (pseudo) náhodných čísel

Užitočné

4. Reálne data a vlastné datové súbory

Užitočné

5. Príprava a programovanie vlastných R skriptov

Užitočné

II. Odhady různych charakteristik rozdělení pravděpodobnosti

1. Výpočet základních popisných charakteristik náhodného výběrů

Užitočné

2. Tabulace náhodného výběru, empirické relativní četnosti

3. Grafické zobrazování náhodného výběru a výběrových charakteristik

4. Nestrannost a konsistence bodových odhadů (pomoci simulace)

a) Nestrannost výběrového rozptylu

b) Konsistence výběrového rozptylu

Užitočné

Samostatná úloha