Phase 01 - Lesson 11

Descomposición en Valores Singulares

This lesson includes a graded coding exercise that runs in your browser, unlocked with lifetime access.

La SVD es la navaja suiza del álgebra lineal. Toda matriz tiene una. Todo científico de datos necesita una.

Tipo: Build Lenguajes: Python, Julia Requisitos previos: Fase 1, Lecciones 01 (Intuición de Álgebra Lineal), 02 (Operaciones de Vectores y Matrices), 03 (Transformaciones de Matrices) Tiempo: ~120 minutos

Objetivos de Aprendizaje

• Implementar SVD via iteración de potencia y explicar el significado geométrico de U, Sigma y V^T
• Aplicar SVD truncada para compresión de imágenes y medir la razón de compresión vs error de reconstrucción
• Calcular la pseudoinversa de Moore-Penrose via SVD para resolver sistemas de mínimos cuadrados sobredeterminados
• Conectar la SVD con PCA, sistemas de recomendación (factores latentes) y Análisis Semántico Latente en NLP

El Problema

Tienes una matriz 1000x2000. Quizas sean calificaciones de usuario-pelicula. Quizas sea una tabla de frecuencia documento-termino. Quizas sean los valores de píxeles de una imagen. Necesitas comprimirla, eliminarle el ruido, encontrar estructura oculta en ella o resolver un sistema de mínimos cuadrados con ella. La descomposición en valores propios solo funciona en matrices cuadradas. Aún así, requiere que la matriz tenga un conjunto completo de vectores propios linealmente independientes.

La SVD funciona en cualquier matriz. Cualquier forma. Cualquier rango. Sin condiciones. Descompone la matriz en tres factores que revelan la geometría de lo que la matriz le hace al espacio. Es la factorización más general y más útil de toda el álgebra lineal.

El Concepto

Que hace la SVD geométricamente

Toda matriz, sin importar su forma, realiza tres operaciones en secuencia: rotar, escalar, rotar. La SVD hace explícita esta descomposición.

A = U * Sigma * V^T

      m x n     m x m    m x n    n x n
     (any)    (rotate)  (scale)  (rotate)

Dada cualquier matriz A, la SVD la factoriza en:

• V^T rota vectores en el espacio de entrada (n-dimensional)
• Sigma escala a lo largo de cada eje (estira o comprime)
• U rota el resultado hacia el espacio de salida (m-dimensional)

graph LR
    A["Input space (n-dim)\nData cloud\n(arbitrary orientation)"] -->|"V^T\n(rotate)"| B["Scaled space\nAligned with axes\nthen scaled by Sigma"]
    B -->|"U\n(rotate)"| C["Output space (m-dim)\nRotated to output\norientation"]

Piénselo así. Le entregas una matriz a la SVD. Te dice: "Esta matriz toma una esfera de entradas, primero la rota por V^T, luego la estira en un elipsoide por Sigma, luego rota el elipsoide por U." Los valores singulares son las longitudes de los ejes del elipsoide.

La descomposición completa

Para una matriz A con forma m x n:

A = U * Sigma * V^T

where:
  U     is m x m, orthogonal (U^T U = I)
  Sigma is m x n, diagonal (singular values on the diagonal)
  V     is n x n, orthogonal (V^T V = I)

The singular values sigma_1 >= sigma_2 >= ... >= sigma_r > 0
where r = rank(A)

Las columnas de U se llaman vectores singulares izquierdos. Las columnas de V se llaman vectores singulares derechos. Las entradas diagonales de Sigma se llaman valores singulares. Siempre son no negativas y convencionalmente se ordenan en orden decreciente.

Vectores singulares izquierdos, valores singulares, vectores singulares derechos

Cada componente de la SVD tiene un significado geométrico distinto.

Vectores singulares derechos (columnas de V): Forman una base ortonormal para el espacio de entrada (R^n). Son las direcciones en el espacio de entrada que la matriz mapea a direcciones ortogonales en el espacio de salida. Piénselos como el sistema de coordenadas natural para el dominio.

Valores singulares (diagonal de Sigma): Estos son los factores de escala. El i-ésimo valor singular te dice cuanto estira la matriz los vectores a lo largo del i-ésimo vector singular derecho. Un valor singular de cero significa que la matriz aplasta esa dirección por completo.

Vectores singulares izquierdos (columnas de U): Forman una base ortonormal para el espacio de salida (R^m). El i-ésimo vector singular izquierdo es la dirección en el espacio de salida donde aterriza el i-ésimo vector singular derecho (después del escalado).

La relación entre ellos:

