<table width = 900 border = 0><tr><td>
# <b>NMST539 | Cvičenie 4</b>

## <b>Wishartovo, Hotellingovo a Wilkovo rozdelenie</b>
### <font color = "#880000"><b>LS 2020/2021 |  23/03/21 | (online výuka)</b></font>
###### <a href = "cvicenie2021_4.Rmd">zdrojový Rmd súbor (kódovanie UTF8)</a>

**************

<b>Outline štvrtého cvičenia:</b> 

<ul>
<li>náhodné matice a niektoré mnohorozmerné rozdelenia -- Wishartove, Hottellingovo a Wilkovo Lambda rozdelenie;</li>
<li>mnohorozmerné testy pre vektor strednej hodnoty (mnohorozmerné  analógie jednovýberových a dvojvýberových t-testov);</li>
</ul><br>


<b>Štatisticky program R</b> je k dispozícii (GNU public licence) na adrese <a href = "https://www.r-project.org">https://www.r-project.org</a><br><br>

<b>RStudio</b> -- "user-friendly" interface (jeden z mnohých, ktoré na internete nájdete): <a href = "https://www.rstudio.com">RStudio</a>.<br><br>

<b>Užitočné návody a manuály pre prácu s programom R:</b><br>

<ul>
<li>Bína, V., Komárek, A. a Komárková, L.: <i>Jak na jazyk R.</i>  (<a href = "NMFM301/Rmanual2.pdf">PDF súbor</a>)</li>
<li>Komárek, A.: <i>Základy práce s R.</i>  (<a href = "http://msekce.karlin.mff.cuni.cz/~komarek/vyuka/dataRko/Rmanual1.pdf">PDF súbor</a>)</li>
<li>Kulich, M.: <i>Velmi stručný úvod do R.</i> (<a href = "NMFM301/uvodr.pdf" >PDF súbor</a>)</li>
<li>De Vries, A. a Meys, J.: <i>R for Dummies.</i> (ISBN-13: 978-1119055808) </li>		
</ul><br>

<b>Odborná literatúra:</b><br>
<ul>
<li>Hardle, W. and Simar.L.: Applied Multivariate Statistical Analysis. Springer, 2015</li>
<li>Mardia, K., Kent, J., and Bibby, J.:Multivariate Analysis, Academic Press, 1979.</li>
</ul>

<br><br>

**************


#### <font color = "#880000"><b>1. Wishartovo pravdepodobnostné rozdelenie</b></font>

Wishartovo rozdelenie je mnohorozmerným zobecnením jednorozmerného  $\chi^2$ rozdelenia, ktoré sa štandardne používa 
pre inferenciu o neznámom parametri rozptylu na základe jednorozmerného náhodného výberu. Rozdelenie nesie názov 
po svojom autorovi,  <a href = "https://en.wikipedia.org/wiki/John_Wishart_(statistician)">John-ovi Wishar</a>t-ovi, ktorý toto rozdelenie formuloval v roku 1928. Wishartovo rozdelenie je základným nástrojom 
pre analýzu/inferenciu ohľadom variančnej-kovariančnej matice na základe mnohorozmerného náhodného výberu  $\boldsymbol{X}_{1}, \dots, \boldsymbol{X}_{n}$, pre  $p$ rozmerné náhodné   vektory $\boldsymbol{X}_{i} = (X_{i 1}, \dots, X_{i p})^\top \in \mathbb{R}^p$ a rozsah náhodného výberu $n \in \mathbb{N}$. 

Wishartovo rozdelenie je mnohorozmerným zobecnením jednorozmerného  $\chi^2$ rozdelenia v nasledujúcom zmysle: pre mnohorozmerný náhodný výber 
 $\boldsymbol{X}_{1}, \dots, \boldsymbol{X}_{n}$ (stĺpcové náhodné vektory) z $p$ rozmerného normálneho rozdelenia $N_{p}(\boldsymbol{0}, \Sigma)$ s nulovým vektorom stredných hodnôt a variančnou-kovariančnou maticou $\Sigma$ (symetrická a pozitívne definitná matica) má príslušná kvadratická forma Wishartovo rozdelenie
 
$$
\mathbb{X}^{\top}\mathbb{X} \sim W_{p}(\Sigma, n), 
$$
s parametrami $p \in \mathbb{N}$ (dimenzia), $n \in \mathbb{N}$ (prozsah náhodného výberu) a $\Sigma$ (variančná-kovariačná matica náhodného výberu $\boldsymbol{X}_{1}, \dots, \boldsymbol{X}_{n}$ kde $\mathbb{X} = (\boldsymbol{X}_{1}, \dots, \boldsymbol{X}_{n})^{\top}$). Podobne ako v jednorozmernom prípade pre náhodný výber $X_{1}, \dots, X_{n}$ (náhodné veličiny) z normálelneho rozdelenia $N(0,1)$, má príslušná kvadratická forma 
$$
\boldsymbol{X}^{\top}\boldsymbol{X} \sim \chi_{n}^2, 
$$
$\chi^2$ rozdelenie, pre $\boldsymbol{X} = (X_{1}, \dots, X_{n})^\top$. Je zrejmé, že pre náhodný výber s rozsahom $n \in \mathbb{N}$ z jednorozmerného normálneho rozdelenia $N(0,1)$ má náhodná veličina $\mathbb{X}^\top\mathbb{X}$  (jednorozmerné) Wishartovo rozdelenie $W_{1}(1, n)$, ktoré je ekvivalentné s $\chi_{n}^2$ rozdelením. 

