Capítulo 2

Mezcla de Gaussianas y EM

2.1 Gaussian Mixture Models

Un GMM modela datos como una suma ponderada de K Gaussianas:

p(x)=k=1KπkN(xμk,Σk)p(\mathbf{x}) = \sum_{k=1}^K \pi_k \, \mathcal{N}(\mathbf{x} \mid \boldsymbol{\mu}_k, \boldsymbol{\Sigma}_k)

Donde πk\pi_k son los pesos (mezcla), con πk=1\sum \pi_k = 1. Es un modelo de variables latentes: cada punto pertenece a un cluster que no observamos.

2.2 El Algoritmo EM

Expectation-Maximization estima parámetros con variables latentes iterando dos pasos:

E-step

Estima la responsabilidad (probabilidad) de cada cluster para cada punto.

rnk=πkN(xnμk,Σk)jπjN(xnμj,Σj)r_{nk} = \frac{\pi_k \mathcal{N}(\mathbf{x}_n \mid \boldsymbol{\mu}_k, \boldsymbol{\Sigma}_k)}{\sum_j \pi_j \mathcal{N}(\mathbf{x}_n \mid \boldsymbol{\mu}_j, \boldsymbol{\Sigma}_j)}

M-step

Maximiza la verosimilitud con las responsabilidades fijas. Reestima μ, Σ, π.

μk=nrnkxnnrnk\boldsymbol{\mu}_k = \frac{\sum_n r_{nk} \mathbf{x}_n}{\sum_n r_{nk}}

2.3 Visualización Interactiva

Haz clic para agregar puntos. Observa cómo EM ajusta las Gaussianas automáticamente.

Iter: 0 | Puntos: 0

Haz clic en el canvas para agregar puntos. Corre EM para ajustar las Gaussianas.

2.4 Implementación desde Cero

python
import numpy as np
from scipy.stats import multivariate_normal

class GMM:
    def __init__(self, K=3, max_iter=100):
        self.K = K
        self.max_iter = max_iter

    def fit(self, X):
        n, d = X.shape
        # Inicializar
        self.pi = np.ones(self.K) / self.K
        idx = np.random.choice(n, self.K, replace=False)
        self.mu = X[idx]
        self.Sigma = np.array([np.eye(d) * X.var() for _ in range(self.K)])
        self.resp = np.zeros((n, self.K))

        for _ in range(self.max_iter):
            # E-step: responsabilidades
            for k in range(self.K):
                self.resp[:, k] = self.pi[k] * multivariate_normal.pdf(
                    X, self.mu[k], self.Sigma[k])
            self.resp /= self.resp.sum(axis=1, keepdims=True) + 1e-12

            # M-step: reestimar
            Nk = self.resp.sum(axis=0)
            self.pi = Nk / n
            self.mu = (self.resp.T @ X) / Nk[:, None]
            for k in range(self.K):
                diff = X - self.mu[k]
                self.Sigma[k] = (diff * self.resp[:, k:k+1]).T @ diff / Nk[k]

    def predict(self, X):
        resp = np.zeros((len(X), self.K))
        for k in range(self.K):
            resp[:, k] = self.pi[k] * multivariate_normal.pdf(
                X, self.mu[k], self.Sigma[k])
        return resp.argmax(axis=1)

2.5 Aplicación: Clustering

python
from sklearn.datasets import make_blobs

X, y_true = make_blobs(n_samples=500, centers=3,
                       cluster_std=1.0, random_state=42)
gmm = GMM(K=3)
gmm.fit(X)
y_pred = gmm.predict(X)

# Comparar con sklearn
from sklearn.mixture import GaussianMixture
skgmm = GaussianMixture(n_components=3)
skgmm.fit(X)
print(f"Sklearn score: {skgmm.score(X):.2f}")

GMM es un modelo generativo: podemos muestrear datos sintéticos del modelo entrenado.

python
# Generar datos sintéticos
z = np.random.choice(self.K, size=100, p=self.pi)
samples = np.array([
    np.random.multivariate_normal(self.mu[k], self.Sigma[k])
    for k in z])