Logo de Pablo MadariagaPablo MadariagaData ScientistAgendar reunión

Capítulo 6

SVD y Aproximación de Matrices

6.1 Descomposición en Valores Singulares

Toda matriz A∈Rm×nA \in \mathbb{R}^{m \times n} se puede factorizar como:

A=UΣVTA = U \Sigma V^T

UU

Matriz ortogonal m×m. Columnas = vectores singulares izquierdos (eigenvectores de AATAA^T).

Σ\Sigma

Diagonal m×n con valores singulares σ1≥σ2≥⋯≥0\sigma_1 \ge \sigma_2 \ge \cdots \ge 0.

VV

Matriz ortogonal n×n. Columnas = vectores singulares derechos (eigenvectores de ATAA^T A).

6.2 SVD en Python

python
import numpy as np

A = np.array([[1, 2, 3],
              [4, 5, 6],
              [7, 8, 9]], dtype=float)

U, S, Vt = np.linalg.svd(A)

print("U:
", U)
print("Valores singulares:", S)
print("V^T:
", Vt)

# Reconstruir
Sigma = np.zeros((3, 3))
np.fill_diagonal(Sigma, S)
A_rec = U @ Sigma @ Vt
print("A reconstruida:
", A_rec)

# Relación con eigendecomposición
# A^T A = V Σ^2 V^T
AtA = A.T @ A
vals, vecs = np.linalg.eigh(AtA)
print("Eigenvalores de A^T A:", np.sort(vals)[::-1])
print("Cuadrados de valores singulares:", S**2)

6.3 Aproximación de Rango k

La mejor aproximación de rango k a una matriz A (en norma de Frobenius) se obtiene tomando los k valores singulares más grandes:

Ak=UkΣkVkTA_k = U_k \Sigma_k V_k^T

Teorema de Eckart-Young: ∥A−Ak∥F=∑i=k+1rσi2\|A - A_k\|_F = \sqrt{\sum_{i=k+1}^r \sigma_i^2}.

python
import numpy as np
from sklearn.datasets import load_sample_image
import matplotlib.pyplot as plt

# Cargar imagen
img = load_sample_image("flower.jpg")
img_gray = img.mean(axis=2)  # a escala de grises

# SVD de la imagen
U, S, Vt = np.linalg.svd(img_gray, full_matrices=False)

# Aproximación con k componentes
k = 20
img_approx = U[:, :k] @ np.diag(S[:k]) @ Vt[:k, :]

print(f"Original: {img_gray.shape}")
print(f"k={k}: compresión {img_gray.size / (k * (img_gray.shape[0] + img_gray.shape[1])):.1f}x")
print(f"Error relativo: {np.linalg.norm(img_gray - img_approx) / np.linalg.norm(img_gray):.4f}")

# Varianza explicada
var_explicada = np.cumsum(S**2) / np.sum(S**2)
print(f"Varianza explicada con k={k}: {var_explicada[k-1]:.2%}")

6.4 SVD y PCA

PCA se puede realizar directamente con SVD sin calcular la matriz de covarianza:

  1. Centrar los datos X (restar media por columna)
  2. SVD: X=UΣVTX = U \Sigma V^T
  3. Componentes principales: T=UΣT = U \Sigma (o XVXV)
  4. Las columnas de V son los loadings (direcciones principales)
python
import numpy as np
from sklearn.decomposition import PCA

# Datos sintéticos
np.random.seed(42)
X = np.random.randn(100, 5)

# PCA vía SVD manual
X_centered = X - X.mean(axis=0)
U, S, Vt = np.linalg.svd(X_centered, full_matrices=False)
scores = U * S  # = X_centered @ Vt.T

# vs sklearn
pca = PCA()
X_pca = pca.fit_transform(X)

print("Componentes (SVD):
", Vt.T[:2])
print("Componentes (sklearn):
", pca.components_[:2])
print("Coinciden?", np.allclose(np.abs(Vt.T), np.abs(pca.components_.T)))

6.5 Usos del SVD en ML

Reducción de dimensionalidad

PCA, compresión de imágenes, reducción de ruido.

Pseudoinversa

A+=VΣ+UTA^+ = V \Sigma^+ U^T

Base de la regresión lineal por mínimos cuadrados.

Sistemas de recomendación

Factorización de matrices en el corazón de SVD en recomendación.

python
# Pseudoinversa vía SVD
A = np.array([[1, 2], [3, 4], [5, 6]])
A_pinv = np.linalg.pinv(A)
print("Pseudoinversa:
", A_pinv)
print("A^+ @ A:
", A_pinv @ A)  # ≈ I (en el sentido de mínimos cuadrados)

# Sistema sobredeterminado: Ax = b
b = np.array([1, 2, 3])
x = A_pinv @ b
print("Solución mínimos cuadrados:", x)