Wishartovo rozdelenie predstavuje rodinu náhodných rozdelení pre definovaných na symetrických, pozitívne semi-definitných .<i>náhodných maticiach</i>.  Príslušná hustota Wishartovho rozdelenia má nasledujäci tvar: 

$$
f(\mathcal{X}) = \frac{1}{2^{np/2} |\Sigma|^{n/2} \Gamma_{p}(\frac{n}{2})} \cdot |\mathcal{X}|^{\frac{n - p - 1}{2}} e^{-(1/2) tr(\Sigma^{-1}\mathcal{X})}, 
$$
kde $\mathcal{X}$ je náhodná matica typu $p\times p$ a  $\Gamma_{p}(\cdot)$ predstavuje mnohorozmerné zobecnenie jednorozmernej Gamma funkcie $\Gamma(\cdot)$. V programe R sú k dispozícii rôzne možnosti (knižnice a príkazy) pre prácu s Wishartovým rozdelením. Niektoré z týchto nástrojov 
využijeme v následujúcej časti. <br><br>

<br><br>

#### <font color = "#880000"><b>Samostatne (teoretické a praktické úlohy)</b></font>
<hr size="3" width="100%" noshade style="color:#440000" align="left" />
<ul>
<li>Uvažujte mnohorozmerný náhodný výber $\boldsymbol{X}_{1}, \dots, \boldsymbol{X}_{n}$ z $p \in \mathbb{N}$ mnohorozmerného normálneho rozdelenia  $N_{p}(\boldsymbol{\mu}, \Sigma)$, s vektorom stredných hodnôt $\boldsymbol{\mu} \in \mathbb{R}^p$ a variančnou-kovariančnou (pozitívne definitnou) maticou $\Sigma$. Ukážte, že náhodná matica definovaná predpisom $n\mathbb{X}^\top \mathcal{H}\mathbb{X}$ má  Wishartovo rozdelenie  $W_{p}(\Sigma, n - 1)$, kde  $\mathcal{H}$ idenpotentná, tzv. centrovacia matica $\mathcal{H} = \mathbb{I}_{n} - \frac{1}{n}\boldsymbol{1}_{n}\boldsymbol{1}_{n}^\top$.</li><br>

<li>V štandardnej inštalácii programu R je pre Wishartovo rozdelenie určená funkcia `rWishart()` -- t.j. generátor náhodných matíc 
s Wishartovým rozdelením. Pomocou helpu v programe R a uvedených príkladov naštudujte implementáciu daného príkazu.</li><br>

<li>Ďalšie užitočné príkazy a funkcie pre prácu s Wishartovým rozdelením sú napr. funkcia `Wishart()` v knižnici 'MCMCpack', funkcia `Wishart()` v knižnici 'mixAK', funkcia `rwishart()` z knižnice 'dlm', alebo `dist.Wishart()` z knižnice 'LaplacesDemon'. </li>
</ul>
<hr size="3" width="100%" noshade style="color:#440000" align="left" />
<br>

Pre jednoduchú ilustráciu použijeme jednorozmerný generátor: náhodný výber z Wishartovho rozdelenia   $W_{1}(\Sigma = 1, n = 10)$ (resp. ekvivalentne z  $\chi^2$ rozdelenia s $n = 10$ stupňami voľnosti) získame pomocou pomocu funkcie `rWishart()` ako

```{r}
S <- as.matrix(1) ### variancna-kovariancna matica (i.e., rozptyl)
sample <- rWishart(10, df = 10, S)
```

Pre obecný rozmer $p \in \mathbb{N}$ použijeme analogický príkaz s vhodne definovanou variančnou-kovariačnou maticou $\Sigma$:

```{r}
Sigma <- cbind(c(1,0,0,0), c(0,1,0,0), c(0,0,1,0), c(0,0,0,1)) ### matica typu 4x4
rWishart(2, df = 100, Sigma) ### sample of size two
```


V prípade, že $p = 1$, môžeme ekvivalentnosť Wishartovho rozdelenia $W_{1}(\Sigma = 1, n)$ a $\chi^2$ rozdelenia s $n$ stupňami voľnosti
jednoducho  overiť aj vizuálne, napr. pomocou histogramu, alebo príslušného neparametrického odhadu hustôt. 

