Introducción

La matriz Hessiana es una herramienta esencial en cálculo multivariable, optimización y aprendizaje automático. Representa todas las segundas derivadas parciales de una función escalar \(f:\mathbb{R}^n\to\mathbb{R}\) y permite analizar la curvatura del espacio de parámetros.

En este artículo encontrarás:

• Definición matemática y propiedades clave.
• Métodos para calcularla: derivación simbólica, diferencias finitas y autodiferenciación.
• Ejemplos de código en Python con NumPy, SciPy, autograd, JAX y PyTorch.
• Comparativas de rendimiento, buenas prácticas y solución de problemas comunes.

Definición matemática

Para una función \(f(\mathbf{x})\) con \(\mathbf{x}= (x_1,\dots,x_n)\), la Hessiana \(H(\mathbf{x})\) es la matriz \(n\times n\) cuyas entradas son:

Propiedades relevantes:

• Es simétrica si las segundas derivadas cruzadas son continuas (teorema de Schwarz).
• Los eigenvalores de \(H\) indican la curvatura: positivos → mínimo local, negativos → máximo local, mixtos → punto de silla.

Métodos para calcular la Hessiana

Ejemplos prácticos en Python

A continuación se presentan fragmentos de código para cada método mencionado.

1️⃣ Derivación simbólica con sympy

import sympy as sp
# Definir variables simbólicas
x, y = sp.symbols('x y')
# Función ejemplo
f = sp.exp(x * y) + sp.sin(x) * sp.cos(y)
# Calcular Hessiana simbólica
H = sp.hessian(f, (x, y))
print(H)

Ventaja: expresión exacta; Desventaja: no escalable a >10 variables.

2️⃣ Diferencias finitas con SciPy

import numpy as np
from scipy.optimize import approx_fprime
def f(v):
    x, y = v
    return np.exp(x * y) + np.sin(x) * np.cos(y)
epsilon = np.sqrt(np.finfo(float).eps)
# Función que devuelve gradiente
grad = lambda v: approx_fprime(v, f, epsilon)
# Hessiana mediante diferencias finitas de segundo orden
def hessian(v):
    n = len(v)
    H = np.zeros((n, n))
    for i in range(n):
        ei = np.zeros(n); ei[i] = epsilon
        H[i, :] = (grad(v + ei) - grad(v - ei)) / (2 * epsilon)
    return H
print(hessian(np.array([0.5, -0.3])))

Consejo: elegir epsilon adecuado para evitar cancelación numérica.

3️⃣ Autodiferenciación con autograd

import autograd.numpy as anp
from autograd import grad, hessian
def f(v):
    x, y = v[0], v[1]
    return anp.exp(x * y) + anp.sin(x) * anp.cos(y)
H = hessian(f)
print(H(anp.array([0.5, -0.3])))

Rápido y sin dependencias pesadas; limitado a CPU.

4️⃣ Autodiferenciación con JAX (GPU/TPU)

import jax.numpy as jnp
from jax import grad, jacfwd, jacrev, hessian
def f(v):
    x, y = v[0], v[1]
    return jnp.exp(x * y) + jnp.sin(x) * jnp.cos(y)
H = hessian(f)
print(H(jnp.array([0.5, -0.3])))

Ventaja: JIT compilation, ejecución en GPU/TPU, soporte para vectores de gran dimensión.

5️⃣ Autodiferenciación con PyTorch

import torch
def f(v):
    x, y = v[0], v[1]
    return torch.exp(x * y) + torch.sin(x) * torch.cos(y)
# Tensor con gradiente habilitado
v = torch.tensor([0.5, -0.3], requires_grad=True)
# Gradiente de primer orden
grad_f = torch.autograd.grad(f(v), v, create_graph=True)[0]
# Hessiana (segunda derivada) mediante grad de grad
H = []
for g in grad_f:
    row = torch.autograd.grad(g, v, retain_graph=True)[0]
    H.append(row)
H = torch.stack(H)
print(H)

Ideal para modelos de deep learning donde la Hessiana se necesita para algoritmos de segunda orden como Newton’s method o Trust‑Region.

Comparativa de rendimiento y escalabilidad

Buenas prácticas y solución de problemas

✅ Validación de resultados

• Comparar la Hessiana obtenida con autodiferenciación contra una versión de diferencias finitas para pequeños vectores.
• Verificar la simetría numérica: \(\|H - H^T\|_F < 1e-8\).

⚠️ Problemas comunes

• Desbordamiento de punto flotante: funciones exponenciales pueden producir valores muy grandes; usar torch.float64 o dtype=jnp.float64 cuando sea necesario.
• Gradientes nulos: ocurre si la función es constante respecto a alguna variable; asegurarse de que todas las variables influyan en la salida.
• Inestabilidad en diferencias finitas: elegir un epsilon demasiado pequeño genera cancelación; usar la regla de thumb \(\epsilon = \sqrt[3]{\text{machine\_epsilon}}\).

🔧 Optimización

• Activar jit en JAX (jax.jit) para compilar la función una sola vez.
• En PyTorch, usar torch.set_grad_enabled(False) cuando solo se necesiten productos Hessiano‑vector.
• Para grandes dimensiones, emplear scipy.sparse.linalg.eigsh sobre la Hessiana esparcida.

Casos de uso en el mundo real

• Optimización de hiperparámetros: Métodos de segunda orden (Newton, BFGS) requieren la Hessiana o su aproximación.
• Modelos de regresión logística: La función de pérdida log‑likelihood tiene Hessiana que define la covarianza del estimador.
• Redes neuronales: Algoritmos de confianza (Laplace Approximation) usan la Hessiana del log‑posterior.
• Visión por computadora: En ajuste de pose y bundle adjustment, la Hessiana estructurada permite resolver sistemas lineales esparcidos.

Conclusión

La matriz Hessiana es una herramienta poderosa que, combinada con los modernos frameworks de autodiferenciación, permite implementar algoritmos de segunda orden de forma eficiente y escalable. Elegir el método adecuado depende del tamaño del problema, del hardware disponible y de la precisión requerida.