import numpy as np #input s0 = 50 r = 0.06 sigma = 0.3 T = 1 M = 50 dt = T/M u = np.exp(sigma*np.sqrt(dt)) d = 1/u a = np.exp(r*dt) p = (a-d)/(u-d) #option price at time t = M s = np.zeros(M+1) min1 = np.zeros(M+1) min2 = np.zeros(M-1) for n in range(M+1): s[n] = s0*(np.power(u,M-n)*np.power(d,n)) min1[n] = s0*np.power(d,n) # minimum from upward path for i in range(M-1): min2[i] = s0*np.power(d,i) #minimum from down path s1 = np.zeros(M-1) for i in range(M-1): s1[i] = s[i+1] #option price at t = M excluding first and last #value of the option vm1 = np.zeros(M+1) vm2 = np.zeros(M-1) for i in range(M+1): vm1[i] = max(s[i]-min1[i],0) for i in range(M-1): vm2[i] = max(s1[i]-min2[i],0) #value at time t = M-1 x = np.zeros(M) for i in range(M): if i == 0: x[i] = np.exp(-r*dt)*(p*vm1[0]+(1-p)*vm2[0]) elif i == M-1: x[i] = np.exp(-r*dt)*(p*vm1[M-1]+(1-p)*vm1[M]) else: x[i] = np.exp(-r*dt)*(p*vm1[i]+(1-p)*vm2[i]) #value at time t = M-1 excluding first and last x1 = np.zeros([M-1,M]) x = np.vstack([x,x1]) for i in range(M-1): for j in range(M-1-i): x[i+1,j] = np.exp(-r*dt)*(p*x[i,j]+(1-p)*x[i,j+1]) print(x[M-1,0]) # result

import numpy as np s0 = 9 r = 0.06 sigma = 0.3 T = 3 M = 1000 dt = T/M u = np.exp(sigma*np.sqrt(dt)) d = 1/u a = np.exp(r*dt) p = (a-d)/(u-d) s = np.zeros(M+1) max1 = np.zeros(M+1) max2 = np.zeros(M-1) for n in range(M+1): s[n] = s0*(np.power(u,M-n)*np.power(d,n)) max1[n] = s0*np.power(u,M-n) for i in range(M-1): max2[i] = s0*np.power(u,i+(M-1)-(2*i+1)) s1 = np.zeros(M-1) for i in range(M-1): s1[i] = s[i+1] #payoff vmax1 = np.zeros(M+1) vmax2 = np.zeros(M-1) for i in range(M+1): vmax1[i] = max(max1[i]-s[i],0) for i in range(M-1): vmax2[i] = max(max2[i]-s1[i],0) #value at time = M-1 x = np.zeros(M) for i in range(M): if i == 0: x[i] = np.exp(-r*dt)*(p*vmax1[0]+(1-p)*vmax1[1]) elif i == M-1: x[i] = np.exp(-r*dt)*(p*vmax2[M-2]+(1-p)*vmax1[M]) else: x[i] = np.exp(-r*dt)*(p*vmax2[i-1]+(1-p)*vmax1[i+1]) #value untill t = 0 x1 = np.zeros([M-1,M]) x = np.vstack([x,x1]) for i in range(M-1): for j in range(M-1-i): x[i+1,j] = np.exp(-r*dt)*(p*x[i,j]+(1-p)*x[i,j+1]) #result print(x[M-1,0])

import numpy as np # input s0 = 9 r = 0.06 k = 10 sigma = 0.3 T = 3 M = 7 dt = T/M u = np.exp(sigma*np.sqrt(dt)) d = 1/u a = np.exp(r*dt) p = (a-d)/(u-d) # option price at time t (S(t)) s = np.zeros(M+1) max1 = np.zeros(M+1) max2 = np.zeros(M-1) for n in range(M+1): s[n] = s0*(np.power(u,M-n)*np.power(d,n)) max1[n] = s0*np.power(u,M-n) for i in range(M-1): max2[i] = s0*np.power(u,i+(M-1)-(2*i+1)) s1 = np.zeros(M-1) for i in range(M-1): s1[i] = s[i+1] #cari max vmax1 = np.zeros(M+1) vmax2 = np.zeros(M-1) for i in range(M+1): vmax1[i] = max(max1[i]-k,0) for i in range(M-1): vmax2[i] = max(max2[i]-k,0) #value at t = M-1 until t = 0 x = np.zeros(M) for i in range(M): if i == 0: x[i] = np.exp(-r*dt)*(p*vmax1[0]+(1-p)*vmax1[1]) elif i == M-1: x[i] = np.exp(-r*dt)*(p*vmax2[M-2]+(1-p)*vmax1[M]) else: x[i] = np.exp(-r*dt)*(p*vmax2[i-1]+(1-p)*vmax1[i+1]) x1 = np.zeros([M-1,M]) x = np.vstack([x,x1]) for i in range(M-1): for j in range(M-1-i): x[i+1,j] = np.exp(-r*dt)*(p*x[i,j]+(1-p)*x[i,j+1]) print(x[M-1,0])

import numpy as np s0 = 9 r = 0.06 k = 10 sigma = 0.3 T = 3 M = 7 dt = T/M u = np.exp(sigma*np.sqrt(dt)) d = 1/u a = np.exp(r*dt) p = (a-d)/(u-d) s = np.zeros(M+1) min1 = np.zeros(M+1) min2 = np.zeros(M-1) for n in range(M+1): s[n] = s0*(np.power(u,M-n)*np.power(d,n)) min1[n] = s0*np.power(d,n) for i in range(M-1): min2[i] = s0*np.power(d,i) s1 = np.zeros(M-1) for i in range(M-1): s1[i] = s[i+1] vm1 = np.zeros(M+1) vm2 = np.zeros(M-1) for i in range(M+1): vm1[i] = max(k-min1[i],0) for i in range(M-1): vm2[i] = max(k-min2[i],0) x = np.zeros(M) for i in range(M): if i == 0: x[i] = np.exp(-r*dt)*(p*vm1[0]+(1-p)*vm2[0]) elif i == M-1: x[i] = np.exp(-r*dt)*(p*vm1[M-1]+(1-p)*vm1[M]) else: x[i] = np.exp(-r*dt)*(p*vm1[i]+(1-p)*vm2[i]) x1 = np.zeros([M-1,M]) x = np.vstack([x,x1]) for i in range(M-1): for j in range(M-1-i): x[i+1,j] = np.exp(-r*dt)*(p*x[i,j]+(1-p)*x[i,j+1]) print(x[M-1,0])