```{r, fig.width = 10}
set.seed(1234)
n <- 10 ### pocet stupnov volnosti
sampleWishart <- rWishart(5000, df = n, S)
sampleChiSq <- rchisq(5000, df = n) 

par(mfrow = c(1,2))
hist(sampleWishart, col = "lightblue", main = expression(paste("Wishart Distribution", sep ="")), xlim = c(0, 40), breaks = 30, freq = F)
lines(density(sampleWishart), col = "red", lwd = 2)
hist(sampleChiSq, col = "lightblue", main = expression(paste(chi^2, "Distribution", sep = "")), xlim = c(0,40), breaks = 30, freq = F)
lines(density(sampleChiSq), col = "red", lwd = 2)
```

<br><br>

Dvojrozmerné Wishartovo rozdelenie (rozdelenie náhodných symetrických matíc typu $2 \times 2$) je už ale výrazne náročnejšie vizualizovať vhodným 
spôsobom. Jedná z možnosti je napr. následujúci graf, ktorý zobrazuje $2 \times 2$ maticu ako usporiadaný vektor - t.j., náhodná matica $\boldsymbol{X}_i$ je tvorená dvoma vektormi $(x_{11}, x_{21})^\top$ a $(x_{12}, x_{22})^\top$. Tieto vektory definujú v dvojrozmernej $xy$ rovine ($xy$ scatterplot) počiatočný a koncový bod úsečky. Úsečky su následne vyzualizované v grafe (počiatočný bod -- t.j. bod určený vektorom $(x_{11}, x_{21})^\top$ je v grafe zvýraznený). <br><br>

Graf možno nie je úplne intuitívny, ale umožňuje napríklad  jednoznačnu a priamočiaru reprodukciu pôvodných hodnôt náhodného výberu. 

```{r, fig.width = 10}
sample2DWishart <- rWishart(100, df = 10, cbind(c(1,1), c(1,3)))

plot(0,0, pch = "", xlim = c(min(sample2DWishart[1,,]), max(sample2DWishart[1,,])), ylim = c(min(sample2DWishart[2,,]), max(sample2DWishart[2,,])),
     xlab = "First column vector", ylab = "Second column vector")
for(i in 1:100){
  lines(sample2DWishart[1,,i], sample2DWishart[2,,i], col = "red")
  points(sample2DWishart[1,1,i], sample2DWishart[2,1,i], pch = 21, bg = "red", cex = 0.8)
}
```

Iné možnosti grafickej reprezentácie náhodného výberu z Wishartovho rozdelenia by boli napríklad vizualizácie založené na vlastných 
čislach, vlastných vektoroch, alebo ďalšie. Väčšinou ale neumožňujú takto jednoduchú spätnú rekonštrukciu náhodných matíc na základe samotného grafu. 

<br><br>



#### <font color = "#880000"><b>Samostatne</b></font>
<hr size="3" width="100%" noshade style="color:#440000" align="left" />
<ul>
<li>V predcházajúcom príklade použijte odlišnú štruktúru pre variančnú-kovariančnú maticu $\Sigma$ alebo voľbu počtu stupňov 
voľnosti a pomocou obrázku sledujte vplyv týchto parametrov na celkovú štruktúru náhodného výberu reprezentovanú pomocou $xy$ scatterplotu. 
Zamerajte sa hlavne na vplyv jednotlivých parametrov variančnej-kovariančnej matice $\Sigma$. </li><br>

<li>Uvažujte obecnú maticu $A \in \mathbb{R}^{p \times q}$ pre nejaké vhodné $q \in \mathbb{N}$, tak, aby $p \neq q$.
Opäť pomocou grafickej vizualizácie oveřte (napríklad pomocou histogramu, alebo KDE), že $A^\top \mathbb{M} A \sim W_{q}(A^\top\Sigma A, n)$, kde $\mathbb{M} \sim W_{p}(\Sigma, n)$; </li>
</ul>
<hr size="3" width="100%" noshade style="color:#440000" align="left" />
<br><br>

#### <font color = "#880000"><b>2. Hotellingovo $\boldsymbol{T^2}$ rozdelnie</b></font>

Podobným spôsobom, ako je v jednorozmernom prípade definované $t$-rozdelenie s $n$ stupňami voľnosti, môžeme defnovať aj jeho mnohorozmerné 
zobecnenie. V jednorozmernom prípade sa jedná o náhodnú veličinu s normálnym $N(0,1)$ rozdelením, ktorá je podelená inou, na nej nezávislou náhodnou veličinou, ktorej kvadrát má $\chi^2$ s $n$ stupňami voľnosti (a navyše je štandardizovaná odmocninou z týchto stupńov voľnosti). 

Mnohorozmerné zobecnenie  (t.j., náhodná veličina s Hotellingovým $T^{2}$ rozdelením) je definované predpisom

