Aceleración de los cálculos: cuándo utilizar Cuasi-Carlo y Monte-Carlo estándar en la fijación de precios

Question

Estoy familiarizado con la teoría de las técnicas de Monte-Carlo en la integración numérica, y recientemente he comenzado mis experimentos con estos métodos aplicados a la fijación de precios de los derivados. Estoy utilizando el libro de Glasserman como referencia.

He empezado por calcular simplemente el precio de la llamada de vainilla de los cajeros automáticos. Mi configuración es muy básica: tengo 100 pasos de tiempo, utilizo el modelo Black-Scholes (para poder comprobar la validez de los resultados) y los esquemas de integración numérica Euler-Maryuama, Runge Kutta y Milstein. Parece que necesito unos $10^6$ muestras para obtener resultados estables, y en Python eso lleva unos 6 minutos. Tengo una máquina bastante rápida, y esta velocidad me parece realmente baja. Inicialmente pensé que la razón está en el muestreador aleatorio. Sin embargo, por lo que he comprobado, el muestreo de variables aleatorias normales sí lleva algo de tiempo, pero $80\%$ del tiempo de cálculo proviene de las operaciones aritméticas (suma y multiplicación), incluso en el caso del esquema básico de Euler-Maryuama. Entiendo que en el caso del Movimiento Browniano Geométrico hay técnicas de simulación mucho más inteligentes disponibles, pero me gustaría mantenerlo genérico ya que en última instancia me gustaría trabajar con difusiones más generales.

Quizá sea sólo un problema de Python, pero ese problema me hizo pensar en acelerar mis simulaciones de Monte-Carlo. En el libro de Glasserman hay dos tipos de procedimientos de aceleración:

1. Reducción de la varianza, incluidas las variantes de control, el muestreo antitético y el muestreo de importancia. Ya he utilizado este último, y mi mejor opción es esa.

2. Quasi Monte Carlo (secuencias de baja discrepancia). Nunca lo he utilizado, y ahora me parece que sólo puede aplicarse para calcular integrales en la forma habitual, digamos que cuando las densidades se dan explícitamente sólo se desconoce la fórmula analítica de la integral.

Hasta ahora, mi objetivo es ser capaz de fijar el precio de las opciones exóticas con Monte-Carlo. Digamos que el modelo es unidimensional, y el pago exótico es extremadamente dependiente de la trayectoria. Tal vez más adelante necesite valorar también una versión americana de dicha opción. Mis preguntas son las siguientes:

1. ¿A cuál de los métodos de reducción de la varianza debo prestar especial atención?

2. ¿Es posible aplicar QMC en un problema de este tipo, o realmente necesito tener una expresión integral (sobre un dominio de dimensión finita) para el precio de la opción?

---

Para aclarar mi pregunta sobre QMC. Veo tres formas de aplicar (Q)MC en la fijación de precios.

• para muestrear incrementos de componentes de movimiento browniano/salto sobre la línea real
• para muestrear trayectorias enteras del espacio de trayectorias
• si el valor de la opción se han encontrado tener forma $$ \tag{1} \int_D p(S_1,\dots,S_n)f(S_1,\dots,S_n)\;\mathrm dS_1\dots \mathrm dS_n $$ donde $p$ es una función de pago, $S_1,\dots, S_n$ son variables de pago (múltiples acciones o múltiples medidas de tiempo de las acciones o cualquier otra cosa), entonces (Q)MC también se puede utilizar para muestrear de $D$ para calcular la integral en $(1)$ numéricamente.

Por lo que he entendido del libro de Glasserman, sólo considera el QMC como una buena solución para el último método, en el que aporta algunas pruebas de que supera al MC aleatorio habitual. En cambio, para los dos primeros métodos habla de MC habitual y de técnicas de reducción de la varianza. Por lo tanto, mi pregunta con respecto a QMC es: ¿puede aplicarse con éxito en los dos primeros métodos, o no está bien diseñado para ellos y no debo esperar muchos beneficios del uso de QMC en los dos primeros métodos?

