Esercizio - k-means a mano e gap statistic su cinque gruppi

Questa pagina non ha ancora la versione ripasso: qui sotto c'è il testo completo.

Testo.

1. Eseguire a mano il k-means con $K=2$ sui punti $x=[1,2,3,11,12,13,25]$ partendo dai centroidi $(1,25)$ e poi da $(2,12)$; confrontare le dispersioni $W$ e commentare.
2. Calcolare $W$ per $K=1,2,3$ e dire quale $K$ indica il metodo del gomito.
3. (Laboratorio) Su $1000$ punti generati da $5$ gruppi gaussiani (centri $(-10,-15)$, $(0,0)$, $(10,-15)$, $(25,-5)$, $(10,15)$, deviazione standard $1{,}5$) scegliere $K$ con il gomito e con la gap statistic.

Teoria usata: Clustering e k-meansIl clustering raggruppa osservazioni simili senza etichette, come preprocessing (un modello per ogni cluster) o come obiettivo (segmentazione clienti, organizzazione di documenti). K-means: si sceglie $K$, si inizializzano $K$ centroidi, si alterna assegnazione di ogni punto al centroide più vicino e aggiornamento di ogni centroide alla media dei suoi punti, fino a convergenza; minimizza $\mathrm{MSE}{\text{within}}=\frac1N\sum_k\sum{x_i\in C_k}|x_i-\mu_k|^2$ ma solo fino a un minimo locale, quindi dipende dall'inizializzazione. Il numero di cluster si sceglie col metodo del gomito (la dispersione cala sempre, si cerca dove rallenta) o con la gap statistic $\mathrm{Gap}(K)=E[\log W_K^{ref}]-\log W_K$ (si prende il più piccolo $K$ con $\mathrm{Gap}(K)\ge\mathrm{Gap}(K+1)-s_{K+1}$). Il clustering gerarchico agglomerativo parte da un cluster per punto e fonde i due più vicini (linkage single, complete, average, Ward) costruendo un dendrogramma; quello divisivo parte da un solo cluster. Programma di Telecomunicazioni: clustering.Clustering e k-means → (algoritmo, dispersione, gomito, gap statistic); la media come minimizzatore: Valore attesoIl valore atteso E[X] = Σ x p_X(x) è la media dei valori di X pesata con le loro probabilità (esiste se la serie converge assolutamente); per una funzione g vale E[g(X)] = Σ g(x) p_X(x) senza trovare la legge di g(X), ed E è lineare: E[aX + bY + c] = aE[X] + bE[Y] + c.Valore atteso →; logaritmi in Esponenziale e logaritmoLa funzione esponenziale a^x (base positiva diversa da 1) e la sua inversa, il logaritmo in base a, con grafici e proprietà.Esponenziale e logaritmo →.

1. k-means a mano ($K=2$)

La dispersione è $W=\sum_k\sum_{x\in C_k}(x-\mu_k)^2$ (in una dimensione la distanza è il valore assoluto della differenza).

Partenza A: centroidi $1$ e $25$.

• Assegnazione. Ogni punto va al centroide più vicino: $1,2,3,11,12,13$ sono più vicini a $1$ (per esempio $13$: distanza $12$ da $1$ e $12$ da $25$, pareggio: si assegna al primo; $12$: $11$ contro $13$), il punto $25$ va al secondo.
• Aggiornamento. $\mu_1=\frac{1+2+3+11+12+13}6=\frac{42}6=7$, $\mu_2=25$.
• Nuova assegnazione. Con centroidi $7$ e $25$: $13$ dista $6$ da $7$ e $12$ da $25$: va al primo; tutto è invariato: convergenza.
• $W=(1-7)^2+(2-7)^2+(3-7)^2+(11-7)^2+(12-7)^2+(13-7)^2+0=36+25+16+16+25+36=154$.

Partenza B: centroidi $2$ e $12$.

• Assegnazione. $1,2,3$ vanno a $2$ (distanze $1,0,1$ contro $11,10,9$); $11$: distanza $9$ da $2$ e $1$ da $12$: va a $12$; così $12$, $13$ e $25$.
• Aggiornamento. $\mu_1=\frac{1+2+3}3=2$, $\mu_2=\frac{11+12+13+25}4=\frac{61}4=15{,}25$.
• Nuova assegnazione. $11$: dista $9$ da $2$ e $4{,}25$ da $15{,}25$: resta nel secondo; $13$: $11$ contro $2{,}25$: secondo. Nulla cambia: convergenza.
• $W=(1-2)^2+0+(3-2)^2+(11-15{,}25)^2+(12-15{,}25)^2+(13-15{,}25)^2+(25-15{,}25)^2=2+18{,}06+10{,}56+5{,}06+95{,}06=130{,}75$.

