5.5 Stime nel caso vettoriale

I risultati della sezione precedente si possono estendere al caso vettoriale, ossia di \(n\) osservazioni di variabili aleatorie \(X_1\), …, \(X_n \in \R^d\), tutte indipendenti tra loro e ciascuna con densità gaussiana vettoriale di parametri comuni \(\mathcal{N}(m, \Sigma)\). Posta per brevità \(X = (X_1, \ldots, X_n) \in \R^{nd}\), funzione di verosimiglianza dei parametri di media e varianza associata all’osservazione di \(X=x = (x_i)_{i=1}^n\), si scrive, a meno di costanti moltiplicative, \[ L(m,\Sigma; x) \propto \exp\bra{ - \frac 1 2 \sum_{i=1}^n (x_i-m)\cdot \Sigma^{-1} (x_i-m)} \frac{1}{\bra{ \det \Sigma}^{n/2}}.\]

I calcoli analitici, già nel caso della stima di massima verosimiglianza, sono meno agevoli e ne riportiamo solamente il risultato. Precisamente, 1. se la varianza \(\Sigma = \Sigma_0\) è nota (rispetto all’informazione prima di osservare le \(X_i\)), la stima di massima verosimiglianza per il parametro di media è la media campionaria \[ m_{{\mle}} = \bar{x} = \frac{1}{n} \sum_{i=1}^n x_i\] (qualsiasi sia \(\Sigma_0\)).

  1. se il parametro di media \(m = m_0\in \R^d\) è noto (rispetto all’informazione a priori), allora la stima di massima verosimiglianza per la covarianza è \[ \Sigma_{{\mle}} = \frac 1 n \sum_{i=1}^n(x_i-m_0)(x_i-m_0)^T,\] dove \(T\) indica l’operazione di trasposizione (quindi il prodotto righe per colonne risulta in una matrice \(d\times d\)); più esplicitamente, la stima della covarianza tra la componente \(j\) e \(k\) è \[ (\Sigma_{{\mle}})_{jk} = \frac 1 n \sum_{i=1}^n(x_i-m_0)_j (x_i-m_0)_k.\] Mettendo insieme i due risultati sopra, si ottiene analogamente al caso reale che la stima (congiunta) di massima verosimiglianza per \((m,\Sigma)\) è data dalla media e dalla covarianza campionarie: \[ m_{{\mle}} = \bar{x} = \frac{1}{n} \sum_{i=1}^n x_i, \quad \Sigma_{{\mle}} = \frac 1 n \sum_{i=1}^n(x_i-\bar{x})(x_i-\bar{x})^T.\] Osserviamo che \(\Sigma_{\mle}\) è una matrice simmetrica e semi-definita positiva. La si può anche interpretare come la matrice di covarianza della variabile aleatoria vettoriale che sceglie uno degli \(n\) valori osservati con probabilità uniforme discreta.

Esempio 5.4 Torniamo al dataset Iris e usiamo le funzioni summary() e cov() per calcolare la media e la covarianza campionaria delle prime \(4\) colonne (escludiamo naturalmente quella contenente il nome della specie). La funzione summary() indica anche mediana e quartili, quindi selezioniamo solamente la media (che corrisponde alla quarta riga).

ind_colonne <- summary(iris[, 1:4])

ind_colonne[4, ]
##      Sepal.Length       Sepal.Width      Petal.Length       Petal.Width 
## "Mean   :5.843  " "Mean   :3.057  " "Mean   :3.758  " "Mean   :1.199  "
cov(iris[, 1:4])
##              Sepal.Length Sepal.Width Petal.Length Petal.Width
## Sepal.Length    0.6856935  -0.0424340    1.2743154   0.5162707
## Sepal.Width    -0.0424340   0.1899794   -0.3296564  -0.1216394
## Petal.Length    1.2743154  -0.3296564    3.1162779   1.2956094
## Petal.Width     0.5162707  -0.1216394    1.2956094   0.5810063

È utile anche considerare la matrice delle correlazioni, in cui al posto delle covarianze è calcolato il coefficiente di correlazione campionario, \[ \bar{\rho}_{jk} = \frac{ \Sigma_{jk}}{ \sqrt{ \Sigma_{jj} \Sigma_{kk}}},\] che è sempre compreso tra \(-1\) e \(1\) (segue dal fatto che la matrice \(\Sigma\) è semi-definita positiva). Il comando in questo caso è cor().

cor(iris[, 1:4])
##              Sepal.Length Sepal.Width Petal.Length Petal.Width
## Sepal.Length    1.0000000  -0.1175698    0.8717538   0.8179411
## Sepal.Width    -0.1175698   1.0000000   -0.4284401  -0.3661259
## Petal.Length    0.8717538  -0.4284401    1.0000000   0.9628654
## Petal.Width     0.8179411  -0.3661259    0.9628654   1.0000000

Per visualizzare la correlazione si può usare un correlogramma, in cui i valori sono accompagnati (o addirittura) sostituiti da colori opportuni. Questo è particolarmente utile se le componenti del vettore sono in gran numero. In R si può usare la funzione corrplot() dalla libreria corrplot (da installare la prima volta tramite il comando install.packages('corrplot')). Poiché la matrice è quadrata basta rappresentarne una parte triangolare, ad esempio superiore.

library(corrplot)

corrplot(cor(iris[, 1:4]), method = "color",
  type = "upper", tl.col = "black", tl.srt = 45)
Correlogramma del dataset Iris

Figura 5.18: Correlogramma del dataset Iris

Il comando plot() applicato direttamente ai dati fornisce invece un diagramma a dispersione (in inglese scatter plot) di tutte le possibili coppie.

plot(iris[, 1:4], col = miei_colori[2], pch = 16)

Tralasciamo invece l’approccio bayesiano, più complesso.

5.5.1 Esercizi

Esercizio 5.8 Ripetere le osservazioni fatte sopra considerando separatamente ciascuna specie (selezionare le prime 50 righe per la prima specie, le ulteriori 50 per la seconda e le ultime 50 per la terza).

Esercizio 5.9 Si scarichino dalla pagina del progetto Pageview stats, https://pageviews.toolforge.org/, i dati relativi alle visualizzazioni di \(4\) pagine di Wikipedia che si possano ritenere correlate (ad esempio si parta da una pagina e si considerino poi il primo collegamento da essa, oppure i primi due, e si ripeta). Si calcoli la correlazione empirica e la si visualizzi graficamente.