9.5: Integración Montecarlo

El esquema final de integración numérica que discutiremos es la integración de Monte Carlo, y conceptualmente es completamente diferente de los esquemas anteriores. En lugar de asignar un conjunto de puntos de discretización (ya sea explícitamente, como en las reglas del punto medio/trapeo/Simpson, o a través de un procedimiento de optimización de máquina, como en el método de cuadratura adaptativa), este método muestrear aleatoriamente los puntos en el dominio de integración. Si los puntos de muestreo son independientes y hay un número suficientemente grande de ellos, la integral se puede estimar tomando un promedio ponderado del integrando sobre los puntos de muestreo.

Para ser precisos, considere una integral 1D sobre un dominio\(x \in [a,b]\). Deje que cada punto de muestreo se dibuje independientemente de una distribución\(p(x)\). Esto significa que la probabilidad de dibujar muestra\(x_{n}\) en el rango\(x_n \in [x, x + dx]\) es\(p(x) dx\). La distribución se normaliza, de manera que

\[\int_a^b p(x) \; dx = 1.\]

Tomemos\(N\) muestras y evaluemos el integrando en esos puntos: esto nos da un conjunto de números\(\{f(x_n)\}\). Luego calculamos la cantidad

\[\mathcal{I}^{\mathrm{mc}} = \frac{1}{N} \sum_{n=0}^{N-1} \frac{f(x_n)}{p(x_n)}.\]

A diferencia de los estimadores que hemos estudiado previamente,\(\mathcal{I}^{\mathrm{mc}}\) es un número aleatorio (porque las variables subyacentes\(\{x_n\}\) son todas aleatorias). Fundamentalmente, su valor promedio es igual a la integral deseada:

\[\begin{align} \left\langle\mathcal{I}^{\mathrm{mc}}\right\rangle &= \frac{1}{N} \sum_{n=0}^{N-1} \left\langle \frac{f(x_n)}{p(x_n)}\right\rangle \\ & = \left\langle \frac{f(x_n)}{p(x_n)}\right\rangle \qquad\mathrm{for}\;\mathrm{each}\; n \\ & = \int_a^b p(x) \, \left[\frac{f(x)}{p(x)}\right]\, dx \\ & = \int_a^b f(x)\, dx \end{align}\]

Para integrales de baja dimensión, normalmente no hay razón para utilizar el método de integración Monte Carlo. Requiere un número mucho mayor de muestras para alcanzar un nivel de precisión numérica comparable a los otros métodos de integración numérica. (Para integrales 1D, la integración de Monte Carlo generalmente requiere millones de muestras, mientras que la regla de Simpson solo requiere cientos o miles de puntos de discretización). Sin embargo, la integración de Monte Carlo supera a los esquemas de integración basados en la discretización cuando la dimensionalidad de la integración se vuelve extremadamente grande. Tales integrales ocurren, por ejemplo, en cálculos mecánicos cuánticos que involucran sistemas de muchos cuerpos, donde la dimensionalidad del espacio Hilbert escala exponencialmente con el número de partículas.