Reducir dimensiones: PCA, t-SNE y UMAP

a) Sobre la recta de 0°, la sombra de cada punto es su coordenada x: 2, −1 y −1. La varianza de la sombra es (4 + 1 + 1) / 3 = 2. Lo que se pierde es la coordenada y (1, 0 y −1), así que el error es (1 + 0 + 1) / 3 = 2/3 ≈ 0,67. Sobre la recta de 90° se intercambian: varianza 2/3 y error 2. En los dos casos la suma es 8/3 ≈ 2,67, la varianza total.

b) Los elementos de la diagonal son las varianzas de cada columna, 2 y 2/3. Fuera de la diagonal va el promedio de los productos x·y: (2·1 + (−1)·0 + (−1)·(−1)) / 3 = 3/3 = 1.

c) Suman lo mismo que la diagonal, 8/3: la varianza total no cambia al girar los ejes, solo se reparte distinto.

P = np.array([[2, 1], [-1, 0], [-1, -1]])
print("0°:  varianza y error", np.round(sombra(P, 0), 2))
print("90°: varianza y error", np.round(sombra(P, 90), 2))
C3 = P.T @ P / 3
valores3, vectores3 = np.linalg.eigh(C3)
print("Valores propios:", valores3.round(3),
      " suma:", valores3.sum().round(3))
v = vectores3[:, 1]
print(f"Ángulo de la primera componente: "
      f"{np.rad2deg(np.arctan2(v[1], v[0])) % 180:.1f}°")

0°:  varianza y error [2.   0.67]
90°: varianza y error [0.67 2.  ]
Valores propios: [0.131 2.535]  suma: 2.667
Ángulo de la primera componente: 28.2°

La primera componente se lleva casi toda la varianza: es una recta inclinada unos 28°, entre el eje horizontal (que ya era bueno) y la diagonal.

Xd_c = Xd - Xd.mean(axis=0)
U, S, Vt = np.linalg.svd(Xd_c, full_matrices=False)
varianzas_svd = S ** 2 / len(Xd)
print("Proporción con SVD:  ",
      (varianzas_svd / varianzas_svd.sum())[:5].round(3))

Z2, comps2, proporcion2 = pca(Xd, 2)
print("Proporción con eigh: ", proporcion2[:5].round(3))
print("Mismas componentes:", np.allclose(np.abs(Vt[:2]), np.abs(comps2)))

Proporción con SVD:   [0.149 0.136 0.118 0.084 0.058]
Proporción con eigh:  [0.149 0.136 0.118 0.084 0.058]
Mismas componentes: True

Dan lo mismo, salvo el signo de cada componente, que es arbitrario. scikit-learn usa SVD por dentro porque es más estable numéricamente y no necesita formar la matriz de covarianza, que con muchas columnas es enorme.

from sklearn.datasets import load_wine
from sklearn.preprocessing import StandardScaler

vinos = load_wine()
nombres = np.array(vinos.feature_names)
desviaciones = vinos.data.std(axis=0)
for j in [0, 4, 12]:
    print(f"Desviación de {nombres[j]}: {desviaciones[j]:.1f}")

for titulo, datos in [("Sin escalar", vinos.data),
                      ("Estandarizado",
                       StandardScaler().fit_transform(vinos.data))]:
    p = PCA().fit(datos)
    principal = nombres[np.argmax(np.abs(p.components_[0]))]
    print(f"{titulo}: la 1.ª componente explica "
          f"{p.explained_variance_ratio_[0]:.1%}; pesa más {principal}")

Desviación de alcohol: 0.8
Desviación de magnesium: 14.2
Desviación de proline: 314.0
Sin escalar: la 1.ª componente explica 99.8%; pesa más proline
Estandarizado: la 1.ª componente explica 36.2%; pesa más flavanoids

Sin escalar, la columna proline tiene una desviación de unos 300, mientras la mayoría de las demás andan entre 0,1 y 15. La varianza es desviación al cuadrado, así que proline se lleva casi toda la varianza y la «primera componente» es básicamente esa columna sola: PCA no descubrió nada, solo midió unidades. Estandarizado, cada columna aporta lo mismo, la primera componente explica alrededor de un tercio y mezcla varias mediciones químicas. Regla: si las columnas están en unidades distintas, estandariza antes de PCA (en los dígitos no hizo falta porque todos los píxeles van de 0 a 16).

from sklearn.model_selection import train_test_split
from sklearn.pipeline import make_pipeline
from sklearn.neighbors import KNeighborsClassifier

X_entr, X_val, y_entr, y_val = train_test_split(
    Xd, yd, test_size=0.25, random_state=0, stratify=yd)