$$
n \boldsymbol{Y}^{\top} \mathbb{M}^{-1} \boldsymbol{Y} \sim T^{2}(n, p),
$$
kde $p \in \mathbb{N}$ je dimenzia/rozmer náhodného vektoru $Y \sim N_{p}(0, \mathbb{I})$ a $n \in \mathbb{N}$ je parameter Wishartovho rozdelenia náhodnej matice $\mathbb{M} \sim W_{p}(\mathbb{I}, n)$. Analogicky ako v jednorozmernom prípade, prepokládame, 
že náhodný vektor $\boldsymbol{Y}$  je nezávsilý od náhodnej matice $\mathbb{M}$.


V špeciálnom prípade, pre $p = 1$, dostaneme jednorozmerné Fisherovo F rozdelenie  s jedným a s $n$ stupňami voľnosti  (čo je vlastne ekvivalentné s rozdelením kvadrátu náhodnej veličiny s $t$ rozdelením s   $n$ stupňami voľnosti).  Hotellingovo $T^2$ rozdelenie s prametrami  $p, n \in \mathbb{N}$ preto možno považovať aj za mnohorozmerné zobecnenie Fisherovho F rozdelenia. Medzi oboma rozdeleniami je dokonca jednoznačná analytická súvislosť,
ktorú je možné vyjadriť pomocou vzťahu

$$
T^{2}(p, n) \equiv \frac{n p}{n - p + 1}F_{p, n - p + 1}.
$$

Náhodný výber z jednorozmerného Fisherovho rozdelenia preto môže byť efektívne použitý aj k tomu, aby sme získali kritické hodnoty 
mnohorozmerného  Hotellingovho $T^2$ rozdelenia. Príslušná transformácia medzi Hotellingovým $T^2$ rozdelením a Fisherovým F rozdelením 
závisí iba na parametroch $n, p \in \mathbb{N}$. <br><br>

Úlohu týchto parametrov (resp. efekt týchto parametrov) v Hotellingovom $T^2$ rozdelení môžeme opäť jednoducho vizualizovať, napríklad 
pomocou náhodného generátora z Fisherovho F rozdelenia (v programe R príkaz `rf()`). 

```{r, fig.width = 10, fig.height=10}
samples <- NULL
samples <- cbind(samples, rf(n = 1000, df1 = 1, df2 = 1))
samples <- cbind(samples, rf(n = 1000, df1 = 1, df2 = 10))

samples <- cbind(samples, rf(n = 1000, df1 = 10, df2 = 10))
samples <- cbind(samples, rf(n = 1000, df1 = 10, df2 = 100))

samples <- cbind(samples, rf(n = 1000, df1 = 100, df2 = 10))
samples <- cbind(samples, rf(n = 1000, df1 = 100, df2 = 100))

samples <- cbind(samples, rf(n = 1000, df1 = 1000, df2 = 100))
samples <- cbind(samples, rf(n = 1000, df1 = 1000, df2 = 1000))

par(mfrow = c(4,2))
for (i in 1:dim(samples)[2]){
  hist(samples[,i], xlab=paste("Sample no.", i, sep = ""), col = "yellow", main = "", freq = F, breaks = 20)
  lines(density(samples[,i]), col = "red")
}
```

Pre pripomenutie, náhodná veličina $X$ s Fisherovým F rozdelením $F_{df_1, df_2}$ má strednú hodnotu a rozptyl 
definované predpismi: 

<ul>
<li> $E X = \frac{df_{2}}{df_{2} - 2}$, pre $df_2 > 2$;</li>
<li>  $Var(X) = \frac{2 df_{2}^2 (df_1 + df_2 - 2)}{df_1 (df_{2} - 2)^2 (df_2 - 4)}$ pre $df_2> 4$;</li>
</ul>

<br><br>
Obecne platí, že pre mnohorozmerný náhodný výber $\boldsymbol{X}_{1}, \dots, \boldsymbol{X}_{n} \sim N_{p}(\boldsymbol{\mu}, \Sigma)$ s neznámym vektorom stredných hodnôt $\boldsymbol{\mu} \in \mathbb{R}^{p}$ a neznámou variančnou-kovariančnou maticou $\Sigma$, dostaneme pre výberový priemer 
$\overline{\boldsymbol{X}}$ rozdelenie 

$$
(n - 1)\Big(\overline{\boldsymbol{X}} - \boldsymbol{\mu}\Big)^\top \mathcal{S}^{-1} \Big(\overline{\boldsymbol{X}} - \boldsymbol{\mu}\Big) \sim T^2(p, n - 1), 
$$

čo lze ekvivalentne zapísať aj do tvaru

$$
\frac{n - p}{p} \Big(\overline{\boldsymbol{X}} - \boldsymbol{\mu}\Big)^\top \mathcal{S}^{-1} \Big(\overline{\boldsymbol{X}} - \boldsymbol{\mu}\Big) \sim F_{p, n - p}. 
$$


