Estimar los parámetros cuando solo se conoce la función característica.

Question

Recientemente estuve trabajando con un proceso llamado Variance Gamma con Llegada Estocástica (VGSA) e intentando ajustar este proceso a unos datos dados.

Para obtener VGSA, como se explica en Carr et al. [2001], tomamos el proceso VG que es un proceso de Lévy homogéneo y agregamos volatilidad estocástica evaluándolo en un cambio de tiempo continuo dado por la integral de un proceso Cox, Ingersoll y Ross [1985] (CIR). Formalmente definimos el proceso CIR $y(t)$ como la solución a la ecuación diferencial estocástica, $$ d y_t=\kappa\left(\eta-y_t\right) d t+\lambda \sqrt{y_t} d W_t $$ donde $W_t$ es una marcha Browniana, $\eta$ es la tasa a largo plazo de cambio de tiempo, $\kappa$ es la tasa de reversión a la media, y $\lambda$ es la volatilidad del cambio de tiempo. El proceso $y(t)$ es la tasa instantánea de cambio de tiempo y por lo tanto el cambio de tiempo está dado por $Y(t)$ donde $$ Y(t)=\int_0^t y(u) d u $$

La función característica para $Y(t)$ está dada por $$ \begin{aligned} \mathbb{E}\left(e^{i u Y(t)}\right) & =\phi(u, t, y(0), \kappa, \eta, \lambda) \\ & =A(t, u) e^{B(t, u) y(0)} \end{aligned} $$ donde $$ \begin{aligned} A(t, u) & =\frac{\exp \left(\frac{\kappa^2 \eta t}{\lambda^2}\right)}{\left(\cosh (\gamma t / 2)+\frac{\kappa}{\gamma} \sinh (\gamma t / 2)\right)^{2 \kappa \eta / \lambda^2}} \\ B(t, u) & =\frac{2 i u}{\kappa+\gamma \operatorname{coth}(\gamma t / 2)} \end{aligned} $$ con $$ \gamma=\sqrt{\kappa^2-2 \lambda^2 i u} $$

El proceso de Lévy de volatilidad estocástica, denominado proceso VGSA, está definido por $$ \begin{aligned} Z(t) & =X_{V G}(Y(t) ; \sigma, \nu, \theta) \\ & =\theta \gamma(Y(t) ; 1, \nu)+\sigma W(\gamma(Y(t) ; 1, \nu)) \end{aligned} $$

Por lo tanto $\sigma, \nu, \theta, \kappa, \eta$, y $\lambda$ son los seis parámetros que definen el proceso. Su función característica está dada por $$ \mathbb{E}\left(e^{i u Z_{V G S A}(t)}\right)=\phi\left(-i \Psi_{V G}(u), t, \frac{1}{\nu}, \kappa, \eta, \lambda\right) $$ donde $\phi$ es la función característica de $Y(t)$ dada en (6.68) y $\Psi_{V G}$ es la función característica logarítmica del proceso de Variance Gamma en tiempo unitario, es decir, $$ \Psi_{V G}(u)=-\frac{1}{\nu} \log \left(1-i u \theta \nu+\sigma^2 \nu u^2 / 2\right) $$

Definimos el proceso de acción en el tiempo $t$ por la variable aleatoria $$ S(t)=S(0) \frac{e^{(r-q) t+Z(t)}}{\mathbb{E}\left[e^{Z(t)}\right]} $$

Notamos que $$ \mathbb{E}\left[e^{Z(t)}\right]=\phi\left(-i \Psi_{V G}(-i), t, \frac{1}{\nu}, \kappa, \eta, \lambda\right) $$ lo cual es equivalente a $e^{-\omega t}$ en el caso de VG. Por lo tanto la función característica del $\log$ del precio de la acción en el tiempo $t$ está dada por $$ \mathbb{E}\left[e^{i u \log S_t}\right]=\exp \left(i u\left(\log S_0+(r-q) t\right)\right) \times \frac{\phi\left(-i \Psi_{V G}(u), t, \frac{1}{\nu}, \kappa, \eta, \lambda\right)}{\phi\left(-i \Psi_{V G}(-i), t, \frac{1}{\nu}, \kappa, \eta, \lambda\right)^{i u}} $$

Para obtener información adicional, consulte la página 236 del libro titulado "Métodos Computacionales en Finanzas" por Ali Hirsa, específicamente en relación al Proceso Estocástico Variance-Gamma (VGSA).

Mi enfoque implica inicialmente transformar por Fourier inversa la función característica y luego utilizar el método de Máxima Verosimilitud. A pesar de usar una suposición inicial y el método de Nelder–Mead para minimización, los resultados obtenidos suelen ser subóptimos, probablemente debido a la convergencia a mínimos locales.

Busco orientación sobre un método más robusto para ajustar el modelo VGSA a mis datos para estimación de parámetros, ya que el enfoque actual parece quedarse atascado en mínimos locales. Cualquier sugerencia e implementaciones de código en Python, R o MATLAB serían muy apreciadas.

