0. pre

在《给你的二叉树期权定价》中就挖了坑要写期权定价的代码，这会有时间来填坑啦。
本文将会用python实现欧式期权定价。具体的定价算法分别是基于BS公式的、蒙特卡洛的以及二叉树的。
对于二叉树和BS公式还不熟悉的小伙伴可以移步至往期关于二叉树期权定价和BS公式的3篇文章先熟悉一下二叉树和BS公式本身，然后再来看代码实现。当然，如果您对于理论推导部分不感兴趣也可以考虑跳过该部分直接看代码。只是这样可能理解起来稍微麻烦些。
接下来我们将先描述基于蒙特卡洛的期权定价方法，然后再分别对三种方法进行python代码的实现。

1.基于蒙特卡洛的期权定价

1.1 蒙特卡洛思想

蒙特卡洛是基于随机抽样的一种统计模拟方法。由于大数定律的成立，此方法可以用于计算一些难搞的期望的数值解或者是进行数值积分的计算。
由于在风险中性的假设下，欧式期权价格等于到期时期权价值的期望在无风险利率下的贴现。因此，这里的问题就成了如何计算时刻期权价值的期望。
由于标的资产服从的既定分布，蒙特卡洛模拟通过生成随机数的方式模拟资产价格的若干条路径，并以此得到这若干条价格路径在到期时的价格。
又由于看涨和看跌期权到期时的价值分别为和。因此我们便可以通过到期时标的资产价格计算期权价值，再由大数定律知样本量足够大时均值趋于期望而计算出时刻期权价值的期望。最后用无风险利率将该期望贴现到0时刻便得到期权价格。

1.2 蒙特卡洛期权定价的数学语言

(1) 标的资产价格服从对数正态分布

首先，资产价格服从几何布朗运动：又有伊藤引理：因此，对于来讲有：所以：

(2) 蒙特卡罗模拟步骤

1、生成若干个标准正态分布随机数。
2、将这n个标准正态分布随机数带入方程得到n个到期时的标的资产价格。
3、求这n个标的资产价格的均值得。
4、用无风险利率对贴现得0时刻的期权价格。

2. 参数初始化

import numpy as npfrom numpy.random import normalfrom scipy.stats import normclass european_option_pricing:"""
    Desc: 此类用于实现各种方式的欧式期权定价
    """def __init__(self, r, sigma, T, K, s0, call):"""
        :param r: 无风险利率
        :param sigma: 标的资产波动率
        :param T: 期权期限
        :param K: 行权价
        :param s0: 标的资产当前价格
        :param call: 是否是看涨期权
        """
        self.r = r
        self.sigma = sigma
        self.T = T
        self.K = K
        self.s0 = s0
        self.call = call

3.基于BS公式的代码实现

基于BS公式的期权定价相对简单，只需要将5个参数填入一下两个看涨和看跌的定价公式。注意：下面的代码都是类european_option_pricing下的方法，每一种方法分别对应一种期权定价算法。

def bs_formula(self):"""
    Desc: 由BS公式计算欧式期权价格
    """
    d1 = (np.log(self.s0/self.K) + (self.r + self.sigma**2 / 2) * self.T) / (self.sigma * np.sqrt(self.T))
    d2 = (np.log(self.s0/self.K) + (self.r - self.sigma**2 / 2) * self.T) / (self.sigma * np.sqrt(self.T))if self.call:
        f0 = self.s0 * norm.cdf(d1) - \
             self.K * np.exp(-self.r * self.T) * norm.cdf(d2)                      # c = s0*N(d1) - Ke^{-rT} * N(d2)else:
        f0 = self.K * np.exp(-self.r * self.T) * norm.cdf(-d2) - \
             self.s0 * norm.cdf(-d1)                                           # p = Ke^{-rT} * N(-d2) - s0 * N(-d1)return f0

4.基于蒙特卡洛的代码实现

def monte_carlo(self, log, path_num, step_num):"""
    Desc: 风险中性假设下利用蒙特卡洛方法计算欧式期权价格
    :param log: 是否根据标的资产对数来模拟路径。
        True: 假设标的资产价格服从对数正态分布[logSt ~ N(logS0 + (r-sigma^2/2)t, sigma^2 *t)]并以此计算sT。
        False: 假设标的资产价格服从几何布朗运动[ds = r*s*dt + sigma*s*dz]并以此计算sT。
    :param path_num: 蒙特卡洛生成的路径数量
    :param step_num: 每条价格路径的步数
    """# (1) 风险中性下模拟标的资产价格路径并计算到期时各路径下的价格if log:# 生成1000条价格路径
        log_sT = np.log(self.s0) + (self.r - self.sigma**2/2) * self.T + \
                 self.sigma * normal(size=path_num) * np.sqrt(self.T)
        sT = np.exp(log_sT)else:# 生成path_num条价格路径，每条路径走1000步
        sT = []for i in range(path_num):
            simu_ret = self.r * self.T/step_num + self.sigma * normal(size=step_num) * np.sqrt(self.T/step_num)     # ds/s = r * dt + sigma * dz
            cum_ret = (1 + simu_ret).prod()
            sT.append(self.s0 * cum_ret)# (2) 基于模拟的到期日标的资产价格计算期权价格if self.call:
        fT = np.array([max(0, s - self.K) for s in sT])                   # 看涨期权到期时各路径下的价格else:
        fT = np.array([max(0, self.K - s) for s in sT])                   # 看跌期权到期时各路径下的价格
    f0 = np.mean(fT * np.exp(-self.r * self.T))return f0