Toto je následne možné aplikovať pri úlohach o inferencii pre neznámy vektor stredných hodnôt -- napr. pri štatistických testoch, alebo pri konštrukcii konfidenčných oblastí pre vektor neznámych stredných hodnôt $\boldsymbol{\mu}$, kde $\boldsymbol{\mu} = (\mu_{1}, \dots, \mu_{p})^{\top}$. 

V praxi sa často využíva namiesto confidenčného regiónu pre $\boldsymbol{\mu}$ (čo môže byť nepraktické hlavne pre vyšší rozmer dimenzie $p$) radšej súbor intervalov spoľahlivosti pre jednotlivé zložky $\boldsymbol{\mu}$, tak aby bola celková pravdepodobnosť pokrytia pod kontrolou -- najčastejšie vyžadujeme pokrytie minimálne  $(1 - \alpha)\times 100~\%$ pre vhodné a dostatočne malé hodnoty $\alpha \in (0,1)$.

Analogicky, v prípade štatistického testu  definujeme dvojicu hypotéz

$$
H_{0}: \boldsymbol{\mu} = \boldsymbol{\mu}_{0} \in \mathbb{R}^{p}
$$

$$
H_{1}: \boldsymbol{\mu} \neq \boldsymbol{\mu}_{0} \in \mathbb{R}^{p}
$$

a využijeme testovú štatistiku

$$
(n - 1)\Big(\overline{\boldsymbol{X}} - \boldsymbol{\mu}_{0}\Big)^\top \mathcal{S}^{-1} \Big(\overline{\boldsymbol{X}} - \boldsymbol{\mu}_{0}\Big), 
$$

ktorá má za platnosti nulovej hypotézy $T^2(p, n - 1)$ rozdelenie.  Ekvivalentne, testová štatistika vyjadrená v tvare

$$
\frac{n - p}{p} \Big(\overline{\boldsymbol{X}} - \boldsymbol{\mu}_{0}\Big)^\top \mathcal{S}^{-1} \Big(\overline{\boldsymbol{X}} - \boldsymbol{\mu}_{0}\Big) 
$$

má za platnosti nulovej hypotézy Fisherovo $F$-rozdelenie s $p$ a $n - p$ stupňami voľnosti.

<br><ve>

V programe R je k dispozícii knižnica `DescTools` -- nutné nainštalovať pomocou príkazu  `install.packages("DescTools")`a následne inicializovať pomocou príkazu `library(DescTools)`. Samotný test nulovej hypotézy vyššie je potom možné uskutočniť  pomocou funkcie `HotellingsT2Test()` (pre konkrétnu implementáciu príkazu použijte help k danej funkcii). 


Analogickým spôsobom je možné skonštruovať aj konfidenčnú oblasť - tzv. konfidenčný elipsoid  $\boldsymbol{\mu} \in \mathbb{R}^p$. 
Keďže platí, že 
$$
\frac{n - p}{p} \Big(\overline{\boldsymbol{X}} - \boldsymbol{\mu}_{0}\Big)^\top \mathcal{S}^{-1} \Big(\overline{\boldsymbol{X}} - \boldsymbol{\mu}_{0}\Big)  \sim F_{p, n - p},
$$
tak nasledujúca množina 
$$
\left\{\boldsymbol{\mu} \in \mathbb{R}^p;~  \Big( \overline{\boldsymbol{X}} - \boldsymbol{\mu}\Big)^\top \mathcal{S}^{-1} \Big( \overline{\boldsymbol{X}} - \boldsymbol{\mu}\Big) \leq \frac{p}{n - p} F_{p, n- p}(1 - \alpha) \right\}
$$
je konfidenčnou oblasťou pre neznámy vektor stredných hodnôt $\boldsymbol{\mu} \in \mathbb{R}^p$ s pravdepodobnosťou pokrytia  $\alpha = 1 - \alpha$. Jedná sa o  ellipsoid v $\mathbb{R}^p$.<br><br>

<b>Ilustračný príklad v 2D: </b><br>
```{r, fig.width=10}
library(mvtnorm)
s=matrix(c(1,-0.5,-0.5,1),2);x=seq(-3,3,by=0.015) ### variance-covariance matrix

contour(x,x,outer(x,x, function(x,y){dmvnorm(cbind(x,y),sigma=s)}))
n <- 20
X <- rmvnorm(n,sigma=s)
m <- apply(X,2,mean)
S <- cov(X)
points(m[1],m[2],pch=8,col="red",cex=2)

S1=solve(S)
contour(x,x,outer(x,x,function(x,y){(n-2)*
                                      apply(t(t(cbind(x,y))-m),1,function(x){t(x)%*%S1%*%x})<
                                      2*qf(0.95,2,n-2)}),col="red",add=TRUE)
bodyx=m[1]+c(-1,1)*sqrt(S[1,1]*2*qf (0.95,2,n-2) /(n-2))
bodyy=m[2]+c(-1,1)*sqrt(S[2,2]*2*qf (0.95,2,n-2) /(n-2))
polygon(x=bodyx[c(1,1,2,2,1)],y=bodyy[c(1,2,2,1,1)],border="blue")
```




