I am trying to use delta-gamma method with montecarlo simulations to calculate the VAR of a portfolio consisting in options and equities.
To use the method I need to compute a gamma matrix, that has gammas in its diagonal and cross-gammas in the i,j components.
My question is how do I compute the cross-gammas? I have looked everywhere but cannot find a closed formula or method.
BR
Consider an instrument value $f(S_0^1, \ldots, S_0^n)$ that depends on $n$ spot levels. Let $$\overrightarrow{S}_0=[S_0^1, \ldots, S_0^n]^T$$ be an $n$-dimensional vector representing the spot levels. We can approximate the cross gamma \begin{align*} \frac{\partial^2 f\big(\overrightarrow{S}_0\big)}{\partial S_0^i \partial S_0^j} \end{align*} by a finite difference scheme of the form \begin{align*} &\frac{1}{4 \varepsilon^2 S_0^i S_0^j}\Big[f\big(\overrightarrow{S}_0 + \varepsilon S_0^i\overrightarrow{1}_i + \varepsilon S_0^j\overrightarrow{1}_j\big) - f\big(\overrightarrow{S}_0 + \varepsilon S_0^i\overrightarrow{1}_i - \varepsilon S_0^j\overrightarrow{1}_j\big)\\ & \qquad\qquad - f\big(\overrightarrow{S}_0 - \varepsilon S_0^i\overrightarrow{1}_i + \varepsilon S_0^j\overrightarrow{1}_j\big) + f\big(\overrightarrow{S}_0 - \varepsilon S_0^i\overrightarrow{1}_i - \varepsilon S_0^j\overrightarrow{1}_j\big)\Big], \end{align*} where $\overrightarrow{1}_i$ is an $n$-dimensional vector with $1$ in the $i^{th}$ element and zeros elsewhere. Here, $\varepsilon$ is a small perturbation, which is usually set to $0.01$.
To justify the above, we use Taylor expansion. For example, \begin{align*} f\big(\overrightarrow{S}_0 + \Delta_i\overrightarrow{1}_i + \Delta_j\overrightarrow{1}_j\big)&\approx f\big(\overrightarrow{S}_0\big)+\Big(\Delta_i\frac{\partial }{\partial S_0^i}+ \Delta_j\frac{\partial }{\partial S_0^j}\Big)f\big(\overrightarrow{S}_0\big)\\ & \ \ + \frac{1}{2}\Big(\Delta_i\frac{\partial }{\partial S_0^i}+ \Delta_j\frac{\partial }{\partial S_0^j}\Big)^2f\big(\overrightarrow{S}_0\big) + o(\Delta_i \Delta_j ). \end{align*}