for k in [1, 2, 3, 5, 10, 20, 40, 64]:
    tubo = make_pipeline(PCA(n_components=k, random_state=0),
                         KNeighborsClassifier(n_neighbors=5))
    tubo.fit(X_entr, y_entr)
    print(f"{k:2d} componentes: exactitud {tubo.score(X_val, y_val):.3f}")

 1 componentes: exactitud 0.356
 2 componentes: exactitud 0.618
 3 componentes: exactitud 0.749
 5 componentes: exactitud 0.913
10 componentes: exactitud 0.967
20 componentes: exactitud 0.976
40 componentes: exactitud 0.980
64 componentes: exactitud 0.980

Con 1 o 2 componentes el k-NN se confunde mucho (es lo que viste en el panel de PCA: los dígitos encimados). Con 10 ya pasa del 96 % y desde unas 20 se queda prácticamente igual que con los 64 píxeles. Es decir: con un tercio de los números se clasifica igual de bien, y el k-NN, que compara distancias con todos los ejemplos, trabaja tres veces menos. El PCA va dentro del pipeline para que se ajuste solo con los datos de entrenamiento (tema 10): ajustarlo con todos sería una pequeña fuga de datos.

La 2 es la que sí vale: t-SNE se esfuerza justamente en que los vecinos cercanos sigan cerca. Lo comprobaste en la teoría: el vecino más cercano en el mapa es del mismo dígito en el 98,8 % de los casos.

Las otras tres son falsas:

1. Las distancias entre grupos no significan mucho. t-SNE cuida los vecindarios, no la distancia entre grupos lejanos. Mídelo en las 64 dimensiones, con el dígito promedio de cada grupo:

centros = np.array([Xd[yd == d].mean(axis=0) for d in range(10)])
dist_64 = np.linalg.norm(centros - centros[0], axis=1)
centros_mapa = np.array([np.median(Z_tsne[yd == d], axis=0)
                         for d in range(10)])
dist_mapa = np.linalg.norm(centros_mapa - centros_mapa[0], axis=1)
print("Del 0 al más cercano, en 64 dimensiones:", np.argsort(dist_64)[1:4])
print("Del 0 al más cercano, en el mapa t-SNE: ",
      np.argsort(dist_mapa)[1:4])

Del 0 al más cercano, en 64 dimensiones: [9 8 5]
Del 0 al más cercano, en el mapa t-SNE:  [5 9 6]

En los píxeles, el dígito más parecido al 0 promedio es el 9, no el 5.

3. Los ejes no tienen nombre. Si giras o reflejas el mapa entero, es igual de bueno: t-SNE solo mira qué puntos quedan cerca de cuáles. En PCA, en cambio, cada eje sí es una dirección concreta de los datos originales.

4. Los dos mapas son vistas válidas. La perplejidad dice cuántos vecinos mira cada punto: con 5, cada punto solo se preocupa por sus vecinos más inmediatos y los grupos se desarman en hilos y grumitos. Por eso se prueban varios valores y solo se confía en lo que se repite (como que el 1 se parte en dos grupos y que unos cuantos unos quedan cerca del 2 o sueltos).

Este código hace los pasos 1 a 4 sin gráficos; tus números deberían coincidir.

from sklearn.datasets import load_sample_image

foto = load_sample_image("china.jpg").mean(axis=2) / 255
alto, ancho = foto.shape
print("Forma:", foto.shape)

for k in [5, 10, 20, 50, 100]:
    p = PCA(n_components=k, random_state=0).fit(foto)
    reconstruida = p.inverse_transform(p.transform(foto))
    guardados = k * alto + k * ancho + ancho
    error = np.mean((foto - reconstruida) ** 2)
    print(f"k = {k:3d}:  tamaño {guardados / foto.size:5.1%}"
          f"   error {error:.4f}")

Forma: (427, 640)
k =   5:  tamaño  2.2%   error 0.0140
k =  10:  tamaño  4.1%   error 0.0107
k =  20:  tamaño  8.0%   error 0.0078
k =  50:  tamaño 19.8%   error 0.0045
k = 100:  tamaño 39.3%   error 0.0023

• La curva baja rápido al principio y después se aplana: las primeras componentes se llevan lo grueso (las zonas grandes de luz y sombra) y las siguientes van agregando detalles cada vez más finos, como bordes y texturas.
• Con k = 50, guardando una quinta parte de los números, la foto ya se reconoce sin problema; con 100, la diferencia cuesta verla, pero el tamaño ya va por el 40 %.
• Un JPG suele comprimir más con la misma calidad, porque trabaja con bloques pequeños y aprovecha cómo ve el ojo. PCA no fue diseñado para fotos; lo interesante es que aprende las «piezas» de la imagen a partir de los propios datos.