Monte Carlo Greeks: Finite Differences, Pathwise Derivatives, and Likelihood Ratios

# How to get greeks using Monte-Carlo for arbitrary option?


# How to get greeks using Monte-Carlo for arbitrary option?












Let's assume I have an arbitrary option that I can price using Monte-Carlo simulation. What is the general approach (i.e. without relying on specific option type) to calculating the greeks in this case?

Edit: I woud like to add a few links on the topic that I found useful:

- Vibrato monte carlo and the calculation of greeks

- Conquering the Greeks in Monte Carlo



## Answer by Brian B (score 12, accepted)

https://quant.stackexchange.com/a/3334

You need to compute your greeks as finite differences, but the full procedure may be pretty tricky. I will use vega $\aleph$ as the example here. Let's begin by designating your Monte Carlo estimator as a function $V(\sigma,s,M)$ where $\sigma$ is the volatility as usual, $s$ is the seed to your random number generator, and $M$ is the sample count.

To begin with, recall that the Monte Carlo estimate of any value converges with the square root of the sample count. In particular, if you choose, say, $M=100$, you can run your estimator $N=500$ times to get estimate $\{V_n\}_{i=1}^{500}$, obtaining the standard deviation $\Sigma_{100}$ of those estimates.

Having done this, we now know the standard error of the estimator for any $M$ to be

$$ e_M \approx \Sigma_{100} \sqrt{\frac{100}{M}} $$

There are three possible cases:

- You can control the random seed $s$, or the set of random samples, used by the Monte Carlo estimator

- You cannot control $s$.

- You cannot even control the sample count $M$.

In the first case, you can use the fact that $s$ has been controlled to get a reasonable estimate of vega with relatively little extra work.

Find an $M$ such that the error in option price $e_M$ is tolerable. Choose a seed $s_0$ and a small increment $\Delta\sigma$ in the volatility, and compute

$$ \aleph^{(1)} = \frac{V(\sigma+\Delta\sigma,s_0,M)-V(\sigma-\Delta\sigma,s_0,M)}{2\Delta\sigma} $$ and $$ \aleph^{(2)} = \frac{V(\sigma+\frac12\Delta\sigma,s_0,M) - V(\sigma-\frac12\Delta\sigma,s_0,M)}{\Delta\sigma} $$

If $\aleph^{(1)} \approx \aleph^{(2)}$ then you have a good estimate and you are done.

The reason this works so nicely is that, by controlling the seed, our difference computations

$$ \delta=V(\sigma+\Delta\sigma,s_0,M)-V(\sigma-\Delta\sigma,s_0,M) $$ are direct Monte Carlo estimators of the vega, since the shared seed implies the samples $x_i$ match in the difference of sums. That is $$ \delta = ( \frac1M \sum_{i=1}^M f(x_i, \sigma+\Delta\sigma) ) -( \frac1M \sum_{i=1}^M f(x_i, \sigma-\Delta\sigma) ) \\ =\frac1M \sum_{i=1}^M f(x_i, \sigma+\Delta\sigma)-f(x_i, \sigma-\Delta\sigma) $$

The second case where you cannot control the seed, on the other hand, is rather more difficult. Here, you will have a different error $e$ to the true value every time you run the function.

For brevity, let's let $$ e_\pm = V(\sigma\pm \Delta\sigma,s_\pm,M). $$

Of course we do not know the value of $e_\pm$ or of $s_\pm$, but we do at least have our estimate of the size of $e_\pm$ as noted above. Therefore, the error in $\delta$ is approximately $e_M \sqrt{2}$. You need to choose $M$ so large that $$ \delta \gg e_M \sqrt{2}. $$

Not knowing the value of $\delta$ a priori makes this difficult, but usually in a trading context one can specify an acceptable absolute error $\epsilon$ in vega. In that case, we can demand $$ \epsilon < \frac{e_M \sqrt{2}}{\Delta\sigma} $$ which translates to $$ M > \Sigma_{100}^2 {\frac{200}{\epsilon^2 \Delta\sigma^2}}. $$

The third case, where you can control neither the random seed $s$ nor the sample count $M$ should be treated as the second case above. You simply treat each run of the algorithm as a single sample.

## Answer by Dirk Eddelbuettel (score 5)

https://quant.stackexchange.com/a/3330

The most general answer is to shift your input to approximate the first derivative. Given that you need Monte Carlo to price this, it may get expensive. But that's the way it goes as when you have no analytical solutions as there aint't no free lunch ...

## Answer by Mark Joshi (score 3)

https://quant.stackexchange.com/a/15436

if the pay-off is continuous, the standard approach is to use the path-wise method also known as IPA. This essentially means that you differentiate along each path. It is the limit as the bump size goes to zero of finite differencing.

The main downside of this method is that the differentiation can be fiddly and slow. The Smoking adjoints paper you mention made the observation that using adjoint/automatic differentiation makes it fast and indeed there are packages that will do the differentiating automatically for you. There is a survey article by Homescu on this topic.

Adjoints and Automatic (Algorithmic) Differentiation in Computational Finance

Cristian Homescu

I have also written too many papers on this. See my SSRN page http://papers.ssrn.com/sol3/cf_dev/AbsByAuth.cfm?per_id=550354

I also devote a chapter of More Mathematical Finance to this.

When the pay-off is discontinuous or you want Hessians life is more complicated. Likelihood ratio is one approach but can lead to large variances. This is particularly the case for vegas.

See Glasserman (2003) Monte Carlo Methods in Financial Engineering for general discussion.

I have done various papers on the discontinuous case and on Hessians. My current favourites are

http://ssrn.com/abstract=2431580

http://ssrn.com/abstract=2011690

Essentially, these work by doing change of variables that remove the discontinuities and so then the pathwise method can be applied. They are designed to yield the most pathwise mix of likelihood ratio and pathwise that avoids the discontinuity.

## Answer by Keith A. Lewis (score 2)

https://quant.stackexchange.com/a/3455

Download the code from http://fmsoption.codeplex.com to see how to do that for vanilla options. You are right, you need implementations for transcendental functions that are written for dual numbers. You will find them in the fmsdual project.

If you just want browse some source code, see http://fmsoption.codeplex.com/SourceControl/changeset/view/10924#145366. Note that `eps` is machine epsilon ~= 2e-16. (!)

## Answer by Keith A. Lewis (score 0)

https://quant.stackexchange.com/a/3363

Another way to do this is to use dual numbers. http://fmsdual.codeplex.com.

They let you calculate an arbitrary number of derivatives while running a single Monte Carlo. Here is an example of how to use it:

```
// Monte Carlo derivatives
void fms_test_monte(size_t N)
{
    ::srand(static_cast<unsigned int>(::time(0)));

    double a = 0.5;
    dual::number<double,3> A(a, 1);
    dual::number<double,3> E(0.,1);

    for (int i = 0; i < N; ++i) {
        double x = 1.0*rand()/RAND_MAX;
        E = E + (x - A)*(x - A);
    }

    E = E/(1.*N);

    // X uniform [0,1]
    // E(X - a)^2 = 1/3 - 2a 1/2 + a^2
    ensure (fabs(E._(0) - (1./3 - a + a*a)) < sqrt(1./N));

    // d/da E(X - a)^2 = -2 E(X - a) = 2a - 1
    ensure (fabs(E._(1) - (2*a - 1)) < sqrt(1./N));

    // d^2/da^2 E(X - a)^2 = 2
    ensure (fabs(E._(2) - 2) < sqrt(1./N));
}
```