Algoritmo de la Pseudoinversa de Moore‑Penrose

La pseudoinversa de Moore‑Penrose es una herramienta fundamental en álgebra lineal numérica. Permite resolver sistemas lineales que no son cuadráticos o que son mal condicionados, y se utiliza ampliamente en regresión, control óptimo y aprendizaje automático.

1. Fundamentos Matemáticos

Sea A ∈ ℝ^{m×n}. La pseudoinversa A⁺ es la única matriz que satisface las cuatro condiciones de Moore‑Penrose:

• 1. AA⁺A = A
• 2. A⁺AA⁺ = A⁺
• 3. (AA⁺)ᵀ = AA⁺
• 4. (A⁺A)ᵀ = A⁺A

Cuando A tiene rango completo, A⁺ coincide con la inversa clásica (A⁻¹) si m=n. En caso contrario, A⁺ brinda la solución de mínima norma a Ax = b.

2. Algoritmo paso a paso (SVD)

El método más estable numéricamente para obtener A⁺ es mediante la descomposición en valores singulares (SVD):

 A = U Σ Vᵀ  →  A⁺ = V Σ⁺ Uᵀ 

1. Computar la SVD: A = U Σ Vᵀ, donde U ∈ ℝ^{m×m}, V ∈ ℝ^{n×n} ortonormales y Σ diagonal con los valores singulares σ_i ≥ 0.
2. Invertir los valores singulares no nulos: crear Σ⁺ poniendo 1/σ_i en la posición diagonal correspondiente; los valores singulares σ_i = 0 permanecen en cero.
3. Formar la pseudoinversa: A⁺ = V Σ⁺ Uᵀ.

Este proceso es robusto frente a matrices mal condicionadas porque los valores singulares pequeños pueden filtrarse mediante un umbral (tol).

3. Implementación en Python

Python ofrece dos librerías de referencia: numpy.linalg y scipy.linalg. Ambas exponen la función pinv que implementa el algoritmo SVD bajo el capó.

import numpy as np
# Matriz de ejemplo (3x2)
A = np.array([[1, 2],
              [3, 4],
              [5, 6]], dtype=float)
# Pseudoinversa
A_pinv = np.linalg.pinv(A)
print("A⁺ =\n", A_pinv)
# Verificación de la condición 1 (AA⁺A ≈ A)
print("AA⁺A =\n", A @ A_pinv @ A)

import numpy as np
from scipy.linalg import svd, diagsvd
A = np.array([[1, 0, 2],
              [0, 1, 0]], dtype=float)  # 2x3
U, s, VT = svd(A, full_matrices=False)
# Umbral: 1e-10 * max(m, n) * max(s)
 tol = 1e-10 * max(A.shape) * np.max(s)
# Inversión segura
s_inv = np.array([1/si if si > tol else 0 for si in s])
Sigma_pinv = diagsvd(s_inv, A.shape[1], A.shape[0])
A_pinv = VT.T @ Sigma_pinv @ U.T
print("A⁺ con control de tolerancia =\n", A_pinv)

4. Casos de Uso Reales

A continuación, tres escenarios donde la pseudoinversa es la herramienta de elección.

1. Regresión lineal de mínimos cuadrados: Cuando el número de observaciones m supera al número de variables n, la solución β = A⁺y minimiza ‖Aβ - y‖₂.
2. Control de sistemas bajo‑determinado: En robótica, la cinemática inversa de un manipulador con más grados de libertad que restricciones se resuelve con θ = J⁺ v, donde J es la Jacobiana.
3. Reducción de dimensionalidad (PCA simplificado): La proyección de datos en el subespacio de rango r se puede obtener mediante X⁺ = V_r Σ_r⁻¹ U_rᵀ, facilitando la reconstrucción de versiones comprimidas.

5. Comparativa con otras técnicas

6. Buenas Prácticas y Troubleshooting

• Control de tolerancia: Ajuste rcond (en np.linalg.pinv) o defina su propio tol para evitar la amplificación de ruido cuando σ_i es muy pequeño.
• Escalado de columnas: Normalizar la matriz antes de la SVD mejora la condición numérica.
• Memoria: Para matrices muy grandes (>10⁴ × 10⁴), considere algoritmos iterativos (p.ej., scipy.sparse.linalg.lsmr) que evitan la SVD completa.
• Validación: Verifique siempre al menos una de las condiciones de Moore‑Penrose (AA⁺A≈A) para confirmar la precisión.
• Seguridad: Cuando la pseudoinversa se usa en APIs públicas, limite la dimensión de entrada para evitar ataques de denegación de servicio (DoS) por matrices gigantes.

7. Rendimiento y Escalabilidad

En entornos de producción (por ejemplo, micro‑servicios de inferencia), la generación única de A⁺ se puede cachear. Para flujos de datos en tiempo real:

• Utilice numpy.linalg.pinv para < 1 000 × 1 000 (≈ 0.1 s).
• Para > 5 000 × 5 000, prefiera scipy.sparse.linalg.lsmr o torch.linalg.pinv con GPU.

8. Conclusiones

La pseudoinversa de Moore‑Penrose es una herramienta esencial que combina rigor matemático con una implementación práctica y estable en Python. Su capacidad para resolver sistemas sobredeterminados y subdeterminados la hace indispensable en ciencia de datos, control y modelado.