A * v_i = sigma_i * u_i

The matrix A takes the i-th right singular vector v_i,
scales it by sigma_i, and maps it to the i-th left singular vector u_i.

Esto te da una imagen coordenada por coordenada de lo que hace cualquier matriz.

Forma de producto externo

La SVD puede escribirse como una suma de matrices de rango 1:

A = sigma_1 * u_1 * v_1^T + sigma_2 * u_2 * v_2^T + ... + sigma_r * u_r * v_r^T

Each term sigma_i * u_i * v_i^T is a rank-1 matrix (an outer product).
The full matrix is the sum of r such matrices, where r is the rank.

Esta forma es el fundamento de la aproximación de bajo rango. Cada término agrega una capa de estructura. El primer término capta el patron más importante. El segundo capta el siguiente más importante. Y así sucesivamente. Truncar esta suma te da la mejor aproximación posible en cualquier rango dado.

Rank-1 approx:    A_1 = sigma_1 * u_1 * v_1^T
                  (captures the dominant pattern)

Rank-2 approx:    A_2 = sigma_1 * u_1 * v_1^T + sigma_2 * u_2 * v_2^T
                  (captures the two most important patterns)

Rank-k approx:    A_k = sum of top k terms
                  (optimal by the Eckart-Young theorem)

Relación con la descomposición en valores propios

La SVD y la descomposición en valores propios están profundamente conectadas. Los valores y vectores singulares de A provienen directamente de los valores y vectores propios de A^T A y A A^T.

A^T A = V * Sigma^T * U^T * U * Sigma * V^T
      = V * Sigma^T * Sigma * V^T
      = V * D * V^T

where D = Sigma^T * Sigma is a diagonal matrix with sigma_i^2 on the diagonal.

So:
- The right singular vectors (V) are eigenvectors of A^T A
- The singular values squared (sigma_i^2) are eigenvalues of A^T A

Similarly:
A A^T = U * Sigma * V^T * V * Sigma^T * U^T
      = U * Sigma * Sigma^T * U^T

So:
- The left singular vectors (U) are eigenvectors of A A^T
- The eigenvalues of A A^T are also sigma_i^2

Esta conexión te dice tres cosas:

1. Los valores singulares siempre son reales y no negativos (son raíces cuadradas de valores propios de una matriz positiva semidefinida).
2. Podrías calcular la SVD via descomposición en valores propios de A^T A, pero esto eleva al cuadrado el número de condición y pierde precisión numérica. Los algoritmos de SVD dedicados evitan esto.
3. Cuando A es cuadrada y simétrica positiva semidefinida, la SVD y la descomposición en valores propios son lo mismo.

SVD truncada: aproximación de bajo rango

El teorema de Eckart-Young-Mirsky afirma que la mejor aproximación de rango k para A (tanto en la norma de Frobenius como en la espectral) se obtiene conservando solo los k mayores valores singulares y sus vectores correspondientes:

A_k = U_k * Sigma_k * V_k^T

where:
  U_k     is m x k  (first k columns of U)
  Sigma_k is k x k  (top-left k x k block of Sigma)
  V_k     is n x k  (first k columns of V)

Approximation error = sigma_{k+1}  (in spectral norm)
                    = sqrt(sigma_{k+1}^2 + ... + sigma_r^2)  (in Frobenius norm)

Esta no es solo "una buena" aproximación. Es demostrablemente la mejor aproximación posible de rango k. Ninguna otra matriz de rango k está más cerca de A.

Conserva los top 3: A_3 capta los tres mayores valores singulares. Error = valores restantes (sigma_4 hasta sigma_8).

Si los valores singulares decaen rápido, un k pequeño capta la mayor parte de la matriz. Si decaen lento, la matriz no tiene estructura de bajo rango.

Compresión de imágenes con SVD

Una imagen en escala de grises es una matriz de intensidades de píxeles. Una imagen 800x600 tiene 480.000 valores. La SVD permite aproximarla con muchos menos.

Original image: 800 x 600 = 480,000 values

SVD with rank k:
  U_k:      800 x k values
  Sigma_k:  k values
  V_k:      600 x k values
  Total:    k * (800 + 600 + 1) = k * 1401 values

  k=10:   14,010 values   (2.9% of original)
  k=50:   70,050 values  (14.6% of original)
  k=100: 140,100 values  (29.2% of original)

  The compression ratio improves as k gets smaller,
  but visual quality degrades.