#### <font color = "#880000"><b>Poznámka</b></font>
<hr size="3" width="100%" noshade style="color:#440000" align="left" />
<ul>
<li>Konfidenčný elipsoid je výhodnejší z geometrického hľadiska, keďže pri stejnej hladine spoľahlivosti definuje menšiu plochu;</li>
<li>Súbor konfidenčných intervalov pre jednotlivé zložky -- t.j. kartezský súčin intervalov je na druhu stranu výhodnejší z
praktického hľadiska a poskytuje omnoho lepšiu intuíciu pre pochopenie a interpretáciu. Pri stejnej hladine spoľahlivosti ale pokrýva 
výrazne vačšiu plochu. </li>
</ul>
<hr size="3" width="100%" noshade style="color:#440000" align="left" />
<br><br>



#### <font color = "#880000"><b>Príklad I: ozdiel dvoch vektorov stredných hodnôt / stejná variančná-kovariančná matica</b></font>
Hotellingove $T^2$ rozdelenie použijeme najprv pre test o rovnosti dvoch neznámych vektorových parametrov stredných hodnôt. 

K dispozícii máme náhodný výber $\boldsymbol{X}_{1}, \dots, \boldsymbol{X}_{n} \sim N_{p}(\boldsymbol{\mu}_{1}, \Sigma)$
a na ňom nezávislý náhodný výber $\boldsymbol{Y}_{1}, \dots, \boldsymbol{Y}_{M} \sim N_{p}(\boldsymbol{\mu}_{2}, \Sigma)$. Parametre stredných hodnôt -- t.j., p rozmerné pevné, ale neznáme vektory -- sú obecne rôzne $\boldsymbol{\mu}_{1} \neq \boldsymbol{\mu}_{2}$, ale variačna-kovariančná matica je pre oba náhodné výbery stejná (hoci neznáma). Rovnaká je samozrejme aj dimenzia náhodných výberov.

Zaujíma nás test nulovej hypotézy 
$$
H_{0}: \boldsymbol{\mu}_1 = \boldsymbol{\mu}_{2}
$$
oproti obecnej alternatíve, že nulová hypotéza neplatí. Je zrejmé, že pre výberové priemery $\overline{X}_{1}$  a $\overline{X}_{2}$ platí, že

$$
(\overline{X}_{1} - \overline{X}_{2}) \sim N_{p}\Bigg(\boldsymbol{\Delta}, \frac{n + m}{n m} \Sigma\Bigg),
 $$
 
a tiež 

$$
n\mathcal{S}_{1} + m\mathcal{S}_{2} \sim W_{p}(\Sigma, n + m - 2), 
$$
kde $\mathcal{S}_{1}$ a $\mathcal{S}_{2}$ sú empriciké odhady variančnej-kovariančnej matice  $\Sigma$ spočítané samostatne na základe prvého a druhého náhodného výberu.  Kritický obor pre daný test je preto definovaný následovne: 

$$
\frac{nm(n + m - p - 1)}{p(n + m)^2}(\overline{\boldsymbol{X}}_{1} - \overline{\boldsymbol{X}}_{2})^{\top} \mathcal{S}^{-1} (\overline{\boldsymbol{X}}_{1} - \overline{\boldsymbol{X}}_{2}) \geq F_{p, n + m - p - 1}(1 - \alpha).
$$
<br>
<hr size="3" width="100%" noshade style="color:#440000" align="left" /><br><br>

#### <font color = "#880000"><b>Príklad II: Rozdiel dvoch vektorov stredných hodnôt / rôzne variančné-kovariančné matice</b></font>

Analogicky ako v predchádzajúcom príklade, opäť uvažujme náhodný výber  $\boldsymbol{X}_{1}, \dots, \boldsymbol{X}_{n} \sim N_{p}(\boldsymbol{\mu}_{1}, \Sigma_{1})$ a na ňom nezávislý druhý náhodný výber $\boldsymbol{Y}_{1}, \dots, \boldsymbol{Y}_{M} \sim N_{p}(\boldsymbol{\mu}_{2}, \Sigma_{2})$, s obecne rôznymi parametrami stredných hodnôt $\boldsymbol{\mu}_{1} \neq \boldsymbol{\mu}_{2}$,
ale tentokrát aj s potenciálne rôznymi (neznámymi) variančnými-kovariančnými maticami. Opäť nás zaujíma test nulovej hypotézy o rovnosti vektorových parametrov stredných hodnôt, t.j., nulová hypotéza v tvare