Answer 1

Por definición, el valor razonable de una opción viene dado por el valor de la expectativa del pago, $\mathbf{E}\left[\textrm{payoff}(\textit{paths})\right]$ . La distribución de probabilidad del paths es la medida de riesgo neutral. Esto es sólo una expresión integral de la forma que escribiste. Esto se aplica a todo precios de las opciones. Muchas opciones son, por supuesto, especiales en el sentido de que sólo dependen del valor terminal de cada camino.

Por lo tanto, siempre se puede aplicar Quasi Monte Carlo en estas situaciones, y el muestreo de trayectorias es esencialmente lo mismo que el cálculo de una integral (es decir, el valor de la expectativa del pago de la opción). Lo más probable es que el comentario de Glasserman se refiera al hecho de que si se utiliza simplemente QMC para generar trayectorias, es posible que no se esté recibiendo todo el beneficio de la baja discrepancia de los números aleatorios generados por QMC. Dicho de otro modo: las muestras tienen una discrepancia baja, pero las trayectorias generadas podrían no tenerla. Así que la reducción de la varianza podría no ser tan grande. Todo esto es relativo, en el sentido de que otros métodos de reducción de la varianza podrían producir mejores resultados (es decir, menor varianza en la estimación del precio).

Ahora, tu problema real no tiene nada que ver con QMC. Es Python el que te está frenando. Esquemas como Euler-Maryuama, Runge Kutta y Milstein son todos esquemas que no pueden ser vectorizados. Por lo tanto, sospecho que puedes tener algunos fragmentos de código de la forma

import numpy as np
... some definitions ...
for i in range(1, N_timesteps):
    for j in range(N_paths):
        S[i, j] = S[i - 1, j] * np.exp((r-0.5*sigma**2)*dt
                 + sigma*np.sqrt(dt) * np.random.randn())

El problema es que esto suele ser muy lento en Python. Cosas que puedes hacer para acelerar esto:

• Muestrear todos sus números aleatorios fuera del bucle.
• Precálculo de tantos factores fuera de los bucles como sea posible.
• Vectorizar los bucles internos, si es posible.
• Evite en lo posible las llamadas a funciones dentro del bucle.

Así que en general, queremos hacer lo menos posible dentro de los bucles for, ya que estas llamadas son muy caras en Python.

Por ejemplo, podríamos escribir lo anterior como

import numpy as np
... some definitions ...
Z = sigma*np.sqrt(dt) *np.random.randn(N_timesteps -1, N_paths)
drift = (r-0.5*sigma**2)*dt
for i in range(1, N_timesteps):
        S[i] = S[i - 1] * np.exp(drift + Z[i-1])

Así, el bucle interno se ha "vectorizado", el muestreo se realiza antes de los bucles, etc.

Si eso no es suficiente (lo que es muy probable si tiene muchos bucles for), entonces tendrá que realizar alguna optimización más "avanzada". Esto normalmente se reduce a tomar la pieza más lenta de su código (casi siempre el bucle for), ponerlo en una función y escribir esta función en, por ejemplo, Cython o Numba.

Por ejemplo, con Numba podríamos definir el bucle for en una función decorada como

from numba import autojit

@autojit
def my_loop(S, Z, drift):
    for i in range(1, S.shape[0]):
        for j in range(S.shape[1]):
            S[i, j] = S[i - 1, j] * np.exp(drift + Z[i-1, j])

... some definitions ...
Z = sigma*np.sqrt(dt) *np.random.randn(N_timesteps -1, N_paths)
drift = (r-0.5*sigma**2)*dt
my_loop(S, Z, drift)

La primera vez que se llame a my_loop, Numba compilará la función (tarda medio segundo). El resultado es una función que se ejecuta mucho, mucho más rápido. Es una gran manera de optimizar ciertas rutinas en su código.

Pero para ayudarte de verdad necesitaremos ver algo de código.