El insight clave: las imágenes naturales tienen valores singulares que decaen rápidamente. Los primeros valores singulares captan la estructura amplia (formas, gradientes). Los posteriores captan detalles finos y ruido. Truncar en el rango 50 a menudo produce una imagen que se ve casi idéntica a la original mientras usa un 85% menos de almacenamiento.

SVD para sistemas de recomendación

El Premio Netflix hizo esto famoso. Tienes una matriz de calificaciones usuario-pelicula donde la mayoría de las entradas faltan.

             Movie1  Movie2  Movie3  Movie4  Movie5
  User1      [  5      ?       3       ?       1  ]
  User2      [  ?      4       ?       2       ?  ]
  User3      [  3      ?       5       ?       ?  ]
  User4      [  ?      ?       ?       4       3  ]

  ? = unknown rating

La idea: esta matriz de calificaciones tiene bajo rango. Los usuarios no tienen gustos completamente independientes. Existe un puñado de factores latentes (acción vs drama, antiguo vs nuevo, cerebral vs visceral) que explican la mayoría de las preferencias.

La SVD sobre la matriz de calificaciones (rellenada) la descompone en:

• U: perfiles de usuario en el espacio de factores latentes
• Sigma: importancia de cada factor latente
• V^T: perfiles de película en el espacio de factores latentes

La calificación predicha de un usuario para una película es el producto punto entre su perfil de usuario y el perfil de la película (ponderado por los valores singulares). La aproximación de bajo rango rellena las entradas faltantes.

En la practica, usas variantes como la SVD incremental de Simon Funk o ALS (alternating least squares) que manejan los datos faltantes directamente. Pero la idea central es la misma: descomposición en factores latentes via SVD.

SVD en NLP: Análisis Semántico Latente

El Análisis Semántico Latente (LSA), también llamado Indexación Semántica Latente (LSI), aplica SVD a una matriz termino-documento.

             Doc1   Doc2   Doc3   Doc4
  "cat"      [  3      0      1      0  ]
  "dog"      [  2      0      0      1  ]
  "fish"     [  0      4      1      0  ]
  "pet"      [  1      1      1      1  ]
  "ocean"    [  0      3      0      0  ]

After SVD with rank k=2:

  Each document becomes a point in 2D "concept space."
  Each term becomes a point in the same 2D space.
  Documents about similar topics cluster together.
  Terms with similar meanings cluster together.

  "cat" and "dog" end up near each other (land pets).
  "fish" and "ocean" end up near each other (water concepts).
  Doc1 and Doc3 cluster if they share similar topics.

El LSA fue uno de los primeros métodos exitosos para captar similitud semántica a partir de texto crudo. Funciona porque los términos sinónimos tienden a aparecer en documentos similares, así que la SVD los agrupa en las mismas dimensiones latentes. Los word embeddings modernos (Word2Vec, GloVe) pueden verse como descendientes de esta idea.

SVD para reducción de ruido

Los datos ruidosos tienen la señal concentrada en los valores singulares principales y el ruido repartido por todos los valores singulares. Truncar elimina el piso de ruido.

Valores singulares de señal limpia:

Valores singulares de señal ruidosa (el ruido se suma a todos):

graph TD
    A["All singular values"] --> B{"Clear gap?"}
    B -->|"Above gap"| C["Signal: keep these (top k)"]
    B -->|"Below gap"| D["Noise: discard these"]
    C --> E["Reconstruct with A_k to get denoised version"]

Esto se usa en procesamiento de señales, medición científica y limpieza de datos. Cada vez que tienes una matriz corrompida por ruido aditivo, la SVD truncada es una forma fundamentada de separar la señal del ruido.

Pseudoinversa via SVD

La pseudoinversa de Moore-Penrose A+ generaliza la inversión de matrices a matrices no cuadradas y singulares. La SVD hace trivial su cálculo.

If A = U * Sigma * V^T, then:

A+ = V * Sigma+ * U^T

where Sigma+ is formed by:
  1. Transpose Sigma (swap rows and columns)
  2. Replace each non-zero diagonal entry sigma_i with 1/sigma_i
  3. Leave zeros as zeros

For A (m x n):      A+ is (n x m)
For Sigma (m x n):  Sigma+ is (n x m)

La pseudoinversa resuelve problemas de mínimos cuadrados. Si Ax = b no tiene solución exacta (sistema sobredeterminado), entonces x = A+ b es la solución de mínimos cuadrados (minimiza ||Ax - b||).

Overdetermined system (more equations than unknowns):

  [1  1]         [3]
  [2  1] x   =   [5]       No exact solution exists.
  [3  1]         [6]

  x_ls = A+ b = V * Sigma+ * U^T * b

  This gives the x that minimizes the sum of squared residuals.
  Same result as the normal equations (A^T A)^(-1) A^T b,
  but numerically more stable.