Commento. Due inizializzazioni portano a due minimi locali diversi: $W=154$ (il punto $25$ da solo) e $W=130{,}75$ (migliore). Il k-means non garantisce il minimo globale: si ripete da più inizializzazioni e si tiene la soluzione con $W$ minore.

2. Gomito

• $K=1$: unico centroide, la media $\mu=\frac{1+2+3+11+12+13+25}7=\frac{67}7=9{,}571$ e $W=\sum(x-\mu)^2=431{,}71$ (la varianza totale per $n$).
• $K=2$ migliore: $W=130{,}75$.
• $K=3$ con centroidi iniziali $(2,12,25)$: cluster $\{1,2,3\}$, $\{11,12,13\}$, $\{25\}$ con $\mu=2,12,25$ e $W=(1+0+1)+(1+0+1)+0=4$.

$W$ scende da $431{,}7$ a $130{,}8$ a $4{,}0$: il calo da $K=2$ a $K=3$ è ancora enorme e un $K=4$ aggiungerebbe quasi nulla (con $K=7=n$, $W=0$ sempre). Il gomito è a $K=3$, che corrisponde ai tre gruppi visibili ($\{1,2,3\}$, $\{11,12,13\}$, $\{25\}$): anche un punto isolato è un cluster.

3. Cinque gruppi: gomito e gap statistic

Gomito. Errore medio entro i cluster $W_K/N$ per $K=1,\dots,8$: $264{,}5;\ 136{,}7;\ 67{,}8;\ 36{,}0;\ \mathbf{4{,}28};\ 3{,}95;\ 3{,}66;\ 3{,}45$. Il calo da $K=4$ a $K=5$ è enorme ($36{,}0\to4{,}28$), poi quasi nulla: gomito a $K=5$, il numero vero di gruppi. (Per ogni $K$ il k-means è stato lanciato più volte tenendo il $W$ minimo.)

Gap statistic.

1. Per ogni $K=1,\dots,7$: $W_K$ sui dati reali.
2. Si generano 10 dataset di riferimento con $1000$ punti uniformi nel rettangolo dei dati ($x\in[\min,\max]$ per ciascuna coordinata) e per ciascuno si calcola $W_K^{\text{ref}}$.
3. $\mathrm{Gap}(K)=\overline{\log W_K^{\text{ref}}}-\log W_K$ e $s_K=$ deviazione standard dei $\log W_K^{\text{ref}}$.

Risultati (una replica): $\mathrm{Gap}=0{,}02;\ 0{,}13;\ 0{,}41;\ 0{,}62;\ \mathbf{2{,}56};\ 2{,}45;\ 2{,}37$ per $K=1,\dots,7$ ($s_K\approx0{,}02$), cioè lo stesso andamento del notebook ($0{,}04;\ 0{,}16;\ 0{,}41;\ 0{,}64;\ 2{,}60;\ 2{,}47$). Scelta: il più piccolo $K$ con $\mathrm{Gap}(K)\ge\mathrm{Gap}(K+1)-s_{K+1}$:

• $K=1$: $0{,}02\ge0{,}13-0{,}02$? No. $K=2$: $0{,}13\ge0{,}41-0{,}02$? No. $K=3$: $0{,}41\ge0{,}62-0{,}02$? No. $K=4$: $0{,}62\ge2{,}56-0{,}02$? No.
• $K=5$: $2{,}56\ge2{,}45-0{,}02$? Sì: $K=5$.

Gomito e gap statistic concordano. Il vantaggio della gap statistic è un criterio oggettivo (il gomito è un criterio visivo): compara con dati senza struttura.

Codice

def wk(X, labels, centroids):                                # dispersione entro i cluster
    return sum(((X[labels == j] - centroids[j]) ** 2).sum() for j in range(len(centroids)))
def gap(X, K, n_ref=10):
    lo, hi = X.min(0), X.max(0)
    logw = np.log(best_of_runs(X, K))                        # k-means ripetuto, tiene W minimo
    ref = [np.log(best_of_runs(np.random.uniform(lo, hi, X.shape), K)) for _ in range(n_ref)]
    return np.mean(ref) - logw, np.std(ref)