Choice of epsilon for numerical calculation of vega in binomial option pricing model

Question

I have a binomial option-pricing model (I don't think the details of how its implemented are relevant). However, when I go to calculate vega, I am essentially running the model a second time with new volatility input we call $\sigma + \epsilon$, where $\sigma$ is the volatility used to calculate the price. I then have two prices $p1$ and $p2$, the first of which is a function of $\sigma$, and the second a function of $\sigma + \epsilon$. I then calculate vega to be $\frac{(p2 - p1)}{\epsilon}$.

I'm having trouble coming up with a good value of $\epsilon$ to avoid floating-point underflows and overflows in the resulting calculation. Any suggestions on how to choose it?

Some things that I have available to me at the time of the choice of $\epsilon$:

• expiration
• strike
• volatility
• risk-free rate
• underlying price
• current time
• option price (as calculated by my model)
• delta (as calculated by my model)
• gamma (as calculated by my model)

Answer 1

Let $\beta$ denote the relative machine precision, usually $\beta = 1E-16$. Assume the you can evaluate the value V up to precision $\alpha$. The best you can get is $\alpha = \beta \cdot V$ if $V$ is not underflow or overflow. Then you can calculate the finite difference up to precision $4 \alpha / \epsilon$ (the 4 might be a rough estimate, but it comes from the fact that there may be an additional cancelation in the finite difference of relative order $\beta$ and we use $\alpha > \beta \cdot V$.

Thus you can calculate the finite difference up to $4 \alpha / \epsilon$. On the other hand, from Taylor expansion, the approximation error of a finite difference is $C \epsilon$, where $C$ is a bound on the second derivative. The best choice of $\epsilon$ is the minimum of $\epsilon \mapsto 4 \alpha / \epsilon + C \epsilon$, which is attained for $-4\alpha/\epsilon^2 + C = 0$, i.e. $\epsilon = 2 \cdot \sqrt{\alpha/C}$. If your tree achieves machine accuracy, i.e., $\alpha = \beta V$, then you should choose $\epsilon = 2 \sqrt{\frac{V}{V''}} \cdot \sqrt{\beta} = 2 \sqrt{\frac{V}{V''}} 1E-8$. In the last equation V denote "the order of magnitude of the value and $V''$ the order of magnitude of the second derivative. Clearly, if the second derivative is small, you can use larger shifts.