Diferencia entre el valor de la opción binomial de forma cerrada y la simulación monte carlo

Question

Estoy tratando de calcular el precio de una opción de compra europea utilizando tanto la expresión de forma cerrada como una simulación monte carlo. Pero los valores que obtengo de ambos métodos no son los mismos:

Expresión de forma cerrada:

$$q = \frac{(1+r)-d}{u-d}$$

$$C \frac{1}{(1+r)^T} * \left [\sum \limits_{i=0}^T \binom{T}{i}*q^i*(1-q)^{T-i}*max(u^i*d^{T-i}*S_0-K, 0) \right ]$$

Implementación en Python de la expresión cerrada desde:

import math
T = 10 # Number of periods
S0 = 8 # Starting price of stock
K = 9 # Strike price of option
r = 0.2 # Risk free interest rate
u = 1.5 # Up factor
d = 0.5 # Down factor

C = 0 #Value of call

risk_free = 1 / (1 + r)**T

q = ((1 + r) - d) / (u - d)

for i in range(T+1):
    prob = math.comb(T, i)*(q**T)*(1-q)**(T-i)
    ST = max(((u**i)*(d**(T-i))*S0)-K, 0) 
    C += ST*prob

print(risk_free*C)

Salida: 4.945275514422904

Implementación en Python de la simulación monte carlo:

import random
T = 10 # Number of periods
S0 = 8 # Starting price of stock
K = 9 # Strike price of option
r = 0.2 # Risk free interest rate
u = 1.5 # Up factor
d = 0.5 # Down factor

n = 20000 # Number of runs
for j in range(n):
    S = S0
    for i in range(T):
        S *= u if random() < q else d
    value += max(S - K, 0)
value /= n * (1 + r) ** T
print("For {} runs the value is {}".format(n, value))

Salida: 6.876698097695621

No entiendo cuál es la causa de esta diferencia, porque el código produce los mismos valores cuando pongo T=2 and S0=10 pero esa entrada tiene un valor p diferente de 0.2 mientras que la entrada actual tiene un valor p de 0.25, pero no entiendo que significa el valor p ya que no se usa en la fórmula.

Answer 1

Bien, he encontrado el problema, mi implementación de la fijación de precios binomial estaba mal.

Esta implementación de Python:

T = 10 # Number of periods
S0 = 8 # Starting price of stock
K = 9 # Strike price of option
r = 0.2 # Risk free interest rate
u = 1.5 # Up factor
d = 0.5 # Down factor

C = 0

q = ((1+r) - d) / (u - d)
risk_free = 1 / ((1 + r)**T)

for i in range(0, T+1):
    prob = math.comb(T, i) * (q**i) * (1-q)**(T-i)
    ST = (u**i) * (d**(T-i)) * S0
    max_value = max(ST - K, 0)
    C += max_value * prob

print(C * risk_free)

Salidas: 6.836045774062984

Que está mucho más cerca de la salida del MC