$$
H_{0}: \boldsymbol{\mu}_1 = \boldsymbol{\mu}_{2}
$$
oproti obecnej alternatíve, že mnulová hypotéza neplati. Je zrejmé,  že platí
$$
(\overline{X}_{1} - \overline{X}_{2}) \sim N_{p}\Bigg(\boldsymbol{\Delta}, \frac{\Sigma_{1}}{n} + \frac{\Sigma_{2}}{m}\Bigg),
 $$
 
a taktiež
$$
(\overline{X}_{1} - \overline{X}_{2})^\top \Big(\frac{\Sigma_{1}}{n} + \frac{\Sigma_{2}}{m}\Big)^{-1} (\overline{X}_{1} - \overline{X}_{2}) \sim \chi_{p}^{2}.
$$

S využitím vyššie uvedeného môžeme zostrojiť testovú štatistiku a analogickym spôsobom 
ako v predchádzajom prípade získame kritický obor na základe ktorého rozhodneme o zamietnuti, resp. nezamietnuti nulovej hypotézy. 

<br><br>

#### <font color = "#880000"><b>Samostatne</b></font>
<hr size="3" width="100%" noshade style="color:#440000" align="left" />
<ul>
<li>Nainštalujte si v programe R knižnice  'Hotelling' a 'DescTools';</li>
<li>Pre jednovýberové a dvojvýberové Hotellingové testy využijte príkaz `HotellingsT2Test()` z knižnice 'DescTools';</li>
<li>Pre dvojvýberové testy môžete využiť aj príkaz `hotelling.test()` z knižnice 'Hotelling';</li>
<li>Pomocou helpu k príkazom `hotelling.stat()` a `hotelling.test()`zistite konkrétnu implementáciu 
jednotlivých funkcii;
</ul>
<hr size="3" width="100%" noshade style="color:#440000" align="left" />
<br><br>


<b>Ilustračný príklad: </b><br>
```{r, message = F, warning = F}
library("Hotelling")
```

Využijeme datový súbor `container.df` z knižnice 'Hotelling'. Data reprezentujú konkrétne merania koncentrácie rôznych chemických 
prvkov (podrobnejšie v helpe `?container.df`). 

```{r}
data(container.df)
summary(container.df)
```

Príslušné empirické odhady neznámych vektorov stredných hodnôt aj variačných-kovariačných matíc sú: 

```{r}
apply(container.df[,2:10], 2, mean) ### mean vector estimate
cov(container.df[,2:10])            ### variance-covariance estimate
```


Zaujíma nás test nulovej hypotézy, či (neznáme) stredné koncentrácie v prvom kontajneri (`gp == 1`) sú rovnaké, 
ako v druhom kontajneri (`gp == 2`). <br><br>

```{r}
split.data = split(container.df[,-1],container.df$gp)
x <- split.data[[1]]
y <- split.data[[2]]
statistic <- hotelling.stat(x, y)
```

Samotný test: 

```{r}
hotellingTest <- hotelling.test(.~gp, data = container.df)
hotellingTest
```

<br><br>

#### <font color = "#880000"><b>Otázky</b></font>
<hr size="3" width="100%" noshade style="color:#440000" align="left" />
<ul>
<li>Ako by ste na základe prezentovaných čísel interpretovali výsledok testu?</li>
<li>Ako by ste výsledky vizuálne reprezentovali pomocou nejakého obrázku? </li>
<li>Pre jednotlivé zložky aplikujte jednorozmerý dvojvýberový  $t$-test(s). Ako sa líšia závery?</li>
<li>Ako by ste zostrojili simultánny interval spoľahlivosti pre jednotlivé zložky neznámych vektorových 
parametrov stredných hodnôt?</li>
</ul>
<hr size="3" width="100%" noshade style="color:#440000" align="left" />
<br><br>


Príslušné simultánne intervaly spoľahlivosti pre (všetky) lineárne kombinácie jednotlivých zložiek $\boldsymbol{\mu} \in \mathbb{R}^{p}$ 
(lineárne kombinácie definované parametrami $\boldsymbol{a}^\top$ získame pomocou vzťahu

$$
P\Big(\forall \boldsymbol{a} \in \mathbb{R}^{p};~ \boldsymbol{a}^\top \boldsymbol{\mu} \in \big(
\boldsymbol{a}^\top\overline{\boldsymbol{X}} - \sqrt{K_{\alpha} \boldsymbol{a}^\top \mathcal{S} \boldsymbol{a}}, \boldsymbol{a}^\top\overline{\boldsymbol{X}} + \sqrt{K_{\alpha} \boldsymbol{a}^\top \mathcal{S} \boldsymbol{a}} \big)
\big)\Big) = 1 - \alpha,
$$
kde $K_{\alpha}$ je príslušný transformácia kvantilu Fisherovho $F$ rozdelenia v tvare $K_{\alpha} = \frac{p}{n - p} F_{p, n - p}(1 - \alpha)$ a $\mathcal{S}$ výberová variančná-kovariančná matica.

