Si quiere hacerlo de manera formal, los textos más básicos sobre el tema son Duan (1995) y Heston y Nandi (2000). En este último caso, en realidad proponen una fórmula cuasi-analítica para fijar el precio de las opciones europeas bajo su modelo GARCH. El modelo tiene esta forma bajo la medida física: \begin{align} ln S_{t+1} - ln S_t &= r_{ft+1} + \lambda h_{t+1} - \frac{1}{2} h_{t+1} + \sqrt{h_{t+1}} z_{t+1}, \; z_{t+1} \sim N(0,1) \\ h_{t+1} &= \sigma^2 + \pi( h_t - \sigma^2 ) + \alpha h_t(z_t^2 - 1 - 2 \gamma z_t ) \end{align} Lo escribo de forma diferente a Heston y Nandi (2000) porque quiero destacar la corrección de convexidad. En concreto, dada la normalidad de las perturbaciones, esta es la elección que se hace para garantizar que los rendimientos brutos tengan la forma que se espera, es decir, el tipo libre de riesgo más alguna compensación por el riesgo de volatilidad: \begin{align} E\left( \frac{S_{t+1}}{S_t} \big| F_t \right) = \exp \left( r_{ft+1} + \lambda h_{t+1} \right). \end{align}
Este no es en absoluto el mejor modelo que se puede encontrar, pero capta algunas cosas. En particular, si se calcula la correlación entre los rendimientos y la varianza, está controlada por $\gamma$ -- introduce una asimetría en la respuesta de la varianza que, a su vez, permite captar el "efecto palanca" (Black, 1976) por el que los rendimientos parecen estar correlacionados negativamente con las medidas de volatilidad y varianza (al menos para los índices bursátiles, parece cierto).
Ahora, para neutralizar el riesgo de este modelo, hay que elegir un núcleo de fijación de precios de forma explícita o implícita. Si se utiliza un núcleo de precios exponencialmente afín, básicamente hay que definir un choque distorsionado $\eta_{t+1}$ de tal manera que el rendimiento de tus acciones sea la tasa libre de riesgo. Entonces, se puede aprovechar el hecho de que este modelo implica la función característica condicional de $S_t$ es exponencialmente afín en las variables de estado ( $S_t$ y $h_t$ ) para obtener una fórmula de precios de forma cerrada. Dicha fórmula aparece en muchos textos. Puede consultar Heston y Nandi (2000) y utilizar la fórmula que dan, teniendo en cuenta $\lambda^{HN200} = \lambda^{Mine} - 1/2$ para que se adapte a la forma en que escribí el modelo anterior. En términos más generales, Christoffersen, Dorion, Jacobs y Wang (2010) detallaron y compararon muchos modelos de valoración de opciones GARCH. Christoffersen, Heston y Jacobs (2013) propusieron una forma de mejorar la fijación de precios mediante un núcleo cuadrático de fijación de precios. Bobaoglu, Christoffersen, Heston y Jacobs (2018) también propusieron combinar características interesantes como una estructura de componentes (se tiene un valor de tendencia a largo plazo para la varianza y un valor a corto plazo que fluctúa en torno a él) con perturbaciones no gaussianas (por lo que se obtienen colas gordas).
Ahora, ¿qué debes hacer? Bueno, mientras su modelo GARCH permanezca en la clase (exponencialmente) afín (es decir, la función característica es exponencialmente lineal en las variables de estado), existe una solución de forma casi cerrada para la fijación de precios. De lo contrario, se verá obligado a realizar una simulación de Montecarlo. En cualquier caso, la calibración del modelo para valorar o cubrir las opciones requiere datos de opciones (idealmente, tanto datos de acciones como de opciones). Se puede utilizar un estimador del método de los momentos para $\lambda$ forzándola a igualar el valor empírico de la prima de riesgo para alguna estimación de la misma y utilizando algún valor para la varianza a largo plazo del subyacente.
