期权定价的有限差分方法.docx
《期权定价的有限差分方法.docx》由会员分享，可在线阅读，更多相关《期权定价的有限差分方法.docx（22页珍藏版）》请在冰豆网上搜索。
期权定价的有限差分方法
第九章期权定价的有限差分方法
在本章中，我们将给出几个简单的例子来说明基于偏微分方程（PDE框架的期权定价方法。
具体的方法的是利用第五章中讲述的有限差分方法来解决Black-scholes偏微分方程。
在9.1节中，我们会回顾衍生品定价的数值解法以及指出如何利用适当的边界条件来模拟一个特定的期权。
在9.2节中我们将会应用简单的显式（差分）方法来求解一个简单的欧式期权。
正如你已熟知的那样，这种方法只能解出一些可以从金融角度来解释的不稳定的数值解。
在9.3节中我们
将可以看到使用完全的隐式方法可以解决这种不稳定问题。
在9.4节
中我们将介绍Crank-Nicolson方法在障碍期权定价中的应用，它可以看做是一种显式与完全隐式方法的混合。
最后，在9.5节中，我们会看到迭代松弛方法可以用于解决使用全隐式方法来解决美式期权定价时由于存在提前执行的可能性而导致的自由边界问题。
9.1使用有限差分法解BS方程
在2.6.2节中，我们给出了一个标的资产在时间t的价格为S（t）的期权，该期权的价格是一个函数f（S,t）,且f（S,t）满足偏微分方程客…翼+慧”（9.1）
通过不同的边界条件可以让这个方程刻画不同的期权的特征。
在
某些地方可能因为假设的改变或者对路径依赖的改变而导致方程式的具体形式改变，但是此处仅仅作为一个起点，帮助读者了解如何应用基于有限差分方法来解决期权定价的问题。
正如我们在第五章中遇到的情况那样，要用有限差分方法来解偏微分方程，在此处我们必须建立资产价格和时间的离散网格。
设T是
期权的到期日，而Sma娓一个足够大的资产价格，在我们所考虑的时间范围内，S⑴的数值不能超过Smax设定SmaxM因为偏微分方程的区域关于资产价格是无边界的。
但是为了达到计算的目的，必须要求它是有界的。
Sma而当于+8。
网格通过点（S,t）取得，其中（S,t）满足
S=0,6S,25S,
t=0,比，2dt2,,NSt=T。
本章中使用网格符号为二/四品）阴，我们回顾一下（9.1）方程式的几种不同解法：
向前差分
dS53
向后差分
dfj
dS
52/=（九厂九_』J一力
爵~"8S8S
"+1厂2力十九.
as?
至于究竟采用哪种方法进行离散化，我们将在后面的实际操作过程中对显式和隐式的方法作出详细的阐述说明
另一个值得我们注意的问题是如何设置边界条件。
对于执行价格
为K的欧式看涨期权，其终值条件为
/GT）=max{S-KQVS
对于执行价格为K的欧式看跌期权，其终值条件为
/⑸T）=max{K—S,0}VS
但是当我们涉及到资产价格的边界问题时，这个问题就变得复杂了。
因为在数学上要求必须在一个有界区域里来解这个方程，而资产价格的区域是无界的。
我们可以通过一些例子来说明这个问题。
例9.1我们首先考虑一支简单的欧式看跌期权。
当资产的价格S（t）非常大时，看跌期权几乎是一文不值，因为我们几乎可以肯定的是其执行价低于现价：
〃4皿/）=。
Smax的值必须足够大才能使得假设的边界条件是合理的。
当资产的价格S⑴=0时，鉴于资产的价格符合几何布朗运动的动力学模型，我们可以认为资产的价格会保持为0;因此到期收益为K,在时间t贴现值为
其边界条件可以写为:
£N=T3&0],i=0,LM
九炉°，，=。
工…，N
例9.2我们可以从9.1推导出处理普通欧式看涨期权的方法。
不妨设资产的价格S（t）=0,那么在任何时刻t,我们都可以知道该期权价值为0:
而对于一个足够大的资产价格S（t）,我们可以肯定的是该期权一定为价内期权，即具有大于零的价值的期权，那么期满行权我们可以得到收益为S（T）_K。
将终值K贴现到t时刻，并且考虑标的资产的价格是一个函数S⑴，那么我们根据无套利定价原理可以得到其边界条件如下：
/（Smax二Sjnax-底卬-。
具体条件可以写为：
工小—…郎-&0Li=QL……M
"尸题-鼓tj=OJ.，…A'
对于资产价格S较大时，需要一个替代的边界条件，但这时候必须要求期权的△为1；在这种情况下，我们所提出的边界条件一般针对未知函数的倒数，而非函数本身。
这便是数学物理方程中所谓的Neumann边界条件。
我们不会采取这种做法，因为这将使得数值解法复杂化了。
当我们面对障碍式期权时，事情可能就变得简单了。
对敲出期权来说，例如一个向下触及失效看跌期权，在标的资产价格触碰到障碍时期权作废，价值为0。
向上触及失效看涨期权的情况也是类似的，它们的优点都在于我们计算必须涉及到的边界即域是自然有界的（显
然存在的）。
美式期权的定价就显得更加复杂难以处理了，因为这涉及到了提前执行边界的问题。
我们必须考虑什么水平的资产价格和在什么时间是期权的最佳实施点。
因此，在我们解决问题的过程中必然产生一个自由边界。
这就需要为不同的奇异期权设定不同的边界条件。
如何计算出正确的边界条件以及使用数值方法逼近正确的结果，这依赖于对具体期权的选择。
利用显式法求解欧式期权
由于初次尝试解方程（9.1）,我们不妨考虑一个简单的欧式看跌期权。
我们分别用中心差分对S进行求导，用向后差分对时间进行求导来逼近求解。
这并不是唯一的方法，但是无论采取什么样的方法都必须保证边界条件是符合特定期权模型的。
例9.1中的边界条件差分后的结果可以用以下的方程式来表示：
（9.2）
需要特别指出说明的是，由于存在终点条件的假设，因此方程必须适当地用回推法求解。
不妨设（9.2）式中有j=N,由于存在给定的终点条件，我们有一个未知量fi，N」，可以认为它是关于三个已知量的函数。
对其他的每个时间层（i）作同样的考虑。
整理方程，我们可以得到：
%~广4：
力tj+号叫+啦&j
（9.3）
/二N—LN—2rLeL
此时有
b；=1-况
functionprice=EuPutExpl（S0,K,r,T,sigma,Smax,dS,dt）
犍立网格，并且在必要时调整增量
M=round（Smax/dS）;
dS=Smax/M;
N=round（T/dt）;
dt=T/N;
matval=zeros（M+1,N+1）;
vetS=linspace（0,Smax,M+1）;
veti=0:
M;
vetj=0:
N;
犍立边界条件
matval（:
N+1）=max（K-vetS,0）;
matval（1,:
）=K*exp（-r*dt*（N-vetj））;
matval（M+1,:
）=0;
犍立三对角矩阵
a=0.5*dt*（sigmaA2*veti-r）.*veti;
b=1-dt*（sigmaA2*veti.A2+r）;
c=0.5*dt*（sigmaA2*veti+r）.*veti;
用解方程
forj=N:
-1:
1
fori=2:
M
matval（i,j）=a（i）*matval（i-1,j+1）+b（i）*matval（i,j+1）+…
c（i）*matval（i+1,j+1）;
end
end
领回价格，有可能在网格外线性插值生成price=interp1（vetS,matval（:
1）,S0）;
图9.1简易欧式期权显式差分方法定价的MATLAB弋码
定价不是直接利用MATLA来完成白图9.1中是定价所使用的MATLA的程序代码，另外计算还需要给出Smax的值以及上文提到的
两个离散步骤。
唯一一点需要注意的是，在数学符号中哥指数可以从0开始，但是在MATLA的矩阵计算中，哥指数是从1开始的。
另外,如果资产的初始价格不在网格上，那么我们只能够寻找两个相邻点进行插值。
一般我们采取比较粗略的线性插值。
但是如果有需要计算分
析期权价格的敏感性（在实际案例中经常需要这样做）则需要进行更加复杂的插值计算。
>>[c,p]=blsprice（50,50,0.1,5/12,0.4）;
>>P
P=4.076
>>EuPutExpl（50,50,0.1,5/12,0.4,100,2,5/1200）ans=4.0669
>>[c,p]=blsprice（50,50,0.1,5/12,0.3）;
>>P
P=2.8446
>>EuPutExpl（50,50,0.1,5/12,0.3,100,2,5/1200）
ans=2.8288
我们可以看到，数值方法能够提供比较精确的结果，为了提高精
度，我们也可以采用更加精细的网格进行计算。
>>EuPutExpl（50,50,0.1,5/12,0.3,100,1.5,5/1200）
ans=2.8597
>>EuPutExpl1（50,50,0.1,5/12,0.3,100,1,5/1200）
ans=-2.8271e+022
在此处我们看到在第五章中分析过的不稳定数值解的例子。
有一种解决方案就是使用隐式方法进行求解。
另外一种方法则是进行稳定性分析，然后在从离散化过程中导出边界。
此处我们不会再使用第二
种方案，这与第五章中对热传导方程进行的简易变换处理方法十分相似。
接着，本章在下一节中将从金融的角度来阐释不稳定性，这暗示我们需要通过改变变量来改变方程式。
9.2显式法不稳定性的金融学解释
在显式法中，我们通过f（S+3S,t+5t）,f（S,t+『）和f（S—gt+Et）三
个值得到f（S,t）。
这看起来与我们在7.4这部分（见图9.2（a））介绍过的三叉树方法有点类似。
这样方便我们从显式法推导出新的方
法。
如下【参考文献1中第18章】，假设在点（i,j）和点（i,j+1）对S求一阶导数和二阶导数是相等的：
4ij+i熬一放
♦2/—+1J+1-2工j+i1Z-L#1
另外一种方法可以取得同样的效果，用品」代替9.2式右侧的fi,jO这种方法引入的误差是有界的，并且这种方法的网格细化至趋于零。
改变过后的替换后的差分方程变成了：
/J+]一—十】J+1———LJ+1
-a2SS
+1/+LJ+12片J+[+/—LJ+1
2市
也可写成（当i=1,2,…,M-1且j=1,2,…,N-1时）：
%二用力-1j+i+，工J+1十弓加j+i
同时有：
&=i।£（_1方拉+2M=1翼
:
227T+rSid
h.=1（l—仃2超阳=1r九
用函*1+rdr0
c一L_（1fiJt+g仃』2§讣二1笈工\+rSi221+rSlu
这个方法仍然是显式的而且仍然受到不稳定性的影响。
系数哥，
白以及己本身也对此作出了很好的解释。
回顾之前我们使用过的二叉树和三叉树方法，我们是在风险中性的假设下将后继节点对期权的预期价值进行贴现得到期权的价格。
事实上，上述系数，包括1/（1+咖）可以看做是长度为黄的时间间隔的贴现系数。
此外我们还有：
71'\Jln\=1
这也说明了上述系数等于概率乘以一个折现系数。
他们是否是风险中性的概率呢？
我们应当首先测试资产价格在时间间隔6t内增长的值：
A]=—SSttd4-十=ri6S3t=rSSt
这正是我们所期望的风险中性世界，对于其方差的增加，有：
此外，当&极小时有：
=即/i2]一E*[A]=这也是服从风险中性假设的几何布朗运动的。
因此我们可以认为
排除一些小问题后，显式法是可以和三叉树方法取得同样的效果的。
即“概率"%和队可能是负数。
细心的读者可能会发现一个循环，在第五章中我们学习了一个线性组合的系数的稳定性条件。
在热传导
方程中，我们必须确定这个组合是一个凸组合，即系数是正的并且相加等于1。
就像一个概率函数。
有一种解决这个问题的方法，就是【参考文献U中通过改变变
量来实现的。
通过Z=lnS来代换重写B-S方程，我们得到了实现稳定性的简单条件。
然而，变量的改变对某些奇异期权来说却也未必是好事。
在下一节中我们将使用隐式法来避免数值不稳定的问题。
9.3使用全隐式法对欧式期权定价
为了克服显式法带来的不稳定性问题，我们采用全隐式法。
首对
时间t采用向前差分求其偏导数。
我们得到一个网格方程：
7+1F+对5sL+L■/
3tJ26S
25s2J7
也可写成（当i=1,2,…,M-1且j=1,2,…,N-1时）：
此时有:
M1+1J7KJ-1
工22
b=1+c.=—JriSt-1-n-^i^St
上22
此时我们有需要使用3个未知量去计算出一个未知量（见图
9.2b）。
首先要注意的是，在每一个时间层，我们有M-1个方程式和
N-1个未知量。
边界条件决定两个未知量，终点条件也给出了最后时间层的值。
正如显式法的案例那样，我们必须用回推法，解决当j=N-1,…/情况下的所有线性方程组。
在时间层j,得到的方程组如
0
值得我们注意的是，上述的矩阵是个三对角矩阵，并且相对于时间层i来说它是一个不变的常量。
因此我们可以使用LU因式分解来
解决问题1。
所有工作的MATLA歌码在图9.3中给出
functionprice=EuPutImpl（S0,K,r,T,sigma,Smax,dS,dt）
犍立网格并在必要时调整增量
M=round（Smax/dS）;
dS=Smax/M;
N=round（T/dt）;
dt=T/N;
matval=zeros（M+1,N+1）;
vetS=linspace（0,Smax,M+1）';
veti=0:
M;
vetj=0:
N;
%建立边界条件
matval（:
N+1）=max（K-vetS,0）;
matval（1,:
）=K*exp（-r*dt*（N-vetj））;
matval（M+1,:
）=0;
%建立三对角矩阵系数
a=0.5*（r*dt*veti-sigmaA2*dt*（veti.A2））;
b=1+sigmaA2*dt*（veti.A2）+r*dt;
c=-0.5*（r*dt*veti+sigmaA2*dt*（veti.A2））;
coeff=diag（a（3:
M）,-1）+diag（b（2:
M））+diag（c（2:
M-1）,1）;
[L,U]=lu（coeff）;
%求解线性方程组
aux=zeros（M-1,1）;aux
（1）=-a
（2）*matval（1,j）;%othertermfromBCiszeromatval（2:
M,j）=U\（L\（matval（2:
M,j+1）+aux））;
end
%返回价格，如果初始价格在网格外则可能进行性线性插值
price=interp1（vetS,matval（:
1）,S0）;
图9.3全隐式方法求解简易欧式期权的MATLAB
>>[c,p]=blsprice（50,50,0.1,5/12,0.4）;
>>P
P=4.0760
>>EuPutlmpl（50,50,0.1,5/12,0.4,100,0.5,5/2400）
ans=4.0718
这个计算结果是非常精确的，这是因为使用了排除不稳定因素存在的更加精细的网格的结果。
另外一种方法就是通过Crank-Nicolson
方法来提高精度，这将在后面的章节中做进一步的介绍。
9.4使用Crank-Nicolson方法对障碍式期权定价
Crank-Nicolson方法已经在5.3.3节中提到过，它通过结合使用显式法和隐式法来提高数值解的精确度。
在B-S模型的求解中应用此方法可以得到：
%—Aj-i.6sf二i-八一i上6s（~/i-1j\
一St-F\US）~〈26S）
4I（㈣2J
♦汽3S）2
rv
=①£，，-1十
上述方程可以改写为:
+（1-A）Aj-l_++（1+
（9.5）
此时有:
-ri）
4（村+『）^（cr2i2+ri）,4
我们首先考虑在2.7,1中介绍过的向下触及失效期权，假设该障碍期权是处于连续监测下的。
在这种情况下我们需要考虑&Smax边界条件是：
/（5maxt0=0,/（Sfc,£）—0.
综合考虑上述边界条件，我们可以将方程（9.5）改写为如下矩阵形式:
（9.6）
=M泻,
同时，有:
一
1一色
—g
i/ij?
/sjj---j—i,j]■
相应的MATLA献码在图9.4中列出
这个结果可以同2.7.1中通过定价分析公式得到的结果进行比
>>DOPut（50,50,0.1,5/12,0.4,40）
ans=
0.5424
>>DOPutCK（50,50,0.1,5/12,0.4,40,100,0.5,1/1200）
ans=
0.5414
障碍期权还有很多种，可以在【参考文献9］中找到更多的关于
障碍式期权定价的偏微分方程应用。
functionprice=DOPutCK（S0,K,r,T,sigma,Sb,Smax,dS,dt）
%建立网格并且在必要的时候调整增量
M=round（（Smax-Sb）/dS）;
dS=（Smax-Sb）/M;
N=round（T/dt）;
dt=T/N;
matval=zeros（M+1,N+1）;
vetS=linspace（Sb,Smax,M+1）';
veti=vetS/dS;
vetj=0:
N;
%建立边界条件
matval（:
N+1）=max（K-vetS,0）;
matval（1,:
）=0;
matval（M+1,:
）=0;
%建立系数矩阵
alpha=0.25*dt*（sigmaA2*（veti.A2）-r*veti）;
beta=-dt*0.5*（sigmaA2*（veti.A2）+r）;
gamma=0.25*dt*（sigmaA2*（veti.A2）+r*veti）;
M1=-diag（alpha（3:
M）,-1）+diag（1-beta（2:
M））-diag（gamma（2:
M-1）,1）;
[L,U]=lu（M1）;
M2=diag（alpha（3:
M）,-1）+diag（1+beta（2:
M））+diag（gamma（2:
M-1）,1）;
%求解线性方程组
forj=N:
-1:
1
matval（2:
M,j）=U\（L\（M2*matval（2:
M,j+1）））;
end
%返回价格，这个结果可能因为初始资产价格在网格外而有线性插值产生。
price=interp1（vetS,matval（:
1）,S0）;
图9.4Crank-Nicolson方法对向下触及失效期权定价的MATLAB码。
9.5美式期权定价
有限差分方法在对简单的欧式期权进行定价时的确是十分方便
的，但是这在美式期权定价中却显得并不是很实用。
我们将这个方法应用于美式期权定价的时候，精确的公式就不大适用了。
这是因为美式期权定价时必须考虑存在提前行权的可能性。
根据无套利原理，在
（S,t）空间中任一点的期权的价格不能低于其内在价值（即立刻行权的所得收益）。
对于一个一般的美式期权，这意味着：
从一个实际操作来看，在显式法中要充分考虑这种情况是很困难的。
但是我们可以做一些修改以便能够简单的应用9.2节中提到的方
法。
在计算得到fi，j后，我们可以检测提前行权的概率。
并且可以像我们处理二叉树模型那样设：
fij=max伉-K-iSS],由于存在着不稳定性问题，我们可能更加倾向于使用隐式法。
但
是这样处理后又导致了另一个复杂的问题，即根据上述关系要求,
必须是已知的，而在隐式法中却并不能已知均要解决这个问题，我
们可以采取一个迭代的方法来求解这个线性方程组，而不是基于LU
分解直接求解。
在3.2.5节中我们介绍了Gauss-Seide松弛迭代算法。
为了方便我们先回顾一下这个方法。
给出一个线性方程组，如下：
Ax=b>
我们将使用下面的迭代方法，初始点为x（0）：
（i-lN\
瓦一^口同川二工国同>i=l，…，N,其中k是迭代计数器，而⑴是松弛变量。
该迭代将会一直进行到
满足收敛条件，如：
|限佳+1）一⑹“r
其中£是容限参数。
现在我们将尝试使用Crank-Nicolson方法对美式期权进行定价。
我们要解决的是和（9.6）式大致相同的一个线性方程组，但是这里
的边界条件有一点不同，价值为0的期权不用再考虑是否有障碍。
我
们要用回推法求解的线性方程组如下：
M107=门］
其中右边有：
/oj-1+/oj
0
*=MaG+ai.
9
0
式子中额外增加的一项是为了将看跌期权独有的边界条件带入
模型方程。
在使用松弛方法的时候应当考虑矩阵M1的三对角性质,
并且这能够适应美式期权提前执行的特征。
设gi为内在价值，i=1,2,…，M-1,同时S=i6S。
对于每个时间层j,我们进行如下迭代:
，严=口皿｛班,居
+E?
卜-a-①）君+"到｝
段*"=max虑
+卜空+山，产）-（1-伪）唐、为石？
］｝
rAf-J+~~Alf-1）皮Li4
当计算从一个时间层进行到下一个时间层的时候，将上一个时间层的计算结果作为下一个时间层迭代计算的初始向量可能是更加合理的。
完成上述步骤需要的代码在图9.5中将列出。
使用该代码有点麻烦，原因在于MATLA的向量募运算是从1开始的，而现实却要求从0开始。
这样我们便不能建立一个包含所有fi，j的向量矩阵，同时系数矩阵M1也无法储存。
上述迭代最好是直接通过口、P和¥计算得
该代码的实际效果可能要比MATLA金融工具箱中附带的使用二叉树方法对美式期权定价的binprice函数（见7.1节）要好。
>>tic,[pr,opt]=
binprice（50,50,0.1,5/12,1/1200,0.4,0）;,toc
Elapsedtimeis0.408484seconds.
>>opt（1,l）
ans=
4.2830
>>
tic,AmPutCK（50,50,0.1,5/12,0.4,100,1,1/600,1.5,0.001）,toc
ans=
4.2815
Elapsedtimeis0.031174seconds.
>>
tic,AmPutCK（50,50,0.1,5/12,0.4,100,1,1/600,1.8,0.001）,toc
ans=
4.2794
Elapsedtimeis0.061365seconds.
>>
tic,AmPutCK（50,50,0.1,5/12,0.4,100,1,1/600,1.2,0.001）,toc
ans=
4.2800
Elapsedtimeis0.023053seconds.
>>
tic,AmPutCK（50,50,0.1,5/12,0.4,100,1,1/1200,1.2,0.001）,toc
ans=
4.2828
Elapsedtimeis0.036693seconds.
>>
tic,AmPutCK（50,50,0.1,5/12,0.4,100,1,1/100,1.2,0.001）,tocans=
4.2778
Elapsedtimeis0.009989seconds.
functionprice
AmPutCK（S0,K,r,T,sigma,Smax,dS,dt,omega,tol）
M=round（Smax/dS）;dS=Smax/M;%建立网格
N=round（T/dt）;dt=T/N;
oldval=zeros（M-1,1）;%Gauss-Seidel更新向量
newval=zeros（M-1,1）;
vetS=linspace（0,Smax,M+1）';
veti=0:
M;vetj=0:
N;
%建立边界条件
payoff=max（K-vetS（2:
M）,0）;
pastval=payoff;%valuesforthelastlayerboundval=K*exp（-r*dt*（N-vetj））;%边界值%建立系数矩阵和式子右边矩阵
升级会员 