Answer 1

Observaciones generales

Un problema difícil pero muy interesante. Algunas ideas:

• Usa GMM en lugar de MLE. MLE es un caso especial del método generalizado de momentos. La función característica te da los momentos a través de la diferenciación ($\mathbb{E}[X^n]=i^{-n}\varphi_X^{(n)}(0)$). Por lo tanto, no necesitas invertir la función característica. Sin embargo, como hay seis parámetros, necesitarías 6-7 momentos que pueden ser difíciles de estimar a partir de una serie temporal.

• ¿Has intentado usar la transformada rápida de Fourier? La velocidad aumentada podría darte más capacidad para enfocarte en la precisión y estabilidad. Carr, Geman, Madan y Yor (2002, p. 320) explican la estimación de su modelo CGMY a partir de datos históricos de la siguiente manera

Para cada activo subyacente, formamos la serie temporal de logaritmos diarios de los precios relativos y luego estimamos los parámetros de la densidad de Lévy $C, G, M, Y,$ y $\eta$ a partir de los datos de retorno ajustados a la media. La estimación directa por máxima verosimilitud es computacionalmente costosa, ya que requiere una inversión de Fourier para cada punto de datos para evaluar la densidad, y estas inversiones deben anidarse en un algoritmo de optimización de búsqueda de gradiente para la estimación de parámetros. La transformada rápida de Fourier se usó para invertir la función característica una vez para cada configuración de parámetros. Este método representa eficientemente el nivel de la densidad de probabilidad en un conjunto predeterminado de valores de retornos. (...) Con la densidad evaluada en estos puntos predeterminados, agrupamos la serie de retornos contando el número de observaciones en cada punto de retorno predeterminado, asignando las observaciones de datos al punto de retorno predeterminado más cercano. Luego buscamos estimaciones de parámetros que maximizaran la verosimilitud de estos datos agrupados. Por lo tanto, las estimaciones informadas son para esta estimación de máxima verosimilitud agrupada utilizando la transformada rápida de Fourier.

• Probablemente, usar el método COS en lugar de FFT sería aún mejor.

• De hecho, hay una literatura sobre cómo disciplinar mejor la estimación por MLE o usar estimadores alternativos (por ejemplo, MCMC). Consulta aquí para una breve descripción.

• No es de mucha ayuda, pero quizás quieras considerar usar un modelo con densidad de forma cerrada (como Variance Gamma).

---

Detalles sobre el método COS

Empiezo citando un conocido lema del análisis usando la transformada de Fourier como aproximación de la función original.

Sea $f:[a,b]\to\mathbb{R}$ una función integrable con transformada de Fourier $\hat{f}$. Entonces, para $N\in\mathbb{N}$, tenemos \begin{align} f(x) &\approx \frac{2}{b-a} \sum_{n=0}^N \tilde{A_n} \cos\left( n\pi\frac{x-a}{b-a}\right), \\ \tilde{A_n} &= \begin{cases} \text{Re}\left(\hat{f}\left(\frac{n\pi}{b-a}\right)\cdot e^{-i\frac{n\pi a}{b-a}}\right) & \text{si } n\geq1, \\ \frac{1}{2} \text{Re}\left(\hat{f}\left(\frac{n\pi}{b-a}\right)\cdot e^{-i\frac{n\pi a}{b-a}}\right) & \text{si } n=0. \end{cases} \nonumber \end{align}

Es fácil de demostrar. Todo el crédito a Fang y Oosterlee (2008). Pero ya puedes ver hacia dónde voy. $\hat{f}$ es la función característica conocida, $f$ es la función de densidad de probabilidad buscada de $\ln(S_t)$. Críticamente, nota que no hay inversión de Fourier, no hay integral, nada. Solo una suma. Eso significa que puedes ser muy rápido y enfocarte en la precisión y estabilidad.

Un problema es la elección de $N$ y los límites $a$ y $b$. Claramente, la mayoría de las funciones de densidad no tienen soporte compacto. Sin embargo, incluso si tienen colas gordas, normalmente se pueden truncar con seguridad. Si un precio de acción logarítmico está en $\ln(10)$, es muy improbable que necesites llegar hasta $\ln(1,000)$. Fang y Oosterlee tienen sugerencias sobre cómo elegir $a$ y $b$ basándose en los cumulantes de $\ln(S_t)$ pero no los conocemos sin los parámetros, así que simplemente tomaría valores elegidos generosamente y haría un análisis de sensibilidad. Elige algunos, estima los parámetros usando MLE y luego cambia $a$ y $b$ y reestima los parámetros. Si no cambia nada, estás bien.

A continuación se muestra un gráfico para la distribución normal estándar con densidad $f:[-10,10]\to\mathbb{R},x\mapsto \frac{1}{\sqrt{2\pi}}e^{-\frac{1}{2}x^2}$ y transformada de Fourier $\hat{f}(u)=e^{-\frac{1}{2}u^2}$. Como puedes ver, converge muy rápido. Para $N=20$, las gráficas son visualmente indistinguibles.