





摘要: 首先, 用原始對偶活躍集方法求解期權(quán)定價正問題, 將相應的數(shù)值解作為監(jiān)督學習的輸出, 然后用訓練好的神經(jīng)網(wǎng)絡替代期權(quán)定價正問題模型. 其次, 結(jié)合Bayes推斷與神經(jīng)網(wǎng)絡進行Metropolis-Hastings采樣, 求解隱含波動率反問題. 該方法減少了采樣過程中正問題計算量龐大的問題, 從而加速了反問題求解過程.
關(guān)鍵詞: 隱含波動率; Bayes推斷; 神經(jīng)網(wǎng)絡; 替代模型
中圖分類號: O241.8""文獻標志碼: A""文章編號: 1671-5489(2024)06-1363-07
Solving "Implied Volatility of American Lookback Options byBayesian Inference and Neural Network
TAO Li1, ZHU Benxi2, QIAN Yiyuan2, XU Jiaqi2
(1. College of International Business, Hainan University, Haikou 570228, China;
2. College of Mathematics, Jilin University, Changchun 130012, China)
Abstract: Firstly, "we used the original dual active set method to solve the forward problem of option pricing, with the corresponding numerical
solutions as the output for supervised learning, and then replaced the "forward problem model of option pricing with a well-trained neural network. Secondly, we combined
Bayesian inference with neural networks for Metropolis-Hastings sampling "to solve the inverse problem of implied volatility. This method reduced the problem of large computational
complexity "of the forward problem during the sampling process, thereby accelerating the solution process for the inverse problem.
Keywords: implied volatility; Bayesian inference; neural network; surrogated model
隱含波動率是期權(quán)的一個重要指標, 它反映了市場對于標的價格的變動幅度. 投資者可以利用它來預測標的價格的未來走勢和供求關(guān)系. 求解期權(quán)的隱含波動率, 可視為解決期權(quán)
定價模型的反問題.Achdou等[1]使用最小二乘公式和下降算法實現(xiàn)了對波動率的求解; Burkovska等[2]用縮基法代替經(jīng)典的有限元法, 從而減少了求解的復雜度; Kutner[3]利用基于二次近似法的數(shù)值方法, 計算了美式期權(quán)的隱含波動率. 傳統(tǒng)數(shù)值方法雖然能保證求解精度, 但當模型為非線性或高維時, 通常需要耗費大量計算時間, 計算效率較低. 近年來, 深度學習的快速發(fā)展提供了解決高維問題的新思路. Han等[4]采用了一種基于深度學習的方法, 處理一般的高維拋物型偏微分方程; Liu等[5]建立了一個校準神經(jīng)網(wǎng)絡框架, 計算了隱含波動率. 上述基于深度學習求解高維偏微分方程的方法, 在保證精度的情況下, 顯著減少了計算所需時間. 本文提出一個利用神經(jīng)網(wǎng)絡替代正問題的模型, 并結(jié)合Bayes推斷的數(shù)值方法, 快速求解隱含波動率.
1"正問題描述和數(shù)值算法
設V是期權(quán)價格, S是當前股票價格, J是股票價格在過去到當前時間的最大值, t是到期日與當前時刻的差值, σ是波動率, r是無風險利率, q
是股息收益率. 美式回望看跌期權(quán)價格V=V(S,J,t;σ)滿足如下自由邊界問題[6]:
LV∶=Vt+12σ2S
22VS2+(r-q)SVS-rV=0,"B(t)lt;S≤Jlt;+∞,"0≤tlt;T,(1)
其中B(t)為未知自由邊界, 收益函數(shù)為F(S,J,T)=J-S且J=max0≤ζ≤t Sζ, 波動率函數(shù)為
σ(t)=∑di=1σiai(t).
這里, {ai(t)}di=1是L2([0,T])有限維逼近空間的一組已知基底(例如{sin(it)}d-1i=0, {eit}d-1i=0), {σi}di=1為參數(shù). 則相應的線性互補模型為
LV·(V-G)=0,JX≤S≤Jlt;+∞,"0≤tlt;T,
LV≤0,"V≥G,JX≤S≤Jlt;+∞,"0≤tlt;T,F(xiàn)(S,J,T)=(J-S),JX≤S≤Jlt;+∞,(2)
其中G=J-S是贖回價格, X滿足:
X=-r+q+σ2/2-(-r+q+σ2/2)2+2rσ2-r+q-σ2/2-(-r+q+σ2/2)2+2rσ2.
利用變換x=ln JS, τ=T-t和u(x,τ)=V(S,J,t)S, 線性互補模型(2)可改寫為
(κ(x)uτ-(γκ(x)ux)x+qκ(x)u)·(u-g)=0,0≤ xlt;L,"0lt;τ≤T,
κ(x)uτ-(γκ(x)ux)x+qκ(x)u≥0,"u-g≥0,0≤xlt;L,"0lt;τ≤T,u(x,0)=g(x),0≤x≤L,
ux(0,τ)=0,"u(L,τ)=g(L),0lt;τ≤T,(3)
其中γ=σ22, κ(x)=e(q-r-γ)x/γ, g(x)=ex-1, L=ln 1X.
設H1L(Ω)∶={v∈H1(Ω); v≥g(x), v(L)=g(L)},H1E0(Ω)∶={v∈H1(Ω); v≥g(x), v(L)=0},""Ω=[0,L].
引理1"設u∈L2(0,T;H2(Ω)), uτ∈L2(0,T;L2(Ω)), 且u(0,τ)gt;g(0), 0lt;τ≤T. 則u是模型(3)的解當且僅當
u是下列變分不等式問題的解:
(κ(x)uτ,v-u)+(γκ(x)ux,vx-ux)+(qκ(x)u,v-u)≥0,"v∈H1E0(Ω),"0lt;τ≤T,
且u(x,0)=g(x).
令時間和空間的剖分分別為
Tt: 0=τ0lt;t1lt;…lt;tM=T,"Δt=TM,"tj=jΔt,"j=0,1,…,M,
Th: 0=x0lt;x1lt;…lt;xN=L,"h=LN,"xi=ih,"i=0,1,…,N,
其中M和N為正整數(shù). 設Vh(Ω)H1(Ω)為分片線性函數(shù)空間, 取有限元試探函數(shù)空間和檢驗函數(shù)空間分別為
VhL(Ω)∶={vh∈Vh(Ω); v(xn)≥g(xn), n=0,1,…,N-1, vh(L)=g(L)}=span{φi(x)}Ni=1,
VhE0(Ω)∶={vh∈Vh(Ω); v(xn)≥g(xn), n=0,1,…,N-1, vh(L)=0}.
則半離散格式為: 求uh(τ)∈VhL(Ω), 使得uh(x,0)=gI(x), 且
(κ(x)uhτ,vh-uh)+(γκ(x)uhx,vhx-uhx)+(qκ(x)uh,vh-uh)≥0,"vh∈VhE0(Ω),"0lt;τ≤T,
其中g(shù)I(x)為g(x)在VhL(Ω)中的插值.
進一步, 在時間方向采取向后Euler格式, 則全離散格式為
(κ(x)((1+qΔτ)umh-um-1h),vh-umh)+Δτ(γκ(x)umhx,vhx-umhx)≥0,"vh∈VhE0(Ω),"m=1,2,…,M,
其相應的矩陣形式為
(DUm,V-Um)≥(Wm,V-Um),
V≥F1,"m=1,2,…,M,(4)
其中D=γΔτA+(1+qΔτ)B, Wm=BUm-1-F2,
V=(v0,v1,2,…,vN-1)T, Um=(um0,um1,…,umN-1)T, A=(Ai,j)N×N,
B=(Bi,j)N×N, i,j=0,1,…,N-1, F1=(g(x0),…g(xN-1))T, F2=(0,…,0,g(L)Δτ(qBN-1,N+γAN-1,N))T,
Ai,j=(κ(x)φ′j,φ′i), Bi,j=(κ(x)φj,φi), i,j=0,1,…,N-1.
算法1[7-8]"原始對偶活躍集方法.
輸入: U0,F(xiàn)1,Wm;
輸出: Um,Xml,Xmr.
for m=1∶M, do
k=0, ρgt;0, εgt;0, U(0)=max{Um-1,F(xiàn)1}, U(1)=0, λ(0)≥0.
while ‖U(k+1)-U(k)‖∞lt;ε, do
① 令活躍集為IS(k+1), 非活躍集為AS(k+1):
IS(k+1)={i∈N, λ(k)i+ρ(Gmi-U(k)i)≤0},
AS(k+1)={i∈N, λ(k)i+ρ(Gmi-U(k)i)gt;0},
N={0,1,…,N-1}.
② 求解下列方程, 得到U(k+1)和λ(k+1):
DU(k+1)-λ(k+1)=Wm,
U(k+1)=Gmi, i∈AS(k+1),
λ(k+1)i=0, i∈IS(k+1).
end while
Um=U(k+1),
lm=min(IS(k+1)), rm=max(IS(k+1)),
xml=xlm, xmr=xrm.
end for
本文采用算法1對式(4)進行數(shù)值求解. 為方便, 上述正問題整體求解過程, 可以抽象為如下算子:
V=T(σ),(5)
其中σ=(σ1,…,σd)T為波動率參數(shù)向量, V=(U0,…,UM)為期權(quán)離散矩陣, T為數(shù)值求解算子.
2"反問題描述和Bayes算法
在式(5)中, 如果已知期權(quán)離散矩陣V和G, 求波動率參數(shù)向量σ, 則稱為反問題. 該類問題是有確定性的, 通常也是病態(tài)的,
一般用Tikhonov正則化方法[9]求解該類問題. 但從統(tǒng)計學的角度, 更關(guān)注的是真解的分布, 而不再是一個點. 式(5)的統(tǒng)計模型版本可寫為
V=T(σ)+δ,(6)
其中δ={δij}j=1,2,…,M+1i=1,2,…,N, δij~δ·N(0,1), N(0,1)是均值為0、 方差為1的標準正態(tài)分布.
2.1"Bayes公式
Bayes公式為
πV(σ)∶=π(σV)=
π(Vσ)π0(σ)∫Zπ(Vσ)π0(σ)dσ,(7)
其中πV(σ)和π0(σ)分別是波動率后驗和先驗的概率密度函數(shù). 于是有
πV(σ)=πη(V-T
(σ))π0(σ)∫Zπη(V-T(σ)π0(σ)
dσ∝πδ(V-T(σ))π0(σ),(8)
其中∝表示兩者成比例的關(guān)系, πδ(·)為最大似然估計, 并通常考慮先驗為正態(tài)分布.
在線性正演模型和高斯先驗的假設下, 最大似然估計等價于Tikhonov正則化[10]. 從確定性的角度可見, 通過優(yōu)化框架只能得到反問題相關(guān)的一個點估計.
而從統(tǒng)計學的角度, 可根據(jù)Bayes公式獲得更多的統(tǒng)計信息. 而且這個后驗分布, 使得反演層參數(shù)的刻畫更全面并包含不確定性[11]
. 本文采用Metropolis-Hastings(MH)算法, 從后驗分布πV(σ)中進行采樣.
2.2"Metropolis-Hastings算法
Markov Chain Monte Carlo (MCMC)算法是一類用于從復雜分布中抽樣的統(tǒng)計方法, 其中MH算法是常用的MCMC方法之一[12].
算法2\V(σ), 一個選定分布q, 最大樣本數(shù) M1和一個初始點σ0;
輸出: 后驗樣本集合{σm}M1m=1.
for m=1∶M1
提議一個候選點: σc~q(·σm-1);
計算接受概率: α(σc,σm-1)=min1,πV
(σc)q(σm-1σc)πV(σm-1)q
(σcσm-1);
生成一個樣本u~U(0,1);
if ult;α(σc,σm-1)
σm=σc;
else
σm=σm-1;
end if
end for
雖然MH采樣算法可以計算真解的后驗分布, 但無論是Tikhonov正則化方法, 還是MH采樣算法, 都需要求解大量的正問題(5). 雖然每個正問題求解只需幾秒鐘, 但計算大量正問
題, 會增加整個算法的計算時間. 本文用神經(jīng)網(wǎng)絡的方法構(gòu)造替代模型, 以減少整體算法的計算量.
2.3"神經(jīng)網(wǎng)絡替代模型
全連接神經(jīng)網(wǎng)絡(FNN)是一種特殊結(jié)構(gòu)的神經(jīng)網(wǎng)絡, 它由一系列全連接層組成, 對非線性函數(shù)可以有效地逼近. 全連接神經(jīng)網(wǎng)絡由(L+1)層網(wǎng)絡組成, 其中前一層的節(jié)點xl-1與
下一層的所有節(jié)點xl都有連接. 在每個節(jié)點上, 有兩個函數(shù)依次復合構(gòu)成. 第一個是一次函數(shù)
Ll(xl-1)=Wlxl-1+bl:nl-1→nl, 第二個是激活函數(shù)σl:nl→nl(例如Sigmoid,tanh,ReLU), 其中Wl表示權(quán)重, bl表示偏置. 圖1為FNN的示意圖. 為簡便,
令Θ={Wl,bl}Ll=1表示FNN的所有網(wǎng)絡參數(shù), 則可以定義一個由Θ參數(shù)化的映射FNN:m→n:
y=FNN(z;Θ)=LLσL-1LL-1…σ1L1(z).(9)
用式(9)求解正問題(5), 即z=σ, y=V, FNN≈T, m=d, n=N×(M+1). 設訓練集為{(σi,Vi)}n1i=1, 測試集為{(σi,Vi)}n2
i=1, n*=n1+n2表示樣本集的總數(shù). 對于監(jiān)督學習, 考慮平方損失函數(shù), 即
Loss(Θ)=12n1∑n1i=1‖Vi-FNN(σi;Θ)‖2.(10)
本文主要目標是通過適當?shù)膬?yōu)化算法, 找到最優(yōu)參數(shù)Θ*, 即Θ*=argminΘ Loss(Θ).(11)
解決問題(10)的常見優(yōu)化算法目前已有很多, 如梯度下降法、 牛頓法和隨機梯度下降法等. 對本文所考慮的問題, 通常n1很大, 導致計算梯度的代價很高. 為減少相應的計算量, 這里使用最小批量梯度下降方法. 該方法將訓練集劃分為小批次, 用于計算損失函數(shù)誤差和梯度, 并更新Θ, 綜合了隨機梯度下降法和梯度下降法的優(yōu)勢.
關(guān)于替代模型, 本文有如下收斂性結(jié)果. 假設:
(H1) 散度函數(shù)f滿足
f(t)≤C1log t,0lt;tlt;1,C2tlog t,t≥1,
其中C1和C2是依賴于f的兩個常數(shù).
引理2[11]"對εgt;0, 存在C=C(ε)gt;0, 使得對任意σ∈d, 有
‖T(σ)-FNN(σ;Θ)‖≤Cexp{ε‖σ‖2}ΦT(Θ),
其中ΦT(·)表示與T相關(guān)的任何函數(shù), 只需滿足lim#{Θ}→∞ ΦT(Θ)=0, #{Θ}表示神經(jīng)網(wǎng)絡參數(shù)的個數(shù).
設μ0為先驗測度, μ為離散算子T(·)對應的真實模型后驗測度, μF為神經(jīng)網(wǎng)絡FNN(·;Θ)替代的后驗測度, 下列定理給出了在f散度下兩者誤差的上界.
定理1"假設(H1)成立, 設測度μ,μFμ0, T(·)和FNN(·;Θ)都是
Lipschitz函數(shù). 則測度μ和μF的誤差在f散度下是收斂的, 即存在常數(shù)Cgt;0, 使得
Df(μ‖μF)∶=∫dfdμ1/dμ0dμ2/dμ0
dμ2dμ0dμ0(σ)≤CΦG(Θ).(12)
定理1的證明類似文獻[11], 首先利用正演模型和替代模型是Lipschitz函數(shù)的性質(zhì), 以及相應逼近的假設, 并基于Fernique定理, 得到歸一化常數(shù)的估計. 然后, 根據(jù)散度函數(shù)f具有分段上界的假設, 以及正演模型和替代模型的假設, 將f散度進行分段估計, 從而可得結(jié)論. 定理1表明, 當本文替代神經(jīng)網(wǎng)絡模型與正演模型足夠逼近時, 在f散度意義下相應的替代后驗測度收斂到真實模型的后驗測度, 故替代模型有一定的理論保障.
3"數(shù)值實例
下面進行數(shù)值實驗(d=1,3)驗證本文提出算法的有效性. 對于正問題(5), 選擇T=1, r=0.03, q=0.01, 并將轉(zhuǎn)換后的截斷長度設為L=1.6, 以確保截斷區(qū)域足夠大. 對于離散化問題, 時間和空間分割分別為M=150和N=100.用MH和MH+FNN算法針對式(6)進行采樣, 并設先驗分布π0(σ)~Nd(-1.8×1d,0.1×Id), 取密度q(·;σ)~Nd(σ,10-6×Id). 此外, 采樣噪聲振幅為δ=0.012. 在產(chǎn)生σ的若干樣本后, 可利用算法1對現(xiàn)性互補問題(4)進行求解, 得到正問題的若干組解(即期權(quán)離散矩陣V), 并作為監(jiān)督學習目標泛函中的真解. 然后, 以σ的若干樣本為輸入?yún)?shù), 正問題的若干組解為輸出, 經(jīng)過優(yōu)化訓練,得到神經(jīng)網(wǎng)絡FNN的參數(shù), 并用該網(wǎng)絡代替正問題的求解.
取神經(jīng)網(wǎng)絡結(jié)構(gòu)為{d,50,50,50,50,151×100}, d=1,3. 使用Adam優(yōu)化器對網(wǎng)絡進行訓練, 每次迭代使用128個樣本的小批量, 并進行5 000次循環(huán)的訓練, 而權(quán)重初始化通過Glorot均勻初始化.這里, 每次使用原始對偶活躍集方法, 對式(6)進行一次求解, 可得一個期權(quán)離散矩陣V, 即為一個樣本. 取總樣本數(shù)n*=20 000, 其中訓練集的樣本數(shù)n1=16 000, 測試集的樣本數(shù)n2=4 000. 所有的數(shù)值實驗都在一臺Intel Core i7的計算機工具包MATLAB R2020b上實現(xiàn).
例1"考慮如下一維的波動率:
σ(t)=σa(t),""a(t)=0.3et.(13)
表1列出了當σ=0.15,0.3,0.45時, Bayes方法(MH)和Bayes+替代模型(MH+FNN)的數(shù)值結(jié)果, 其中時間包括線下時間和線上時間. 這里, FNN的實驗輸出的是期權(quán)價
格V, 但MH和MH+FNN實驗輸出的是波動率σ. 由表1可見, 兩種方法的精度近似, 但后者花費的時間只有前者的1/3左右, 故神經(jīng)網(wǎng)絡作為替代模型, 加速效果顯著.
表2列出了MH+FNN算法固定每層神經(jīng)元為50時, 不同隱藏層的數(shù)值結(jié)果. 表3列出了MH+FNN算法固定隱藏層為3時, 不同神經(jīng)元個數(shù)的數(shù)值結(jié)果. 表2和表3的
數(shù)值結(jié)果表明, 隨著隱藏層和神經(jīng)元個數(shù)的增加, 精度會有些提高, 但會增加一定的時間. 對于增加樣本, 有類似的結(jié)果, 這里不再列出.
例2"考慮含3個參數(shù)的波動率, 即σ(t)=∑3i=1σiai(t),其中a1(t)=1, a2(t)=t, a3(t)=3-t, 波動率的真解為σ1=0.1, σ2=0.15, σ3=0.2. 表4列出了(σ1,σ2,σ3)=(0.1,0.15,0.2)時, Bayes方法(MH)和Bayes+替代模型(MH+FNN)的數(shù)值結(jié)果. 與一維的結(jié)果類似, MH+FNN的計算結(jié)果和MH類似, 但時間上要快2/3左右.
圖2給出了MH和MH+FNN算法采樣的相關(guān)圖. 由圖2可見: (A),(C),(F)是自相關(guān)圖, 其兩條線基本重合, 說明自相關(guān)度非常高; (B),(D),(E)是協(xié)相關(guān)圖, 其各自的兩個環(huán)
線也基本重合, 說明協(xié)相關(guān)度很高. 實驗結(jié)果表明, 兩種采樣效果一致.
參考文獻
[1]"ACHDOU Y, INDRAGOBY G, PIRONNEAU O.
Volatility Calibration with American Options [J]. Methods Appl Anal, 2004, 11(4): 533-556.
[2]"BURKOVSKA O, GLAU K, MAHLSTEDT M, et al. Complexity Reduction for Calibration
to American Options [J]. J Comput Finance, 2019, 23(1): 25-60.
[3]"KUTNER G W. Determining the Implied Volatility for A
merican Options Using the QAM [J]. "Financ Rev, 1998, 33(1): 119-130.
[4]"HAN J Q, JENTZEN A, E W N. Solving High-Dimensional Partial Differential Equations Using Deep L
earning [J]. Proc Natl Acad Sci USA, 2018, 115(34): 8505-8510.
[5]"LIU S Q, LEITAO , BOROVYKH A, et al. On a Calibration Neural Network to Extract Implied Information fr
om American Options [J]. Appl Math Finance, 2020, 28(5): 449-475.
[6]"BLACK F, SCHOLES M. The Pricing of Options and Corporate Liabilities [J]. J Polit Econ, 1973, 81(3): 637-654.
[7]"GAO Y, SONG H M, WANG X S, et al. Primal-Dual Activ
e Set Method for Pricing American Better-of Option on Two Assets [J]. Commun Nonlinear Sci Numer Simulat, 2020, 80: 104976-1-104976-15.
[8]"SONG H M, WANG X S, ZHANG K, et al. Primal-Dual Active Set Method for American Lookback Put Option Pricing [J]. "East Asian J Appl Math, 2017, 7(3): 603-614.
[9]"BOUCHOUEV I, ISAKOV V. The Inverse Problem of Option Pricing [J]. Inverse Probl, 1997, 13(5): L11-L17.
[10]"LEHTINEN M S, PIVRINTA L, SOMERSALO E. Linear Inverse Problems for Generalised Random Variables [J].
Inverse Probl, 1989, 5(4): 599-612.
[11]"GAO Y, LIU H Y, WANG X C, et al. A Bayesian Scheme
for Reconstructing Obstacles in Acoustic Waveguides [J]. J Sci Comput, 2023, 97(3): 53-1-53-45.
[12]"HASTINGS W K. Monte Carlo Sampling Methods Using Markov Chains and Their Applications [J]. Biometrika, 1970, 57(1): 97-109.
(責任編輯: 李"琦)