I have written the following python code to price a european call option. I get reasonable agreement with an online option calculator for the price, but the delta is quite a bit off. As far as I understand, I simulate a random walk, calculate the pay-off, present value the payoff, and then do lots of such simulations and take the average to get the price of the option. In order to calculate the delta, I simulate two paths, where the second path starts at a price [$] S + \Delta S [$], but the random number has to be the same. Can anyone confirm the correctness of the following pseudo python code? This is not the exact code, it just contains the essence of the algorithm. strike = 105r = 0.05sigma = 0.2time_to_expiry = 2time_step = 0.01N = time_to_expiry/time_stepdS = 0.5S[0] = 100S2[0] = S[0] + 2*dSfor i in xrange(1, N): S_current = S[i-1] S2_current = S2[i-1] #Draw a random number from a normal distribution with mean 0 and variance 1 rand = np.random.normal(loc=0.0, scale=1.0) #Euler discretisation of risk-neutral random walk S = S_current*(1 + dt*r + sigma*np.sqrt(dt)*rand) S2 = S2_current*(1 + dt*r + *sigma*np.sqrt(dt)*rand) S_call_payoff = max(S[N-1] - strike, 0) S_put_payoff = max(strike - S[N-1], 0) #Calculate greeks S2_call_payoff = max(S2[N-1] - strike,0)S2_put_payoff = max(strike - S2[N-1],0)call_delta = (S2_call_payoff - S_call_payoff)/(2*dS)put_delta = (S2_put_payoff - S_put_payoff)/(2*dS) call_price = np.exp(-r*(time_to_expiry))*S_call_payoffput_price = np.exp(-r*(time_to_expiry))*S_call_payoffThen one must average over many such simulations. Have I calculated the delta correctly? For the above parameters, an online option calculator gives:Call price: 13.6537, Put price: 8.6486Call delta: 0.6266, Put Delta: -0.3734In contrast, my Monte-Carlo code with 10,000 samples gives:Call price: 14.05, Put price: 8.57Call delta: 0.71, Put Delta: -0.4They're in the right ball-park, but I'm not confident enough to say the code is working correctly. Any help much appreciated.

For the MC, your naïve approach won't work. Kienitz and Wetterau book discuss how to it. The ds factor is crucial based on Jaeckel 2005 "more likely than not" article.However, your fd would better using lattice The Euler method is awful in general, so fd (awful) == awful. This folktale can be made mathematically precise. edit: K&W provide Matlab code, so a port to Python could be a good way to learn it.

There are three issues here. First, I understand that you may want to do time stepping for some other purpose, but it adds noise in the way that you are doing it. You should simulate the log stock price process (and then exponentiate it). This way there is no discretization error, regardless of the length of your time steps. Second, you are taking a relatively large spatial step to estimate your deltas, picking up a bit of gamma in the process. Third, 10,000 paths is not a lot, given that you are not employing any variance reduction techniques. All of that being said, your basic approach should work in the limit as you take shorter time steps, employ more paths and shrink the bump that you are using to calculate deltas.

QuoteOriginally posted by: CuchulainnFor the MC, your naïve approach won't work. Kienitz and Wetterau book discuss how to it. The ds factor is crucial based on Jaeckel 2005 "more likely than not" article.However, your fd would better using lattice The Euler method is awful in general, so fd (awful) == awful. This folktale can be made mathematically precise. edit: K&W provide Matlab code, so a port to Python could be a good way to learn it.Thanks Cuch, do you know of any free resources that explain it? K&W is quite pricey.

QuoteOriginally posted by: bearishThere are three issues here. First, I understand that you may want to do time stepping for some other purpose, but it adds noise in the way that you are doing it. You should simulate the log stock price process (and then exponentiate it). This way there is no discretization error, regardless of the length of your time steps. Second, you are taking a relatively large spatial step to estimate your deltas, picking up a bit of gamma in the process. Third, 10,000 paths is not a lot, given that you are not employing any variance reduction techniques. All of that being said, your basic approach should work in the limit as you take shorter time steps, employ more paths and shrink the bump that you are using to calculate deltas.Thanks, I did the Euler discretisation precisely because I want to learn about the downfalls of the method compared with others. I was also of the understanding that although a risk-neutral random walk has an exact solution, most SDE's don't, and therefore you'd need to do something like Euler-M anyway, hence I wanted to get some intuition for it this way.