```{r}
Fquant <- (9/(20 - 9)) * qf(1 - 0.05, 9, 11)
meanVec <- apply(container.df[,2:10], 2, mean)
LB <- diag(9) %*% meanVec - sqrt(Fquant * diag(diag(9) %*% cov(container.df[,2:10])  %*% diag(9)))
UB <- diag(9) %*% meanVec + sqrt(Fquant * diag(diag(9) %*% cov(container.df[,2:10])  %*% diag(9)))

confInt <- data.frame(LB, UB, meanVec,  row.names = names(container.df)[2:10])
names(confInt) <- c("Lower Boundary", "Upper Boundary", "Mean Estimate")

confInt
```


<br><br>




#### <font color = "#880000"><b>3. Wilkovo Lambda rozdelnie</b></font>
Toto rozdelenie je odvodené z dvoch nezávislých náhodných matíc s  Wishartovým rozdelenim.
Vpodstate sa jedná o mnohorozmerné zobecnenie jednorozmerného Fisherovho F rozdelenia a používa sa hlavne 
k inferencii dvoch variačných-kovariačných matíc  (analogicky, ako sa jednorozmerné Fisherovo F rozdelenie používa 
k inferencii o dvoch parametroch rozptylu). 

Pre dve náhodné a vzájomné nezávislé náhodné matice 
$$
\mathbb{A} \sim W_{p}(\mathbb{I}, n) \quad \textrm{a} \quad \mathbb{B} \sim W_{p}(\mathbb{I}, m) 
$$

môžeme definovať náhodnú veličinu $\frac{|\mathbb{A}|}{|\mathbb{A} + \mathbb{B}|}$ ktorej rozdelenie je  Wilkovo Lambda, t.j., $\Lambda(p, n, m)$. 
Toto rozdelenie sa tiež používa v súvislosti s testom pomerom vierohodnosti. Analogickym spôsobom, ako funguje metóda analýzy rozptylu (ANOVA) na základe Fisherovho F rozdelenia, v mnohorozmernom prípade používame mnohorozmerné zobecnenie  (MANOVA), ktorá vychádza z Wilkovho Lambda rozdelenia. V programe R príkaz `manova()`. 



<br><br>

#### <font color = "#880000"><b>Ilustračný príklad</b></font>
<hr size="3" width="100%" noshade style="color:#440000" align="left" />
<ul>
<li>K ilustrácii Wilkovho Lambda rozdelenia využijeme priamo definíciu náhodnej veličiny s týmto rozdelením. 
Generovať budeme náhodné matice z Wishartovho rozdelenia a vhodnou transformáciou získavame požadované Wilkovo Lambda rozdelenie. <br><br>

```{r, fig.width = 10, fig.height=5}
set.seed(1234)
I <- cbind(c(1,0,0), c(0,1,0), c(0,0,1))
A <- rWishart(100, df = 10, I)
B <- rWishart(100, df = 4, I)

WilksLambda <- apply(A, 3, det)/apply(A + B, 3, det)

hist(WilksLambda, col = "lightblue", freq = F, xlab = "Wilk's Lambda", ylab = "Density", ylim = c(0, 4.5))
lines(density(WilksLambda), col = "red", lwd = 2)
```

</li>

</ul>
<hr size="3" width="100%" noshade style="color:#440000" align="left" />
<br><br>







<br><br>
<hr size="3" width="100%" noshade style="color:#440000" align="left" />
<hr size="3" width="100%" noshade style="color:#440000" align="left" /><br><br> 

## <font color = "#880000"><b>Domáca (samostatná) úloha </b></font>
### (Deadline: Piate cvičenie / 30.03.2021)

V R knižnici `SMSdata` sú k dispozícii rôzne datové súbory (viď podrobnejšie inštrukcie na stránke <a href = "https://dl1.cuni.cz/course/view.php?id=11882">Moodle UK</a>). Na vhodných datach (napr. z knižnice `SMSdata`, ale nie nutne, 
je možné využiť aj vlastné data) aplikujte jednovýberový, alebo dvojvýberový Hotellingov $T^2$ test pre vektor(y) stredných hodnôt. 

<ul>
<li>Nulovú a alternatívnu hypotézu v prvom rade vhodne motivujte pomocou obrázku (t.j. stručná exploratívna analýza).</li>

<li>Pomocou dostupných nástrojov v programe R, alebo pomocou vlastnej procedúry urobte príslušný štatistický test. Výsledky testu 
stručne interpretujte -- vysvetlite neštatistikovi/nematematikovi.</li>

<li>Doplňte vhodné intervaly spoľahlivosti (napr. konfidenčná oblasť, alebo súbor intervalov).</li> 

<li> Riešenie umiestnite na svoju webovú stránku, <font color = "red">najneskôr v utorok, 30.03.2021, do 14:00</font>.</li>
</ul>




</td></tr></table>
<br><br><br>