He buscado a través de docenas de documentos y no encontré en ninguno de ellos respuestas satisfactorias y suficientemente teóricas a mis preocupaciones. Así que he combinado todo lo que he encontrado a continuación. Por favor, indicad si mi comprensión del tema es adecuada y corregidme si es necesario. Además de la parte teórica he encontrado también un obstáculo en la aplicación práctica.
Duan (1995) en su trabajo desarrolló un modelo para la valoración de opciones europeas con el modelo GARCH. Introdujo la medida de riesgo neutral local (LRNVR) $\mathbb{Q} $ que equivale a la medida del mercado físico $\mathbb{P} $ . Dado que el precio de la opción de compra o de la opción de venta no tiene una solución analítica en su marco, es necesario realizar simulaciones de Monte Carlo. Mi primera preocupación para la que no he encontrado una respuesta explícita es la siguiente:
(1) Tenemos una serie temporal con algunos datos de mercado y ajustamos los parámetros del modelo a esta serie temporal bajo la medida física $ \mathbb{P} $ .
(2) A continuación, utilizamos los parámetros estimados para el proceso transformado bajo la medida LRNVR $\mathbb{Q} $ y realizamos simulaciones de Monte Carlo para estimar el precio de una opción.
Dado que la simulación MC bajo medida $\mathbb{Q} $ no es un problema para mí, centrémonos ahora en la derivación teórica de la MLE bajo la medida física $\mathbb{P} $ .
Supongamos que tenemos una muestra de $T$ log-returnos de algún activo financiero. Sea $X = (X_1, \dots , X_T) $ denotan nuestra muestra y suponen que $t$ -El logaritmo de retorno sigue el proceso GARCH(1,1)-M propuesto por Duan (1995): $$ X_t = \ln \frac{S_t}{S_{t-1}} = r + \lambda \sqrt{h_t} - 0.5h_t + \sqrt{h_t} z_t $$ bajo medida física $\mathbb{P} $ , donde $z_t \overset{iid}{=} \mathcal{N} (0,1) $ y $h_t = \omega + \alpha h_{t-1} z_{t-1}^2 + \beta h_{t-1} $ . Aquí asumimos que $ \omega >0 $ y $ \alpha, \beta \geq 0 $ . También requerimos $ \alpha + \beta < 1$ para garantizar la estacionariedad. Parámetro $r$ es un tipo de interés de mercado sin riesgo (parámetro conocido) y $\lambda$ es una prima de riesgo asociada a un determinado activo financiero (parámetro a estimar).
Ya que tenemos $ X_t \sim \mathcal{N} (r + \lambda \sqrt{h_t} - 0.5h_t, h_t )$ entonces la función de probabilidad para $t$ -la observación es: $$ l_t (X_t ; \theta ) = \frac{1}{\sqrt{2 \pi h_t}} \exp \left( - \frac{ ( X_t - r - \lambda \sqrt{h_t} + \frac{1}{2} h_t )^2 }{2h_t} \right) \text{,} $$
donde $\theta = (\omega, \alpha, \beta, \lambda ) $ es un vector de parámetros a estimar. Función de verosimilitud para un vector $X$ es: $$ \mathcal{L} (X; \theta ) = \prod_{t=1}^T l_t(X_t; \theta) \text{.} $$ Como es más fácil calcular los logaritmos naturales, entonces: $$ \ln \mathcal{L} (X; \theta ) = -\frac{T}{2} \ln \left( 2 \pi \right) - \frac{1}{2} \sum_{t=1}^T \ln \left( h_t \right) - \frac{1}{2} \sum_{t=1}^T \frac{ ( X_t - r - \lambda \sqrt{h_t} + \frac{1}{2} h_t )^2 }{2h_t} $$ para ser maximizado. Estamos buscando $ \hat{\theta} $ que es: $$ \hat{\theta} = \arg \max_{\theta} \ln \mathcal{L} (X; \theta ) $$ con restricciones para $\omega , \alpha, \beta $ como se ha indicado anteriormente.
Pasemos ahora a la aplicación práctica de lo anterior:
Dejemos que nuestros datos de mercado sean las rentabilidades logarítmicas diarias de AAPL para el periodo 2016-2019 (o cualquier otro dato, porque el siguiente problema no desaparece con el cambio de datos de origen). Supongamos que el tipo de interés libre de riesgo es $r=0$ . Como varianza inicial $h_1$ asumimos la varianza de nuestra muestra, es decir $h_1 = Var(X) $ . La función a minimizar es $ - \ln \mathcal{L} (X; \theta ) $ y se define como sigue:
loglike <- function(params, log_returns){
omega <- params[1]
alpha <- params[2]
beta <- params[3]
lambda <- params[4]
bigT <- length(log_returns)
h <- c(var(log_returns))
for (i in 2:bigT) {
h[i] <- omega + alpha * (log_returns[i-1] - lambda * sqrt(h[i-1]) + 0.5 * h[i-1] )^2 + beta * h[i-1]
}
likelihood <- 0.5 * bigT * log(2*pi) + 0.5 * sum(log(h) + ((log_returns - lambda * sqrt(h) + 0.5 * h)^2) / h )
return(likelihood)
}
Los parámetros de partida y las restricciones del problema de optimización son los siguientes:
params <- rep(0.01, 4)
lb <- c(0,0,0, -Inf)
A = cbind(0, 1, 1, 0)
b = 1
Y si trato de usar fmincon de practica para la optimización tengo el siguiente error y advertencias:
> fmincon(loglike, log_returns = log_returns, x0 = params, lb = lb, A = A, b = b)
Error in if (f < 0) { : missing value where TRUE/FALSE needed
In addition: Warning messages:
1: In sqrt(h[i - 1]) : NaNs produced
2: In log(h) : NaNs produced
3: In sqrt(h) : NaNs produced
En primer lugar, no sé por qué recibo estos avisos. Según la función objetivo definida anteriormente $ h \geq 0 $ por lo que tomar raíces cuadradas y logaritmos naturales no debería producir NaNs.
En segundo lugar no entiendo el error devuelto por fmincon . ¿Qué ocurre con mi función objetivo?