QuoteOriginally posted by: barnyQuoteOriginally posted by: CuchulainnFor the MC, your naïve approach won't work. Kienitz and Wetterau book discuss how to it. The ds factor is crucial based on Jaeckel 2005 "more likely than not" article.However, your fd would better using lattice The Euler method is awful in general, so fd (awful) == awful. This folktale can be made mathematically precise. edit: K&W provide Matlab code, so a port to Python could be a good way to learn it.Thanks Cuch, do you know of any free resources that explain it? K&W is quite pricey.I believe the code is freely available. AFAIK And Jaeckel "more... " is on the Web.

Yep I already found it but code without context is pretty meaningless Either way, I assume they have assimilated this knowledge from other books or papers to put into their book. What is the name of the method?

K&W call it the finite difference method using difference estimator. And Glasserman 2004, of course discusses all this stuff.BTW JK == @lapsilago is a Wilmott member.

Another way to do the MC delta (after switching to X = log S as bearish suggested) is to simply treat it as another option value with its own payoff.This avoids the complication of trying to do a numerical differentiation within a MC.

QuoteOriginally posted by: AlanAnother way to do the MC delta (after switching to X = log S as bearish suggested) is to simply treat it as another option value with its own payoff.This avoids the complication of trying to do a numerical differentiation within a MC.Bit confused by this, can you elaborate? Do I simulate the same risk-neutral random walk, but alter the pay-off somehow? In which case, what is the pay-off of the delta of an option?

QuoteOriginally posted by: barnyQuoteOriginally posted by: AlanAnother way to do the MC delta (after switching to X = log S as bearish suggested) is to simply treat it as another option value with its own payoff.This avoids the complication of trying to do a numerical differentiation within a MC.Bit confused by this, can you elaborate? Do I simulate the same risk-neutral random walk, but alter the pay-off somehow? In which case, what is the pay-off of the delta of an option?Good question: I was wondering about that myself Like, how does it work out?QuoteThis avoids the complication of trying to do a numerical differentiation within a MC.ND is an ill-conditioned process. Taking smaller dS, h gives badder results.// As a follow-on, Can you differentiate the BS PDE in S to get a PDE for the delta? It might be mathematical heresy

QuoteOriginally posted by: barnyQuoteOriginally posted by: AlanAnother way to do the MC delta (after switching to X = log S as bearish suggested) is to simply treat it as another option value with its own payoff.This avoids the complication of trying to do a numerical differentiation within a MC.Bit confused by this, can you elaborate? Do I simulate the same risk-neutral random walk, but alter the pay-off somehow? In which case, what is the pay-off of the delta of an option?There are a couple different ways to do it. Assuming no dividend on the stock: 1. From the BS formula, S0 delta(S0) = C(S0) + K e^(-r T) Phi(d2(S0)), so you can simulate the same random walk and valueC(S0) and e^(-r T) Phi(d2(S0)), which both satisfy the BS pde and have obvious payoff values. (See my follow-up postfor a simpler version of the same thing, where the payoff is displayed).or2. Differentiate the BS pde to get the pde followed by delta (This answer Cuch's question, too). For that one, you will need tosimulate a different random walk as the drift differs. Again, the payoff value should be obvious: it's just the delta at expiration.

QuoteOriginally posted by: outrunIf the MC price is $$P = \sum_{i=1..n} p_i / n$$then$$dP/dS = \sum_{i=1..n} (dp_i/dS_i) / n$$ for a vanilla call where $$p_i = max(S_i - K,0)e^{-rt}$$it translates to$$dp_i/dS_i = e^{-rt} \text{ if } S_i >= K\text{, else }0$$Max(S-K,0) is it not differentiable at S = K.