Estimating Option Gamma by Monte Carlo Finite Differences
Summary
The document explains why a finite difference estimate of a vanilla call option’s gamma can fail when implemented path by path in a Monte Carlo simulation. The estimator compares call payoffs under upward-shifted, unshifted, and downward-shifted initial prices, using the same simulated random draw for each path. The questioner incorrectly assumes that all three terminal prices remain on the same side of the strike whenever the unshifted price is there, causing the payoff second difference to collapse to zero in many paths.
The answer gives the condition under which the shifted paths cross the strike and shows why paths near the strike matter. As the price bump shrinks, the pathwise second difference approaches a Dirac delta expression; almost every individual simulation then contributes zero, even though the expectation is nonzero and depends on the terminal-price density at the strike. This explains why direct pathwise Monte Carlo can be unreliable for gamma. The discussion points to other resources but does not compare alternative estimators or provide a complete numerical procedure.
Key ideas
- A finite difference gamma estimator compares option payoffs across bumped and unbumped initial prices using common random draws.
- The shifted terminal prices can cross the strike even when the unshifted terminal price lies on one side.
- Paths near the strike drive the second difference in payoff.
- As the bump shrinks, the pathwise limit is concentrated at the strike, making ordinary Monte Carlo estimates problematic.
- The expected limiting quantity can be nonzero because it depends on the terminal-price density at the strike.
Tags
Full text
# Greeks: Estimate gamma by Monte Carlo finite difference
# Greeks: Estimate gamma by Monte Carlo finite difference
When I was using Monte Carlo to calculate the gamma of a vanilla call option by finite difference method, I stuck in this weird situation as below.
Consider this, $$ Gamma = \frac{CallPrice(S^{up}_{T}) - 2 * CallPrice(S_{T}) + CallPrice(S^{down}_{T})}{dS^2} $$ And we can choose dS small enough such that when $$ S_{T}>K \text{ then } S^{down}_{T}>K $$ and $$ S_{T}<K \text{ then } S^{up}_{T}<K $$
That is, we can write the above Gamma formula as $$ Gamma = \frac{(S^{up}_{T}-K)I(S_{T}>K) - 2 * (S_{T}-K)*I(S_{T}>K) + (S^{down}_{T}-K)I(S_{T}>K)}{dS^2} $$ $$ = \frac{(S^{up}_{T} - 2 * S_{T} + S^{down}_{T})I(S_{T}>K)}{dS^2} $$ $$ = \frac{(S_{0}+dS)*exp(...) - 2*S_{0}*exp(...) + (S_{0}-dS)*exp(...)}{dS^2} = 0 $$
So every time I run the simulation, I always get a correct delta but wrong gamma. (Not equal to zero maybe because of rounding error?)
I know gamma is nonzero, but I can't find where I did it wrong. Any help?
Note: This question is a bit similar to this one "Greeks: Why does my Monte Carlo give correct delta but incorrect gamma?", but slightly different.
## Answer by ir7 (score 3, accepted)
https://quant.stackexchange.com/a/57713
Pathwise finite difference Gamma formula is indeed:
$$\Gamma(S_0,T, dS; Z) = (dS)^{-2} \left[ (S_T^{up} (Z) - K)^+ -2 (S_T (Z) - K)^+ + (S_T^{dn} (Z) - K)^+ \right], $$
where $Z$ is a standard normal rv, and
$$ S_T (Z) = S_0\eta (Z), $$
$$ S_T^{up} (Z) = (S_0+dS)\eta(Z) = (S_0+dS)S_0^{-1} S_T (Z) $$
$$ S_T^{dn} (Z) = (S_0-dS)\eta(Z) = (S_0-dS)S_0^{-1} S_T (Z), $$
and
$$ \eta(Z) = \exp \left( (r-0.5\sigma^2)T + \sigma \sqrt{T} Z \right). $$
For realizations of $Z$, if $ S_T (Z) > K$, then $ S_T^{up} (Z) > K$. We also have $ S_T^{dn} (Z) > K$, but only if:
$$ dS < (S_T(Z) - K)\eta(Z)^{-1}.$$
If $ S_T (Z) <K$, then $ S_T^{dn} (Z) <K$. And $ S_T^{up} (Z) <K$, if
$$ dS < (K -S_T(Z) )\eta(Z)^{-1}. $$
So, for such $Z$ realization and $dS$, we do have:
$$ I(S_T(Z)-K) = I(S_T^{up}(Z)-K) = I(S_T^{up}(Z)-K). $$
If $ S_T (Z) = K$, then $ S_T^{up} (Z) > K$, and $ S_T^{dn} (Z) < K$ no matter what $dS>0$ we choose. So, in this case:
$$ I(S_T(Z)-K) = I(S_T^{dn}(Z)-K) =0 \not= 1 = I(S_T^{up}(Z)-K). $$
The link that you provided offers answers and resources on how to deal with Gamma calculations in Monte Carlo framework.
Note that, with notations above,
$$ \lim_{dS\rightarrow 0}\Gamma(S_0,T, dS; Z) = \delta (S_T(Z)-K) S_T(Z)S_0^{-1}, $$
where $\delta$ is the Dirac delta function. So, computing its expectation via Monte Carlo is bound to give meaningless results. The probability of giving a non-zero result is zero, while in reality we know that:
$$ E[\delta (S_T(Z)-K) S_T(Z)S_0^{-1}] = KS_0^{-1}f(K), $$
where $f$ is the pdf of $S_T(Z)$.Shown in full with attribution under the source's licence. Licence: CC BY-SA 4.0 (Stack Exchange)
This summary was written by Stratmill's research agent from the original; it is not a copy of the source.