前言:这篇笔记是根据姜禄彬老板在公众号上发布的笔记复刻的,不同的是原作者用的是股票数据,我用的是比特币期权数据。
这份笔记里主要是如何用python代码来计算BS模型、如何求隐含波动率、什么是波动率微笑、greeks等,整体还是有点乱,之后有时间再整理下。
另外由于公式用markdown写的,可能显示不出来,可以在nbviewer上查看,网址为:https://nbviewer.org/github/xuelixunhua/options_thing/blob/main/options.ipynb#。因为生成了很多图片,懒得上传,所以还是推荐大家去nbviewer上看。
在这篇文章中,我开始介绍第三种根据市场情绪得出的统计概率模型。这篇文章主要参考了斯文的《Python金融实战案例精粹》P362- P385,和斯文的《基于 Python 的金融分析与风险管理》P351-P362,吴清的《期权交易策略十讲》P68-P78,同时参考了我所写的证券研究系列连载第 145 篇《证券研究中的宏观对冲(5):股市波动性数学推导》, 第 147 篇《证券研究中的宏观对冲(7):股价波动的微分方程,对数正态分布与希腊值》的内容。
用 Python 分析期权交易的方法
看涨期权与看跌期权的模块:期权定价方法
斯文的《python 金融实战案例精粹》P362 到 367 页讲述了期权定价运用 Python 的编程方法。我们先来回顾一下,在系列连载第 144 篇中我们讲述了欧式看涨期权与看跌期权的数学公式。在期权投资当中已经知道了期权的价格是由股票价格 S0、行权价 K、波动率 σ \sigma σ、离到期的时间 T 和短期利率r 决定的。在实际的投资研究中影响期权最重要的变量是方向、波动率和时间,但是相对来说,我认为方向才是最重要的。金融学家推导出一个公式,从这些变量中计算出期权的价格。这个最基础的模型就是布莱克斯科尔斯模型( Black- Scholes Model)。这个模型是在 1973 年早期面世的,这个公式相当容易使用,公式很短,变量不多。c 为看涨期权定价公式,p 为看跌期权定价公式,实际的公式是:
KaTeX parse error: No such environment: align at position 8: \begin{̲a̲l̲i̲g̲n̲}̲ &c = S_0N(d_1)…
这些变量是:
- S0=基础资产在定价日的股票价格;
- K=行权价;
- T=期权合约的剩余期限(年),用 1 年的百分比表示;
- r=目前的无风险利率;
- σ \sigma σ=基础资产价格百分比变化(收益率)的年化波动率;
- ln=自然对数;
- N(x)= 标准正态分布的累积正态密度函数。
为了能够便利的计算欧式看涨期权与欧式看跌期权的定价,通过 Python 自定义一个运用 BSM 模型计算期权价的函数,并且在函数中能够输入包括基础资产价格(S)、期权执行价格(K)、波动率(sigma)、无风险利率(r)、期权合约定价日(即今天的交易日,T0)以及期权到期日(T1)和期权类型 Types(call, put)。
欧式看涨期权与欧式看跌期权的 Python 源代码如下,在以上自定义的函数 Value_BSM 中,输入基础资产价格、期权执行价格、波动率、无风险利率、定价日以及期权到期日等参数,就可以快速计算基于 BSM 模型的期权价格。
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from scipy.stats import norm
from pylab import mpl
import seaborn
seaborn.set()
mpl.rcParams['font.sans-serif'] = ['KaiTi']
mpl.rcParams['axes.unicode_minus'] = False
def value_BSM(S, K, sigma, r, T0, T1, types):
try:
S, K, types = float(S), float(K), types.tolist()[0]
except:
pass
T = (T1 - T0).apply(lambda x: x.days) / 365 # 计算期权的剩余期限
d1 = (np.log(S/K) + (r + pow(sigma, 2) / 2) * T) / sigma * np.sqrt(T)
d2 = d1 - sigma * np.sqrt(T)
if types == 'C':
value_call = S * norm.cdf(d1) - K * np.exp(-r*T) * norm.cdf(d2) # norm.cdf(d1) 是指 N(d1)
return value_call
if types == 'P':
value_put = K * np.exp(-r*T) * norm.cdf(-d2) - S * norm.cdf(-d1)
return value_put
比特币期权来计算期权价格
原文是计算腾讯股票的,当时使用的 6 个月期的 HIBOR 利率是 2.31393% 并作为定价中的无风险利率,这里直接沿用。
整理期权数据
该部分内容借鉴对股票数据的处理,先将所有的期权数据合并到一个 DataFrame,对时间和期权名称两列进行排序。
之后可以直接读取pkl。
import os
rule_type = '1H' # 根据需求修改时间,如5T、30T、1H、1D等
# 导入某文件夹下所有期权的代码
def get_option_code_list_in_one_dir(path):
"""
从指定文件夹下,导入所有csv文件的文件名
:param path:
:return:
"""
option_list = []
# 系统自带函数os.walk,用于遍历文件夹中的所有文件
for root, dirs, files in os.walk(path):
if files: # 当files不为空的时候
for f in files:
if f.endswith('.csv'):
option_list.append(f[:-4])
return sorted(option_list)
# ===读取所有期权代码的列表
path = r'C:\Users\xxx\okex\1H\BTC-USD'
option_code_list = get_option_code_list_in_one_dir(path)
# ===读取比特币价格
df_btc = pd.read_pickle(r'xxx\BTC-USDT.pkl')
df_btc = df_btc[['candle_begin_time', 'close']]
df_btc.columns = ['candle_begin_time', 'btc_price']
df_btc = df_btc.resample(rule_type, on='candle_begin_time', label='left', closed='left').agg({
'btc_price': 'last'})
df_btc['btc_return'] = np.log(df_btc['btc_price'] / df_btc['btc_price'].shift())
df_btc['btc_volatility'] = np.sqrt(365) * df_btc['btc_return'].rolling(365).std() # 计算年化波动率
# 循环读取期权数据
all_option_data = pd.DataFrame() # 用于存储数据
for code in option_code_list:
print(code)
# =读入数据
df = pd.read_csv(path + '/%s.csv' % code, encoding='gbk', skiprows=1, parse_dates=['candle_begin_time'])
agg_dict = {
'open': 'first',
'high': 'max',
'low': 'min',
'close': 'last',
'volume': 'sum',
'currency_volume': 'sum',
}
df = df.resample(rule_type, on='candle_begin_time', label='left', closed='left').agg(agg_dict)
df['symbol'], df['T1'], df['K'], df['types'] = code, pd.to_datetime('20' + code.split('-')[2]), code.split('-')[3], code.split('-')[-1]
df = pd.merge(left=df, right=df_btc, on='candle_begin_time', how='left', sort=True, indicator=True)
df.reset_index(inplace=True)
# 合并数据
all_option_data = all_option_data.append(df, ignore_index=False)
all_option_data.sort_values(['candle_begin_time', 'symbol'], inplace=True)
all_option_data.reset_index(inplace=True, drop=True)
all_option_data.to_pickle('option_btc_%s.pkl' % rule_type)
all_option_data
rule_type = '1D' # 根据需求修改时间,如5T、30T、1H、1D等
df = pd.read_pickle(r'C:\Users\xueli\data\coin\option\okex\option_btc_%s.pkl' % rule_type)
df_btc = pd.read_pickle(r'C:\Users\xueli\python_file\coin_quant\coin_alpha\data\backtest\output\pickle_data\spot\BTC-USDT.pkl')
df_btc = df_btc[['candle_begin_time', 'close']]
df_btc.columns = ['candle_begin_time', 'btc_price']
df_btc = df_btc.resample(rule_type, on='candle_begin_time').agg({
'btc_price': 'last'})
df_btc.reset_index(inplace=True)
计算收益率的历史年化波动率
在约翰∙赫尔《期权、期货及其他衍生产品》P256 页中列出了年化波动率的计算公式,如下:年波动率=每个交易日的波动率*\sqrt{每年的交易日天数}波动率*\sqrt{每年的交易日天数},每年的交易日天数为 365 天。这里采用的是滚动计算。
该部分内容已经在上述整理数据时解决。
df_btc['btc_return'] = np.log(df_btc['btc_price'] / df_btc['btc_price'].shift())
df_btc['btc_volatility'] = np.sqrt(365) * df_btc['btc_return'].rolling(365).std()
运用年化波动率计算某日看涨期权价格
我们随便找个日期,然后选定一个期权合约
df[df['candle_begin_time'] == pd.to_datetime('2021-07-02')]
option_thing = df[(df['symbol'] == 'BTC-USD-210702-35000-C') & (df['candle_begin_time'] == pd.to_datetime('2021-06-22'))]
call_btc = value_BSM(S=option_thing['btc_price'], K=option_thing['K'], sigma=option_thing['btc_volatility'], r=0.0231393, T0=(option_thing['candle_begin_time']), T1=(option_thing['T1']), types=option_thing['types'])
# 原始数据开高收低单位为BTC,通过乘以BTC价格获得当时期权价格
option_thing[['open', 'high', 'low', 'close']] = option_thing.loc[:, ('open', 'high', 'low', 'close')].apply(lambda x: x * option_thing['btc_price'])
模拟(看涨/看跌)期权到期日收益
针对上述的欧式看涨期权,模拟计算当日价格取值是区间[32000-38000]的等差数列所对应的期权价格数列:同时,模拟计算在期权合约到期日价格取值是区间[32000-38000]的等差数列所对应的期权合约收益情况(不考虑期权费并且是多头),将模拟得到的当前期权价格与到期日期权收益进行可视化(要求看涨期权)。假设K代表期权的执行价格,ST时基础资产在期权合约到期时的价格,在期权到期时且不考虑期权费的情况下期权收益(不考虑期权费)
from functools import partial
btc_price_list = np.linspace(32000, 38000, 1000) # 1000等分
retrun_call_list = np.maximum(btc_price_list - 35000, 0) # 看涨期权到期日收益序列,35000是上述的
value_call_list = []
for i in btc_price_list:
value_call_list.append(value_BSM(S=i, K=35000, sigma=option_thing['btc_volatility'], r=0.0231393, T0=option_thing['candle_begin_time'], T1=(option_thing['T1']), types='C').tolist()[0])
value_call_list
plt.figure(figsize=(12,10))
plt.plot(btc_price_list, value_call_list, 'r-', label=u'price of call option', lw=2.5)
plt.plot(btc_price_list, retrun_call_list, 'b-', label=u'return of call option', lw=2.5)
plt.legend(fontsize=13)
plt.show()
上述图不是真实世界的图,真实世界的看涨期权价格不会一直是纯线性的,而是带有曲线性质的。
期权价格都高于期权到期日的收益,这是因为期权价格包含了两部分:一部分是期权的内在价值( intrinsic value),它表示假定期权在定价日就能行权所获得收益,金额等于期权到期日的收益;另一部分是期权的时间价值( time value),它反映了期权在存续期内实现更高收益的可能性。同时,当期权处于平值时(即基础资产价格等于执行价格),期权的时间价值达到最大。
计算期权隐含波动率
隐含波动率:是市场参与者对未来波动率的预期值,当合约的隐含波动率显著高于其他合约时,则该合约被高估,反之则被低估。
随着期权执行价格的上升,隐含波动率不断下降。看涨期权执行价格越低,隐含波动率越高,看跌期权执行价格越低,隐含波动率越高,但这并不是绝对的。
计算隐含波动率的方法:二分法和 Newton Raphson 方法
计算隐含波动率的二分法
举一个简单的例子,假定初始猜测某支期权的波动率是 20%,对应该波动率数值估计得到的欧式看涨期权价格是 0.1035 元,显然,比市场价格 0.1566 元更小。由于期权价格是波动率的增函数,因此合理地估计正确的波动率应该会比 20%更大。然后假定波动率是 30%,对应的期权价格 0.336 元,这个结果又比 0.1566 元高,则可以肯定波动率是介于 20%-30%的区间中。接下来,取上两次波动率数值的均值,也就是波动率 25%,对应的期权价为 0.1662 元,这个值又比 0.1566 元高,但是合理的波动率所处的区间范围收窄至 20%与 25%之间,然后取均值 22.5%继续计算,每次迭代都使波动率所处的区间减半,最终就可以计算出满足较高精确度的隐含波动率近似值。
在系列连载第 145 篇中,我们讲到法博齐等人在《金融经济学》中提出了两种计算隐含波动率的方法,分别是二分法和 Newton Raphson 算法,这两种方法在 Timothy Sauer 的《数值分析》中有详细介绍。
首先举例说明二分法,假设某日微软的股票价格为\$80.8375,并且执行价格为\$85 的看涨期权售价\$2. 875。另外还有额外信息:S= \$80.375,X= \$85,T= 0.1945,r= 2.79%,d = 0.0。
对任何给定的 σ,我们就可以利用布莱克斯科尔斯公式来计算看涨期权的价格。尽管我们现在知道 C=$2.875,但是隐含波动率 σ 是未知的需要倒推才能得到。
为了找到 σ,我们可以从一个低σ值和一个高 σ 值开始。一个 σ= 10%大致说明,在下一年,股票价格有可能上涨 10%,也有可能下跌 10%,并且有可能是大多数股票波动率的低的估计值。另一方面,σ=100%就是一个高的估计值。计算 σ=10%时的看涨期权价格,标记为 C(σ=0.10),同样计算 σ= 100%时的看涨期权价格,我们有 C(σ=0.10) = \$0.2278,C(σ= 1.00) = \$12. 4018。
由于观测到的价格 C= \$2. 875 在 C(σ=0.10)和 C(σ=1.00)之间,由于看涨期权价格是波动率的增函数,波动率必须在 0.10 到 1.00 之间。
现在计算 0.10 和 1.00 中点的看涨期权价格,我们有 C(σ= 0.55) = \$6.0472,由于观测到的价格 C= \$2. 875 在C(σ=0.10)和 C(σ=0.55)之间,波动率就必须在 0. 10 和 0.55 之间。计算 0.10 和 0.55 中点的看涨期权价格,我们有 C(σ= 0.325) = \$2.9168。
由于观测到的价格 C= \$2. 875 在 C(σ=0.10)和 C(σ=0.325)之间,波动率就必须在 0.10 和 0.325 之间。计算 0. 10和 0.325 中点的看涨期权价格,我们有 C(σ= 0.2125)= \$1.4361。
由于观测到的价格 c= \$2. 875 在 C(σ=0.2125)和 C(σ= 0.325)之间,波动率就必须在 0.2125 和 0.325 之间。继续上面的过程,我们将发现
二分法第一种代码
因此,隐含波动率为 σ=0.3219=32.19%。可以非常容易地利用计算机对上述过程编程,从而迅速找到隐含波动率。上述方法称为二分法。总体来说,这里编程的思想是不断给与更高与更低的波动率数值,以找到所对应的期权所对应的波动率。
def implied_volatility(P, S, K, r, T0, T1, optype):
T = (T1 - T0).apply(lambda x: x.days) / 365 # 计算期权的剩余期限
sigma_min = 0.00001
sigma_max = 1.000
sigma_mid = (sigma_min + sigma_max) / 2
if optype == 'C':
def call_bs(S, K, sigma, r, T):
d1 = (np.log(S/K) + (r + pow(sigma, 2) / 2) * T) / (sigma * np.sqrt(T))
d2 = d1 - sigma * np.sqrt(T)
call_value = S * norm.cdf(d1) - K * np.exp(-r*T) * norm.cdf(d2)
return call_value
call_min = call_bs(S, K, sigma_min, r, T)
call_max = call_bs(S, K, sigma_max, r, T)
call_mid = call_bs(S, K, sigma_mid, r, T)
diff = P - call_mid
if P < call_min or P > call_max:
print('error, the price of option is beyond the limit')
while abs(diff) > 1e-6:
diff = P - call_bs(S, K, sigma_mid, r, T)
sigma_mid = (sigma_min + sigma_max) / 2
call_mid = call_bs(S, K, sigma_mid, r, T)
if P > call_mid:
sigma_min = sigma_mid
else:
sigma_max = sigma_mid
else:
def put_bs(S, K, sigma, r, T):
d1 = (np.log(S/K) + (r + pow(sigma, 2) / 2) * T) / (sigma * np.sqrt(T))
d2 = d1 - sigma * np.sqrt(T)
put_value = K * np.exp(-r*T) * norm.cdf(-d2) - S * norm.cdf(-d1)
return put_value
put_min = put_bs(S, K, sigma_min, r, T)
put_max = put_bs(S, K, sigma_max, r, T)
put_mid = put_bs(S, K, sigma_mid, r, T)
diff = P - put_mid
if P < put_min or P > put_max:
print('error, the price of option is beyond the limit')
while abs(diff) > 1e-6:
diff = P - put_bs(S, K, sigma_mid, r, T)
sigma_mid = (sigma_min + sigma_max) / 2
put_mid = put_bs(S, K, sigma_mid, r, T)
if P > put_mid:
sigma_min = sigma_mid
else:
sigma_max = sigma_mid
return sigma_mid
上面代码中,当误差 diff 小于 0.000001,返回 sigma_mid 的值,sigma_mid 的值即为我们所求的期权的隐含波动率的值。
二分法第二种代码
第二种方法是利用二分查找法并运用 Python 构建分别计算欧式看涨、看跌期权隐含波动率的自定义函数,具体的代码如下,这种方法是分别定义了看涨期权隐含波动率和看跌期权隐含波动率的函数算法,且少输入了一个参数
def impvol_call_binary(C, S, K, r, T):
try:
C, S, K = float(C), float(S), float(K)
except:
pass
def call_bs(S, K, sigma, r, T):
d1 = (np.log(S/K) + (r + pow(sigma, 2) / 2) * T) / (sigma * np.sqrt(T))
d2 = d1 - sigma * np.sqrt(T)
return float(S * norm.cdf(d1) - K * np.exp(-r * T) * norm.cdf(d2))
sigma_min = 0.001
sigma_max = 1.000
sigma_mid = (sigma_min + sigma_max) / 2
call_min = call_bs(S, K, sigma_min, r, T)
call_max = call_bs(S, K, sigma_max, r, T)
call_mid = call_bs(S, K, sigma_mid, r, T)
diff = C - call_mid
if (C < call_min) or (C > call_max):
print('error, the price of option is beyond the limit')
while abs(diff) > 0.001:
diff = C - call_bs(S, K, sigma_mid, r, T)
sigma_mid = (sigma_max + sigma_min) / 2
call_mid = call_bs(S, K, sigma_mid, r, T)
if C > call_mid:
sigma_min = sigma_mid
else:
sigma_max = sigma_mid
return sigma_mid
def impvol_put_binary(P, S, K, r, T):
try:
P, S, K = float(P), float(S), float(K)
except:
pass
def put_bs(S, K, sigma, r, T):
d1 = (np.log(S/K) + (r + pow(sigma, 2) / 2) * T) / (sigma * np.sqrt(T))
d2 = d1 - sigma * np.sqrt(T)
return float(K * np.exp(-r * T) * norm.cdf(-d2) - S * norm.cdf(-d1))
sigma_min = 0.001
sigma_max = 1.000
sigma_mid = (sigma_min + sigma_max) / 2
put_min = put_bs(S, K, sigma_min, r, T)
put_max = put_bs(S, K, sigma_max, r, T)
put_mid = put_bs(S, K, sigma_mid, r, T)
diff = P - put_mid
if P < put_min or P > put_max:
print('error, the price of option is beyond the limit')
while abs(diff) > 0.001:
diff = P - put_bs(S, K, sigma_mid, r, T)
sigma_mid = (sigma_max + sigma_min) / 2
put_mid = put_bs(S, K, sigma_mid, r, T)
if P > put_mid:
sigma_min = sigma_mid
else:
sigma_max = sigma_mid
return sigma_mid
依然沿用BTC期权作为分析对象
option_thing_2 = df[(df['symbol'] == 'BTC-USD-210702-35000-P') & (df['candle_begin_time'] == pd.to_datetime('2021-06-24'))]
# 原始数据开高收低单位为BTC,通过乘以BTC价格获得当时期权价格
option_thing_2[['open', 'high', 'low', 'close']] = option_thing_2.loc[:, ('open', 'high', 'low', 'close')].apply(lambda x: x * option_thing_2['btc_price'])
imp_vol_1 = impvol_call_binary(C=option_thing['close'], S=option_thing['btc_price'], K=option_thing['K'], r=0.04, T=(option_thing['T1'] - option_thing['candle_begin_time']).apply(lambda x: x.days) / 365)
imp_vol_2 = impvol_put_binary(P=option_thing_2['close'], S=option_thing_2['btc_price'], K=option_thing_2['K'], r=0.04, T=(option_thing_2['T1'] - option_thing_2['candle_begin_time']).apply(lambda x: x.days) / 365)
这里在计算中出现问题,有的数据会出现error的原因:波动率为0,看涨期权的价格下限为0.65左右,大于市场价0.6248,这个时候期权实际价格已经低于理论价格,BS公式超出其计算范围,用隐波已无法衡量实际价格的合理性了。期权客户端报价隐含波动率一栏有些显示“-”或者“0.0001”也是这个道理
此外,隐含波动率也正在被制作成指数。比如,芝加哥期权交易所(CBOE)于 1993 年对外发布了隐含波动率的指数 VIX( Volatility Index),VIX 也被称作"恐惧指数"(fear factor)。此后,德国、法国、英国、瑞士、韩国等国家及我国的香港、台湾地区也相继推出了波动率指数国际经验表明,波动率指数能够前瞻性地反映市场情绪与风险。例如,在 2008 年国际全融危机中,波动率指数及时准确地为全球各家金融监管机构提供了掌握市场压力、监控市场情绪的指标,有效提升了监管能力与决策水平。
2015 年 6 月 26 日,上海证券交易所发布了首只基于真实期权交易数据编制的波动率数中国波指(iVIX)中国波指是用于衡量上证 50FTF 基金未来 30 日的预期波动,它的推出一方面为市场提供了高效灵敏的风险监测指标,便于市场实时衡量市场风险,增强技术分析手段,提升策略交易能力;另一方面也对进一步丰富上海证券交易所的衍生产品种类、实现行生品发展战略产生了积极意义。
计算隐含波动率的牛顿迭代法
牛顿迭代法的数学推导
简单来说,这一算法是找到一个方程的根,即从任意初始值x0开始,根据下式不断迭代,获得 f(x)=0 的解。
牛顿方法对应的几何图如下图所示。为了找到函数 f(x)=0 的根,给定一个初始估计x0,画出函数 f 在x0点的切线。用切线来近似函数 f,求出其与 x 轴的交点作为函数 f 的根,但是由于函数的弯曲,该交点可能并不是精确解。因而,该步骤要迭代进行。
从x0开始,画出函数 y=f(x)的切线和 x 轴的交点记做x1,这是对于函数根的下一个近似。并且可以看到,再将x1代入 f(x)中求其切线,就可以得到一个更近的数值。从几何图像中我们可以推出牛顿方法的公式,x0点的切线斜率可由导数 f(x)给出,并且切线上的一点是(x0,f(x0))。
我们先回顾一下点斜式公式,设 Q(x,y)是直线 l 上不同于点 P 的任一点,由于点 P,Q 都在 l 上,所以,可以用点 P(x0,y0)的坐标来表示直线 l 的斜率,于是可得以下方程: y − y + 0 = k ( x − x 0 ) y-y+_0=k(x-x_0) y−y+0=k(x−x0)
这就是所求的过点 P(x0,y0),斜率为 k 的直线 l 的方程,这个方程是由直线上的一点和斜率(一个方向)所确定的,称为直线方程的点斜式。
一般地,如果一条直线 l 上任意点的坐标(x,y)都满足一个方程,满足该方程的每一个数对(x,y)所确定的点都在直线 l 上,我们就把这个方程称为直线 l 的方程。如果已知直线 l 上一点 P(x0,y0)及斜率 k,可用上述方法求出直线 l 的方程。
我们来看一下 Timothy Sauer 所写的《数值分析》,由于在 Newton Raphson 中,一条直线的点斜率方程是y − f(x0) = f’(x0)(x − x0),因而切线和 x 轴的交点等价于在直线中令 y=0:
KaTeX parse error: No such environment: align at position 8: \begin{̲a̲l̲i̲g̲n̲}̲ &f'(x_0)-(x-x_…
求解 x 得到根的近似,称为x1,然后,重复整个过程,从x1开始,得到x2,等等,进而得到如下的迭代公式:牛顿方法是:x0=初始估计,xi+1 = xi −f(xi) / f’(xi),其中 i=0,1,2,3,。。。,直到收敛,具体可以见上图。
下面将牛顿迭代方法应用到期权定价中。构造函数f(σ) = C(σ = σi) − C,其中C是当前已知即时期权价格,那么波动率就是对f(σ) =0 求解,并且可以知道隐含波动率的数值是唯一的。
首先需要算出f(σ)的导数 f ′ ( σ ) = e d 1 2 / 2 2 π ( S T ) f'(\sigma)=e^{d_1^2/2}\sqrt{2\pi}(S\sqrt T) f′(σ)=ed12/22π(ST