為什麼我的 CRR 模型實現不收斂？

回想一下，期權定價的 CRR（Cox-Ross-Rubinstein）模型是通常的二叉樹模型 $ u $ （上因素）和 $ p $ （風險中性機率之一）定義如下： $$ u = e^{\sigma\sqrt{\Delta t}}, $$ $$ p = \frac{e^{r\Delta t} - e^{-\sigma\sqrt{\Delta t}}}{e^{\sigma\sqrt{\Delta t}} - e^{-\sigma\sqrt{\Delta t}}}, $$ 在哪裡 $ \sigma $ 是波動率， $ r $ 是利率， $ \Delta t = \frac{T}{M} $ （時間步長），然後 $ d = \frac{1}{u} $ （下降因素）和 $ q = 1- p $ ，其他一切照常。

下面是我使用 Python 3 實現的 CRR 模型：

# Implementation of Cox-Ross-Rubenstein option's pricing model.

import math

T = 0.25 # time horizon
M = 2 # quantity of steps
t = T/M # step
sigma = 0.1391*math.sqrt(0.25) # volatility
r = 0.0214*0.25 # interest rate
u = math.exp(sigma*math.sqrt(t)) # up-factor
d = 1.0/u # down-factor
S0 = 2890.30 # initial underlying stock price
K = 2850 # strike

# compute risk-neutral probabilities
p = (math.exp(r*t)-math.exp(-sigma*math.sqrt(t)))/(math.exp(sigma*math.sqrt(t))-math.exp(-sigma*math.sqrt(t))) # up
q = 1 - p # down

# profit from call option
def call(stock_price, K):
   price = max(stock_price - K, 0)
   return price

# profit from put option
def put(stock_price, K):
   price = max(K - stock_price, 0)
   return price

# price for European style
def european():
   price = 1.0/(1+r)*(p*option_prices[i+1][j+1]+q*option_prices[i+1][j])
   return price

# price for American style, specify call or put in argument
def american(style):
   price = max(style, european())
   return price

stock_final_prices = []
option_final_prices = []

# create dictionary, containing lists of options prices at every time step
option_prices = {}
for i in range(0,M+1):
   option_prices[i] = [None] * (i + 1)

# calculate possible final stock prices
for i in range(0,M+1):
   stock_final_prices.append(S0*math.pow(u,i)*math.pow(d,M-i))

# calculate possible option final prices -- choose call or put function
for i in range(0,M+1):
   option_final_prices.append(put(stock_final_prices[i], K))

option_prices[M] = option_final_prices

# going backwards -- uncomment european or american function, choose call or put for american style
for i in range(M-1,-1,-1):
   for j in range(0,i+1):
       option_prices[i][j] = european()
       #option_prices[i][j] = american(call(S0*math.pow(u,j)*math.pow(d,i-j), K))

print('The price is ${0} for {1} steps.'.format(option_prices[0][0], M))

您可以嘗試使用時間步數變數 $ M $ 看看什麼時候 $ M $ 增長期權的價格變為零，這沒有任何意義。但是，如果您指定的數值 $ p $ 和 $ u $ 手動使用 $ M = T $ ，它將成為完美執行的通常的二叉樹模型（Black-Scholes-Merton）。

那麼，為什麼我的 CRR 實現沒有收斂到一些有意義的非零價格？我在哪裡做錯了？我真的卡住了，找不到它。任何有關程式碼審查的幫助將不勝感激。

1. 從外觀上看，您的貼現是不正確的，因為當您增加 M 時，您應該貼現 1/(1+r0t)（假設 r0=0.0214 是您似乎貼現 1/(1+r0) 的年利率T))
2. 似乎您的“r”變數中有錯字，可能您希望它是 r=r0=0.0214 而不是 r=0.0214*0.25。

現在對樹程式碼的工作重寫是

# Implementation of Cox-Ross-Rubenstein option's pricing model.

from scipy.stats import norm
import numpy as np

#OP inputs as i understand them
T = 0.25 # time horizon
M = 2 # quantity of steps
sigma = 0.1391*np.sqrt(0.25) # volatility
r0 = 0.0214
S0 = 2890.30 # initial underlying stock price
K = 2850 # strike

#size M+1 grid of stock prices simulated at time T
def stock_prices(S0,T,sigma,M):    
   res = np.zeros(M+1)    
   t = T*1.0/M # step
   u = np.exp(sigma*np.sqrt(t)) # up-factor
   d = 1.0/u
   dn = d/u
   res[0] = S0*np.power(u,M) 
   for i in range(1,M+1):
       res[i] = res[i-1] * dn
   return res

# terminal payoff from call option
def payoff(stock_price, K,kind='call'):
   epsilon = 1.0 if kind == 'call' else -1.0
   price = np.maximum(epsilon*(stock_price - K), 0)
   return price

# price for European style option using CRR
def european_crr(S0,K,T,r,sigma,M,kind='call'):
   #terminal payoff
   option_price = payoff(stock_prices(S0,T,sigma,M),K,kind)
   t = T*1.0/M # time_step
   df = np.exp(-r*t) #discount factor
   u = np.exp(sigma * np.sqrt(t))
   d = 1/u
   p = (np.exp(r*t)-d)/(u-d) # risk neutral probability for up-move        
   q=1-p    
   for time_idx in range(M): #move backward in time
       for j in range(M-time_idx):
           option_price[j] = df*(p*option_price[j]+q*option_price[j+1])
   return option_price[0]

#analytical check Black-Scholes formula (no dividend nor repo)
def european_bs(S0,K,T,r,sigma,kind='call'):
   df = np.exp(-r*T)
   F = np.exp(r*T)*S0                   # forward price with no dividend or repo
   m = np.log(F/K)/(sigma * np.sqrt(T)) #moneyness
   epsilon = 1.0 if kind == 'call' else -1.0
   d1 = m + 0.5*sigma*np.sqrt(T)
   d2 = d1 - sigma * np.sqrt(T)
   Nd1 = norm.cdf(epsilon*d1)
   Nd2 = norm.cdf(epsilon*d2)
   return epsilon*df*(F*Nd1-K*Nd2)

• 您可以檢查您的輸入，使用您提供的輸入，分析 (european_bs) 和樹 (european_crr) 方法收斂到看跌期權的 17.97 和看漲期權的 73.48
• 另請注意，對於樹，不需要保存二維數組。使用大小為 M+1 的單個數組來保存臨時結果就足夠了
• 當我對 1+2 中解釋的程式碼執行兩項更改時，您的結果與我的樹和分析版本的結果匹配

引用自：https://quant.stackexchange.com/questions/42157