Ventajas de estabilidad numérica

Calcular la descomposición en valores propios de A^T A eleva al cuadrado los valores singulares (los valores propios de A^T A son sigma_i^2). Esto eleva al cuadrado el número de condición, amplificando los errores numéricos.

Example:
  A has singular values [1000, 1, 0.001]
  Condition number of A: 1000 / 0.001 = 10^6

  A^T A has eigenvalues [10^6, 1, 10^{-6}]
  Condition number of A^T A: 10^6 / 10^{-6} = 10^{12}

  Computing SVD directly: works with condition number 10^6
  Computing via A^T A:     works with condition number 10^{12}
                           (6 extra digits of precision lost)

Los algoritmos modernos de SVD (bidiagonalización de Golub-Kahan) trabajan directamente sobre A, sin formar nunca A^T A. Por eso siempre deberias preferir np.linalg.svd(A) en lugar de np.linalg.eig(A.T @ A).

Conexión con PCA

PCA ES SVD sobre datos centrados. Esto no es una analogía. Es literalmente el mismo cálculo.

Given data matrix X (n_samples x n_features), centered (mean subtracted):

Covariance matrix: C = (1/(n-1)) * X^T X

PCA finds eigenvectors of C. But:

  X = U * Sigma * V^T    (SVD of X)

  X^T X = V * Sigma^2 * V^T

  C = (1/(n-1)) * V * Sigma^2 * V^T

So the principal components are exactly the right singular vectors V.
The explained variance for each component is sigma_i^2 / (n-1).

In sklearn, PCA is implemented using SVD, not eigendecomposition.
It is faster and more numerically stable.

Esto significa que todo lo que aprendiste sobre reducción de dimensionalidad en la Lección 10 es SVD por debajo. El PCA es la aplicación más común de la SVD en machine learning.

Construye

Paso 1: SVD desde cero usando iteración de potencia

La idea: para encontrar el mayor valor singular y sus vectores, usa iteración de potencia sobre A^T A (o A A^T). Luego deflaciona la matriz y repite para el siguiente valor singular.

import numpy as np

def power_iteration(M, num_iters=100):
    n = M.shape[1]
    v = np.random.randn(n)
    v = v / np.linalg.norm(v)

    for _ in range(num_iters):
        Mv = M @ v
        v = Mv / np.linalg.norm(Mv)

    eigenvalue = v @ M @ v
    return eigenvalue, v

def svd_from_scratch(A, k=None):
    m, n = A.shape
    if k is None:
        k = min(m, n)

    sigmas = []
    us = []
    vs = []

    A_residual = A.copy().astype(float)

    for _ in range(k):
        AtA = A_residual.T @ A_residual
        eigenvalue, v = power_iteration(AtA, num_iters=200)

        if eigenvalue < 1e-10:
            break

        sigma = np.sqrt(eigenvalue)
        u = A_residual @ v / sigma

        sigmas.append(sigma)
        us.append(u)
        vs.append(v)

        A_residual = A_residual - sigma * np.outer(u, v)

    U = np.column_stack(us) if us else np.empty((m, 0))
    S = np.array(sigmas)
    V = np.column_stack(vs) if vs else np.empty((n, 0))

    return U, S, V

Paso 2: Prueba y compara con NumPy

np.random.seed(42)
A = np.random.randn(5, 4)

U_ours, S_ours, V_ours = svd_from_scratch(A)
U_np, S_np, Vt_np = np.linalg.svd(A, full_matrices=False)

print("Our singular values:", np.round(S_ours, 4))
print("NumPy singular values:", np.round(S_np, 4))

A_reconstructed = U_ours @ np.diag(S_ours) @ V_ours.T
print(f"Reconstruction error: {np.linalg.norm(A - A_reconstructed):.8f}")

Paso 3: Demo de compresión de imagen

def compress_image_svd(image_matrix, k):
    U, S, Vt = np.linalg.svd(image_matrix, full_matrices=False)
    compressed = U[:, :k] @ np.diag(S[:k]) @ Vt[:k, :]
    return compressed

image = np.random.seed(42)
rows, cols = 200, 300
image = np.random.randn(rows, cols)

for k in [1, 5, 10, 20, 50]:
    compressed = compress_image_svd(image, k)
    error = np.linalg.norm(image - compressed) / np.linalg.norm(image)
    original_size = rows * cols
    compressed_size = k * (rows + cols + 1)
    ratio = compressed_size / original_size
    print(f"k={k:>3d}  error={error:.4f}  storage={ratio:.1%}")

