En el análisis de datos, a menudo nos enfrentamos a la tarea de representar un fenómeno mediante un modelo matemático. La diferencia entre el fenómeno real y nuestro modelo se denomina error de representación.
Supongamos que tenemos \(m\) puntos de datos \((x_i, y_i)\). Si buscamos un polinomio de grado \(m-1\) que pase exactamente por todos los puntos, estamos realizando una interpolación (Armentano 2026 clase 17). Esto requiere resolver un sistema lineal \(Vc = y\), donde \(V\) es la matriz de Vandermonde:
\[V = \begin{pmatrix} 1 & x_1 & \dots & x_1^{m-1} \\ \vdots & \vdots & \ddots & \vdots \\ 1 & x_m & \dots & x_m^{m-1} \end{pmatrix}\]
Si los valores de \(x_i\) son distintos, \(V\) es invertible y existe una solución única. Sin embargo, la interpolación exacta tiene riesgos: si los datos tienen ruido, el polinomio puede presentar oscilaciones salvajes (como el fenómeno de Runge), lo que invalida su capacidad predictiva.
Cuando el número de datos \(m\) es mucho mayor que el número de parámetros \(n\) del modelo (sistema sobredeterminado), generalmente no existe una solución exacta para \(Ax = b\) (Strang 2018 Lec. 9). En lugar de forzar un ajuste exacto en pocos puntos, buscamos una solución aproximada \(\hat{x}\) que sea “lo más cercana posible” a todos los datos simultáneamente.
El criterio estándar, propuesto originalmente por Gauss, consiste en minimizar la norma \(\ell_2\) del residuo:
\[\min \left\| Ax - b \right\| _2^2 = \min \sum_{i=1}^m (ax_i - b_i)^2\]
Este cambio de paradigma transforma un problema sin solución (en el sentido estricto) en un problema de optimización con una solución robusta y única, sentando las bases de la regresión lineal y el aprendizaje automático moderno.
17.1 Las Ecuaciones Normales
Cuando un sistema \(Ax = b\) es sobredeterminado (\(m > n\)), el vector \(b\) generalmente no pertenece al espacio columna de \(A\), por lo que no existe una solución exacta. El objetivo de los mínimos cuadrados es encontrar un vector \(\hat{x}\) que minimice el error cuadrático total:
\[\hat{x} = \mathop{\mathrm{argmin}}_{x \in C(A)} \left\| Ax - b \right\| _2^2\]
Geométricamente, la distancia más corta desde el vector \(b\) hasta el subespacio \(C(A)\) se encuentra mediante una proyección ortogonal. Para que el error \(e = b - A\hat{x}\) sea mínimo, este debe ser perpendicular a todas las columnas de \(A\). En otras palabras, buscamos el vector \(\hat{x}\) en \(C(A)\) más cercano a \(b\) en términos de distancia euclídea.
Para que el error sea mínimo, el vector de error (\(b - A\hat{x}\)) debe ser perpendicular a todas las columnas de \(A\). En lenguaje matricial, esto significa que el producto de \( A^\top \) por el error \(e\) debe ser cero: \( A^\top \cdot e = 0\). Entonces \( A^\top \cdot (b - A\hat{x}) = 0\). Distribuyendo los términos, obtenemos las Ecuaciones Normales:
\[A^\top A \hat{x} = A^\top b \tag{17.1}\]
La matriz \(A^\top A\) es fundamental en este proceso. Es una matriz cuadrada (\(n \times n\)), simétrica y semidefinida positiva. Si las columnas de \(A\) son linealmente independientes, entonces la matriz \(A^\top A\) es invertible y la solución \(\hat{x}\) es única:
\[\hat{x} = (A^\top A)^{-1} A^\top b\]
Si las columnas son dependientes, existen infinitas soluciones que minimizan el error, y suele seleccionarse la de norma mínima mediante la pseudoinversa.
Este método transforma un problema de aproximación en un sistema lineal cuadrado que puede resolverse mediante algoritmos estándar como la eliminación gaussiana.
Resolución como Problema de Optimización
El problema de mínimos cuadrados, tambien puede pensarse como un problema de minimización donde queremos encontrar \(\hat{x}\) que minimice la función de pérdida cuadrática:
\[f(x) = \left\| Ax - b \right\| ^2 = (Ax - b)^\top (Ax - b) = x^T A^T A x - 2 b^T A x + b^T b\] Calculamos la derivada (gradiente) de \(f(x)\) respecto a \(x\) e gualamos a cero:
\[\nabla f(x) = 2 A^\top A x - 2 A^\top b\]
\[\nabla f(x) = 0 \quad \Rightarrow \quad A^\top A x = A^\top b\]
Llegamos nuevamente a la ecuación normal de mínimos cuadrados. La solución \(\hat{x}\) es el punto crítico (mínimo global, ya que \(f\) es convexa).
17.2 Proyeccion Ortogonal
La solución del problema de mínimos cuadrados tiene una interpretación geométrica elegante: estamos buscando el punto en el subespacio \(C(A)\) que se encuentra a la distancia mínima del vector de datos \(b\) (Strang 2018 Lec. 9; Strang 2019, II.2; Armentano 2026 Clase 17).
Cuando un sistema \(Ax = b\) es incompatible, es decir, \(b \notin C(A)\), la solución numérica óptima no busca resolver la igualdad, sino minimizar el tamaño del error: \(\|b - Ax\|\). En otras palabras, buscamos el vector \(\hat{x}\) en \(C(A)\) más cercano a \(b\) en términos de distancia euclídea.
Si las columnas de \(A\) se transforman en una base ortonormal \(Q\), la “mejor aproximación” de \(b\) en el subespacio es simplemente su proyección ortogonal \(Pb = QQ^\top b\).
Pero estamos interesados en encontrar una solucion para cualquier \(A\).
Dado que generalmente el sistema es incompatible ( \(b \notin C(A)\) ), no existe un vector \(x\) tal que \(Ax = b\). La mejor aproximación es el vector \(A\hat{x}\), que es la proyección ortogonal de \(b\) sobre el hiperplano formado por las columnas de \(A\).
El vector error (o residuo) \(e = b - A\hat{x}\) representa la parte de los datos que el modelo no puede explicar. Para que la distancia \(\|e\|\) sea mínima, el vector de error debe ser perpendicular a todas las columnas de \(A\). Esto significa que \(e\) pertenece a \( \mathcal{N}( A^\top ) \), espacio nulo de \( A^\top \), lo que nos lleva nuevamente a las ecuaciones normales.
Operador de Proyeccion Ortogonal
Podemos expresar la proyección directamente a partir de \(b\) mediante una matriz \(P\) tal que \(Pb = A\hat{x}\). Al sustituir la solución \(\hat{x}\) de las ecuaciones normales (Ecuación 17.1), obtenemos \(Pb = A(A^\top A)^{-1}A^\top b\)
Lo que hace visible la estructura de la matriz de proyección: \[P = A(A^\top A)^{-1} A^\top \tag{17.2}\]
- Observaciones
-
La matriz de proyección ortogonal \(P\) es:
- Simétrica (\( P^\top =P\))
- Idempotente (\(P^2=P\)).
- Si A es ortogonal, entonces \(P=A A^\top \).
- Simetria
- \( P^\top = \left( A ( A^\top A)^{-1} A^\top \right) ^\top = A \left( ( A^\top A)^\top \right) ^{-1} A^\top = A ( A^\top A )^{-1} A^\top = P\)
- Idempotencia
- \(P^2 = ( A ( A^\top A)^{-1} A^\top )( A ( A^\top A)^{-1} A^\top ) = \\ A ( A^\top A)^{-1} ( A^\top A) ( A^\top A)^{-1} A^\top = A ( A^\top A)^{-1} A^\top = P\)
- Caso Ortogonal
- Si las columnas de \(A\) son ortonormales (\(A = Q\)), la matriz de proyección se simplifica drásticamente a \(P = A A^\top \), ya que \(P =Q ( Q^\top Q)^{-1} Q^\top = Q \mathbb{I} Q^\top = Q Q^\top \). (Armentano 2026 Clase 5; Strang 2018 Lec. 3)
17.3 Solución vía SVD - La Pseudoinversa \(A^+\)
La Descomposición en Valores Singulares proporciona la herramienta definitiva para resolver el problema de mínimos cuadrados, especialmente en casos donde la matriz \(A\) es singular o tiene columnas dependientes. En estos escenarios, las ecuaciones normales fallan porque \(A^\top A\) no es invertible. La solución es la pseudoinversa \(A^+\).
Si conocemos la SVD de una matriz, \(A = U \Sigma V^\top\), su pseudoinversa se construye invirtiendo los componentes de la descomposición:
\[A^+ = V \Sigma^+ U^\top\]
Donde \(\Sigma^+\) es una matriz de dimensiones \(n \times m\) que contiene los recíprocos de los valores singulares no nulos (\(1/\sigma_i\)) en su diagonal, y ceros en todas las demás entradas.
La importancia de la pseudoinversa radica en que siempre proporciona una solución válida al problema \(\min \|Ax - b\|_2^2\), denotada como \(x^+ = A^+ b\):
Si las columnas de \(A\) son independientes, \(A^+\) coincide con el operador de las ecuaciones normales: \(A^+ = (A^\top A)^{-1} A^\top\). En este caso, \(x^+\) es la única solución de mínimos cuadrados.
Si \(A\) tiene columnas dependientes, existen infinitos vectores \(x\) que minimizan el error. La pseudoinversa selecciona automáticamente el vector de norma mínima (\(\min \|x\|_2\)).
Geométricamente, \(A^+\) actúa como un puente entre los cuatro subespacios fundamentales:
- Lleva el espacio columna \(C(A)\) de vuelta al espacio fila \(C( A^\top )\) de forma perfecta.
- Cualquier componente del vector \(b\) que se encuentre en \( \mathcal{N}( A^\top ) \), espacio nulo de \(A^\top\), es enviado directamente a cero por \(A^+\), eliminando el ruido que no puede ser explicado por el modelo.
17.4 Ortogonalización: El método de Gram-Schmidt
Una tercera vía para resolver el problema de mínimos cuadrados consiste en transformar las columnas de la matriz \(A\) en una base ortonormal. Este proceso, conocido como Gram-Schmidt, permite factorizar la matriz como \(A = QR\). (Strang 2018 Lec. 11; Armentano 2026 Clase 17) 1
\(Q\) Es una matriz de Stiefel \(m \times n\) (columnas ortonormales entre sí, es decir, \( Q^\top Q = \mathbb{I} \)). Geométricamente, \(Q\) abarca el mismo espacio columna que \(A\), pero sin redundancias ni problemas de escala.
\(R\) Es una matriz triangular superior que contiene los coeficientes de las combinaciones lineales que vinculan a \(A\) con \(Q\). Al ser triangular, facilita enormemente la resolución de sistemas mediante sustitución hacia atrás.
Sustituyendo \(A = QR\) en las ecuaciones normales (\( A^\top A \hat{x} = A^\top b\)):
\[ (QR)^\top (QR) \hat{x} = (QR)^\top b\] \[ R^\top ( Q^\top Q) R \hat{x} = R^\top Q^\top b\]
Como \( Q^\top Q = \mathbb{I} \) y asumiendo que \(R\) es invertible (columnas de \(A\) independientes), podemos simplificar el sistema a: \(R\hat{x} = Q^\top b\) y por lo tanto: \[\hat{x} = R^{-1} Q^\top b\]
A diferencia de las ecuaciones normales, que requieren calcular \(A^\top A\) (lo cual puede duplicar el error de redondeo si la matriz está mal condicionada), el método \(QR\) es numéricamente más estable. Al operar directamente sobre \(Q\) y \(R\), se evita la “pérdida de precisión” asociada al cuadrado de los valores singulares, convirtiéndolo en el estándar de la industria para sistemas sobredeterminados de tamaño moderado.↩︎