26. La verosimiglianza#
Obiettivi di apprendimento
Dopo aver completato questo capitolo, sarai in grado di:
Comprendere il concetto di verosimiglianza e il suo ruolo nella dei parametri.
Generare grafici della funzione di verosimiglianza binomiale.
Generare grafici della funzione di verosimiglianza del modello gaussiano.
Interpretare i grafici della funzione di verosimiglianza.
Comprendere il concetto di stima di massima verosimiglianza.
import numpy as np
import pandas as pd
import matplotlib as mpl
import matplotlib.pyplot as plt
import seaborn as sns
# from scipy.stats import binom
import scipy.stats as st
import math
import arviz as az
%matplotlib inline
%config InlineBackend.figure_format = 'retina'
# Initialize random number generator
RANDOM_SEED = 8927
rng = np.random.default_rng(RANDOM_SEED)
plt.style.use("bmh")
plt.rcParams["figure.figsize"] = [10, 6]
plt.rcParams["figure.dpi"] = 100
plt.rcParams["figure.facecolor"] = "white"
sns.set_theme(palette="colorblind")
%load_ext autoreload
%autoreload 2
%config InlineBackend.figure_format = "svg"
26.1. La funzione di verosimiglianza#
Come spesso accade in statistica, il nostro obiettivo è stimare un parametro del modello utilizzando i dati a nostra disposizione. Nel contesto frequentista, una delle metodologie più comuni per stimare i parametri è la stima di massima verosimiglianza. Questo approccio si basa sulla definizione della funzione di verosimiglianza, che indica quanto i nostri dati siano verosimili per un determinato valore del parametro. La funzione di verosimiglianza è utilizzata sia nel paradigma bayesiano che in quello frequentista.
La funzione di probabilità (o densità) e la funzione di verosimiglianza sono formalmente identiche ma hanno scopi diversi. Nella funzione di probabilità, i parametri
La funzione di probabilità è indicata come
Per la funzione di verosimiglianza, è consuetudine utilizzare la lettera
In questa formula,
La funzione di verosimiglianza, a differenza di una funzione di densità di probabilità, non è normalizzata per avere un’area unitaria. Ciò significa che la funzione di verosimiglianza fornisce solo informazioni relative sulla plausibilità dei diversi valori del parametro
26.2. Modello binomiale#
Facciamo un esempio relativo alla distribuzione binomiale. Consideriamo un esperimento ripetuto
dove
La funzione di verosimiglianza rappresenta la plausibilità relativa di osservare i dati
dove possiamo trascurare il coefficiente binomiale in quanto non dipende da
Per fare un esempio pratico, consideriamo la ricerca di Zetsche et al. [ZBR19]. Questi ricercatori hanno trovato che, su 30 pazienti clinicamente depressi, 23 manifestavano delle aspettative distorte negativamente relativamente al loro umore futuro. Per i dati di Zetsche et al. [ZBR19], la funzione di verosimiglianza corrisponde dunque alla funzione binomiale di parametro
Per costruire la funzione di verosimiglianza, è necessario applicare l’eq. (26.1) diverse volte, variando ogni volta il valore di
n = 30
y = 23
Creiamo ora i possibili valori del parametro
theta = np.linspace(0.0, 1.0, num=100)
print(theta)
[0. 0.01010101 0.02020202 0.03030303 0.04040404 0.05050505
0.06060606 0.07070707 0.08080808 0.09090909 0.1010101 0.11111111
0.12121212 0.13131313 0.14141414 0.15151515 0.16161616 0.17171717
0.18181818 0.19191919 0.2020202 0.21212121 0.22222222 0.23232323
0.24242424 0.25252525 0.26262626 0.27272727 0.28282828 0.29292929
0.3030303 0.31313131 0.32323232 0.33333333 0.34343434 0.35353535
0.36363636 0.37373737 0.38383838 0.39393939 0.4040404 0.41414141
0.42424242 0.43434343 0.44444444 0.45454545 0.46464646 0.47474747
0.48484848 0.49494949 0.50505051 0.51515152 0.52525253 0.53535354
0.54545455 0.55555556 0.56565657 0.57575758 0.58585859 0.5959596
0.60606061 0.61616162 0.62626263 0.63636364 0.64646465 0.65656566
0.66666667 0.67676768 0.68686869 0.6969697 0.70707071 0.71717172
0.72727273 0.73737374 0.74747475 0.75757576 0.76767677 0.77777778
0.78787879 0.7979798 0.80808081 0.81818182 0.82828283 0.83838384
0.84848485 0.85858586 0.86868687 0.87878788 0.88888889 0.8989899
0.90909091 0.91919192 0.92929293 0.93939394 0.94949495 0.95959596
0.96969697 0.97979798 0.98989899 1. ]
Per esempio, ponendo
st.binom.pmf(y, n, 0.1)
9.7371682902e-18
Ponendo
st.binom.pmf(y, n, 0.2)
3.58141723492221e-11
Se ripetiamo questo processo 100 volte, una volta per ciascuno dei valori
def like(r, n, theta):
return math.comb(n, r) * theta**r * (1 - theta) ** (n - r)
La curva che interpola i punti ottenuti è la funzione di verosimiglianza, come indicato dalla figura seguente.
plt.plot(theta, like(r=y, n=n, theta=theta), "r-")
plt.title("Funzione di verosimiglianza")
plt.xlabel("Valore della variabile casuale theta [0, 1]")
plt.ylabel("Verosimiglianza")
Text(0, 0.5, 'Verosimiglianza')
26.2.1. Interpretazione#
Si osservi che la verosimiglianza,
Per determinare il valore di argmax
di NumPy. Una volta ottenuto l’indice corrispondente, è possibile recuperare il valore di theta
. In questo modo, otterremo il valore di
l = like(r=y, n=n, theta=theta)
l.argmax()
76
theta[76]
0.7676767676767677
È importante notare che, invece di utilizzare la funzione like()
che abbiamo definito precedentemente per motivi didattici, è possibile ottenere lo stesso risultato utilizzando in modo equivalente la funzione binom.pmf()
.
plt.plot(theta, st.binom.pmf(y, n, theta), "r-")
plt.title("Funzione di verosimiglianza")
plt.xlabel("Valore della variabile casuale theta [0, 1]")
plt.ylabel("Verosimiglianza")
Text(0, 0.5, 'Verosimiglianza')
26.2.2. La log-verosimiglianza#
Dal punto di vista pratico risulta più conveniente utilizzare, al posto della funzione di verosimiglianza, il suo logaritmo naturale, ovvero la funzione di log-verosimiglianza:
Poiché il logaritmo è una funzione strettamente crescente (usualmente si considera il logaritmo naturale), allora
Per le proprietà del logaritmo, la funzione nucleo di log-verosimiglianza della binomiale è
È importante sottolineare che non è strettamente necessario lavorare con i logaritmi, ma è altamente raccomandato. La ragione principale è che i valori di verosimiglianza, che rappresentano il prodotto di numeri di probabilità molto piccoli, possono diventare estremamente piccoli, anche dell’ordine di
Svolgiamo nuovamente il problema precedente usando la log-verosimiglianza per trovare il massimo della funzione di log-verosimiglianza. Ora utilizziamo la funzione binom.logpmf()
.
La funzione di log-verosimiglianza è rappresentata nella figura successiva.
n = 30
r = 23
plt.plot(theta, st.binom.logpmf(y, n, theta), "r-")
plt.title("Funzione di log-verosimiglianza")
plt.xlabel("Valore della variabile casuale theta [0, 1]")
plt.ylabel("Log-verosimiglianza")
Text(0, 0.5, 'Log-verosimiglianza')
Il risultato replica quello trovato in precedenza con la funzione di verosimiglianza.
ll = st.binom.logpmf(y, n, theta)
ll.argmax()
76
theta[76]
0.7676767676767677
26.2.3. Verosimiglianza congiunta di n campioni iid da una distribuzione binomiale#
Nell’esempio precedente, abbiamo considerato la verosimiglianza per una singola osservazione binomiale. Per generalizzare quanto abbiamo finora, consideriamo nuovamente la distribuzione di probabilità binomiale
Sostituendo i dati
26.3. Modello gaussiano#
Estendiamo ora la discussione precedente al caso gaussiano. Inizieremo calcolando la verosimiglianza per il modello gaussiano dato un singolo valore puntuale
Definiamo la funzione di verosimiglianza per la distribuzione gaussiana con i parametri
26.3.1. Una singola osservazione#
Per ottenere la funzione di verosimiglianza dell’eq. (26.3), esaminiamo inizialmente il caso in cui i dati corrispondono ad una singola osservazione
Supponiamo che la variabile casuale
y = 114
L’eq. (26.3) dipende dai parametri
mu = np.linspace(70.0, 160.0, num=1000)
Poiché stiamo esaminando 1000 potenziali valori per il parametro
In ciascun passo dell’esercizio inseriremo nell’eq. (26.3)
il singolo valore
considerato (che viene mantenuto costante),il valore
assunto noto (anch’esso costante),uno alla volta ciascuno dei valori
che abbiamo definito.
Quindi, nelle 1000 applicazioni dell’eq. (26.3), il valore
In Python, la distribuzione gaussiana può essere implementata attraverso la funzione norm.pdf()
. Tale funzione richiede tre argomenti: il valore o il vettore di valori norm.pdf()
consente di calcolare la densità di probabilità gaussiana per ogni valore di
Applicando la funzione norm.pdf()
1000 volte, una volta per ciascuno dei valori
f_mu = st.norm.pdf(y, loc=mu, scale=15)
La funzione di verosimiglianza è la curva che interpola i punti
plt.plot(mu, f_mu, "r-")
plt.title("Funzione di verosimiglianza")
plt.xlabel("Valore della variabile casuale mu [70, 160]")
plt.ylabel("Verosimiglianza")
plt.xlim([70, 160])
(70.0, 160.0)
La funzione di verosimiglianza così trovata ha la forma della distribuzione Gaussiana. Nel caso di una singola osservazione, ma solo in questo caso, ha anche un’area unitaria. Per l’esempio presente, la moda della funzione di verosimiglianza è 114.
l = st.norm.pdf(y, loc=mu, scale=15)
mu[l.argmax()]
113.96396396396396
26.3.2. Un campione casuale indipendente da una distribuzione gaussiana#
Consideriamo ora il caso più generale di un campione casuale indipendente di
Se le variabili casuali
laddove
Per chiarire l’eq. (26.3), consideriamo un esempio che utilizza come dati i valori BDI-II dei trenta soggetti del campione clinico di Zetsche et al. [ZBR19].
y = [
26,
35,
30,
25,
44,
30,
33,
43,
22,
43,
24,
19,
39,
31,
25,
28,
35,
30,
26,
31,
41,
36,
26,
35,
33,
28,
27,
34,
27,
22,
]
Ci poniamo l’obiettivo di creare la funzione di verosimiglianza per questi dati, supponendo di sapere (in base ai risultati di ricerche precedenti) che i punteggi BDI-II si distribuiscono secondo la legge Normale e supponendo
true_sigma = np.std(y)
true_sigma
6.495810615739622
Per la prima osservazione del campione (
Se consideriamo tutte le osservazioni del campione, la densità congiunta è il prodotto delle densità delle singole osservazioni:
Utilizzando i dati del campione e assumendo
Il valore
È più conveniente svolgere i calcoli usando il logaritmo della verosimiglianza. In Python definiamo la funzione di log-verosimiglianza, log_likelihood()
, che prende come argomenti y
, mu
e sigma = 6.50
:
def log_likelihood(y, mu, sigma=true_sigma):
return np.sum(st.norm.logpdf(y, loc=mu, scale=true_sigma))
Consideriamo, ad esempio, il valore
bar_y = np.mean(y)
print(bar_y)
30.933333333333334
L’ordinata della funzione di log-verosimiglianza in corrispondenza di
print(log_likelihood(y, 30.93, sigma=true_sigma))
-98.70288339960591
Troviamo ora i valori della log-verosimiglianza per ciascuno dei 1000 valori mu
.
mu = np.linspace(np.mean(y) - 2 * np.std(y), np.mean(y) + 2 * np.std(y), num=1000)
Troviamo il valore dell’ordinata della funzione di log-verosimiglianza in corrispondenza di ciascuno dei 1000 valori mu
che abbiamo definito.
ll = [log_likelihood(y, mu_val, true_sigma) for mu_val in mu]
Nel caso di un solo parametro sconosciuto (nel caso presente, mu
, ll
). Tale funzione descrive la credibilità relativa che può essere attribuita ai valori del parametro
plt.plot(mu, ll, "r-")
plt.title("Funzione di log-verosimiglianza")
plt.xlabel("Valore della variabile casuale mu")
plt.ylabel("Log-verosimiglianza")
plt.axvline(x=np.mean(y), color="k", alpha=0.2, ls="--")
<matplotlib.lines.Line2D at 0x123f6bed0>
Il valore
Il massimo della funzione di log-verosimiglianza, ovvero 30.93 per l’esempio in discussione, è identico alla media dei dati campionari.
26.4. Derivazione formale#
Le stime di massima verosimiglianza per i parametri
La massimizzazione della probabilità dei nostri dati può essere scritta come:
Espandendo i nostri parametri, abbiamo
Impostando l’espressione sopra uguale a zero, otteniamo
Quindi
In altre parole, la stima di massima verosimiglianza del parametro
26.5. Commenti e considerazioni finali#
In conclusione, la funzione di verosimiglianza costituisce un ponte fondamentale tra i parametri di un modello statistico e i dati osservati. Attraverso questa funzione, otteniamo informazioni preziose sulla plausibilità dei dati alla luce dei parametri del modello. La sua formulazione coinvolge tre elementi chiave: il modello statistico che descrive la generazione dei dati, i possibili valori dei parametri del modello e i dati effettivamente osservati. La funzione di verosimiglianza svolge un ruolo centrale nell’inferenza statistica, consentendoci di valutare la compatibilità tra i dati osservati e i diversi valori dei parametri del modello. La sua comprensione e il suo corretto utilizzo sono essenziali per il processo di analisi dei dati e l’interpretazione dei risultati.
26.6. Watermark#
%load_ext watermark
%watermark -n -u -v -iv -w
Last updated: Sat Jun 17 2023
Python implementation: CPython
Python version : 3.11.3
IPython version : 8.12.0
numpy : 1.24.3
arviz : 0.15.1
pandas : 1.5.3
seaborn : 0.12.2
matplotlib: 3.7.1
scipy : 1.10.1
Watermark: 2.3.1