5.基于二叉树的代码实现

def binary_tree(self, step_num):"""
    Desc: 利用二叉树，从期权到期时的叶节点倒推地计算前面各个节点的期权价值，直到0时刻即得出期权价格。
    Note: 欧式期权只需要倒推地计算每一层各个节点期权的隐含价值但美式期权还要比较每个节点期权隐含价值与直接行权价值
          的相对大小，并最终取价值大的那一个作为该节点的期权价值。
    :param step_num:  二叉树步数
    """# (1) 计算 u, d, p, q and dt
    dt = self.T/step_num
    u = np.exp(self.sigma * np.sqrt(dt))
    d = np.exp(-self.sigma * np.sqrt(dt))
    p = (np.exp(self.r * dt) - d)/(u - d)
    q = 1 - p# (2) 构造二叉树价格路径的矩阵(n步二叉树需要一个n+1*n+1矩阵来装价格路径)
    s = np.array([[0] * (1 + step_num)] * (1 + step_num), dtype=float)for i in range(0, 1 + step_num):                           # 遍历每一步for j in range(i+1):                                   # 遍历每一步中标的资产价格的每种可能
            s[j, i] = self.s0 * u**(i-j) * d**j# (3) 构造二叉树下期权内涵价值路径的矩阵
    f = np.array([[0] * (1 + step_num)] * (1 + step_num), dtype=float)for i in range(step_num, -1, -1):# 倒推地计算期权价值if i == step_num:# 最后一步的期权价格if self.call:
                f[:, i] = [max(sT - self.K, 0) for sT in s[:, i]]else:
                f[:, i] = [max(self.K - sT, 0) for sT in s[:, i]]else:# 其他时候的期权价格(因为基于风险中性，所以第t时刻的期权价值等于t+1时刻期权价值的贴现)for j in range(i+1):
                f[j, i] = np.exp(-self.r * dt) * (p * f[j, i+1] + q * f[j+1, i+1])
    f0 = f[0, 0]return f0

6.测试比较

(1)首先，计算BS公式下的期权价格。

europe_opt_price = european_option_pricing(r=0.05, sigma=0.2, T=1, K=10, s0=10, call=True)
bs_price = europe_opt_price.bs_formula()
print('BS公式下的期权价格：', bs_price)
BS公式下的期权价格： 1.045058357218557

(2)接着，计算蒙特卡洛下的期权价格。

for num in [10, 100, 1000, 10000, 100000, 100000]:
    monte_carlo_price = europe_opt_price.monte_carlo(log=True, path_num=num, step_num=100)
    print('{}条路径下的蒙特卡洛下的期权价格：'.format(num), monte_carlo_price)10条路径下的蒙特卡洛期权定价： 1.28783643643895100条路径下的蒙特卡洛期权定价： 0.87416519605453771000条路径下的蒙特卡洛期权定价： 0.987155982802952410000条路径下的蒙特卡洛期权定价： 1.062215087666785100000条路径下的蒙特卡洛期权定价： 1.0439738276791009500000条路径下的蒙特卡洛期权定价： 1.0445916779825162

(3)最后，计算二叉树下的期权价格。

for num in [10, 50, 100, 500, 1000]:
    binary_tree_price = europe_opt_price.binary_tree(step_num=num)
    print('{}步二叉树下的期权价格：'.format(num), binary_tree_price)10步二叉树下的期权价格： 1.025340904487193750步二叉树下的期权价格： 1.041069154073265100步二叉树下的期权价格： 1.043061166224914500步二叉树下的期权价格： 1.04465851364465381000步二叉树下的期权价格： 1.0448584103764602

从上面三种方法的计算结果可知BS公式、二叉树和蒙特卡洛模拟相互间是一致的。当二叉树步数和蒙特卡洛路径数量越来越大时他们的定价结果与BS公式越接近。而这一点也是与理论预期相符的。
话说欧式搞定了，美式还会远么

坑继续挖着，心情好的时候又来填
参考文献：《期权期货及其他衍生品》 by John Hull