Paso 4: Reducción de ruido

np.random.seed(42)
clean = np.outer(np.sin(np.linspace(0, 4*np.pi, 100)),
                 np.cos(np.linspace(0, 2*np.pi, 80)))
noise = 0.3 * np.random.randn(100, 80)
noisy = clean + noise

U, S, Vt = np.linalg.svd(noisy, full_matrices=False)
denoised = U[:, :5] @ np.diag(S[:5]) @ Vt[:5, :]

print(f"Noisy error:    {np.linalg.norm(noisy - clean):.4f}")
print(f"Denoised error: {np.linalg.norm(denoised - clean):.4f}")
print(f"Improvement:    {(1 - np.linalg.norm(denoised - clean) / np.linalg.norm(noisy - clean)):.1%}")

Paso 5: Pseudoinversa

A = np.array([[1, 1], [2, 1], [3, 1]], dtype=float)
b = np.array([3, 5, 6], dtype=float)

U, S, Vt = np.linalg.svd(A, full_matrices=False)
S_inv = np.diag(1.0 / S)
A_pinv = Vt.T @ S_inv @ U.T

x_svd = A_pinv @ b
x_lstsq = np.linalg.lstsq(A, b, rcond=None)[0]
x_pinv = np.linalg.pinv(A) @ b

print(f"SVD pseudoinverse solution:  {x_svd}")
print(f"np.linalg.lstsq solution:   {x_lstsq}")
print(f"np.linalg.pinv solution:    {x_pinv}")

Usalo

Las demos completas y funcionales están en code/svd.py. Ejecutalo para ver la SVD aplicada a compresión de imágenes, sistemas de recomendación, análisis semántico latente y reducción de ruido.

python svd.py

La versión en Julia en code/svd.jl demuestra los mismos conceptos usando la función nativa svd() de Julia y el paquete LinearAlgebra.

julia svd.jl

Entregalo

Esta lección produce:

• outputs/skill-svd.md - una skill para saber cuando y como aplicar SVD en proyectos reales

Ejercicios

1. Implementa la SVD completa desde cero sin usar iteración de potencia. En su lugar, calcula la descomposición en valores propios de A^T A para obtener V y los valores singulares, luego calcula U = A V Sigma^{-1}. Compara la precisión numérica con tu versión de iteración de potencia y con NumPy.

2. Carga una imagen real en escala de grises (o convierte una a escala de grises). Comprímela en los rangos 1, 5, 10, 25, 50, 100. Para cada rango, calcula la razón de compresión y el error relativo. Encuentra el rango donde la imagen se vuelve visualmente aceptable.

3. Construye un pequeño sistema de recomendación. Crea una matriz de calificaciones usuario-pelicula 10x8 con algunas entradas conocidas. Rellena las entradas faltantes con las medias de las filas. Calcula la SVD y reconstruye una aproximación de rango 3. Usa la matriz reconstruida para predecir las calificaciones faltantes. Verifica que las predicciones sean razonables.

4. Crea una matriz documento-termino 100x50 con 3 tópicos sintéticos. Cada tópico tiene 5 términos asociados. Agrega ruido. Aplica SVD y verifica que los 3 mayores valores singulares sean mucho mayores que el resto. Proyecta los documentos en el espacio latente 3D y comprueba que los documentos del mismo tópico se agrupen.

5. Genera una matriz limpia de bajo rango (rango 3, tamaño 50x40) y agrega ruido gaussiano en distintos niveles (sigma = 0.1, 0.5, 1.0, 2.0). Para cada nivel de ruido, encuentra el rango de truncamiento óptimo barriendo k de 1 a 40 y midiendo el error de reconstrucción respecto a la matriz limpia. Grafica como cambia el k óptimo con el nivel de ruido.

Términos Clave

Lectura Adicional

• Gilbert Strang: Linear Álgebra and Its Applications, Chapter 7 - tratamiento completo de la SVD con aplicaciones
• 3Blue1Brown: But what is the SVD? - intuición geométrica para la SVD
• We Recommend a Singular Value Decomposition - visión general accesible de la American Mathematical Society
• Netflix Prize and Matrix Factorization - entrada de blog original de Simon Funk sobre SVD para recomendaciones
• Latent Semantic Analysis - la aplicación original de NLP de la SVD
• Numerical Linear Álgebra by Trefethen and Bau - el estándar de oro para entender los algoritmos de SVD y sus propiedades numéricas