1973年Black和Scholes推导出著名的Black-Scholes期权定价公式^[1],是金融数学领域的一个里程碑.然而, Black-Scholes模型的一系列假设(如波动率为常数,价格过程服从对数正态分布、市场是完备的)使得期权价格与市场观察价格存在较大差异.尤其1987年金融危机的爆发, Black-Scholes模型中期权价格与市场报价的偏差更加显著.从那时起,人们不断的探寻能够精确复制股票价格动态的定价模型.指数Lévy模型通常被当作Black-Scholes模型的替代模型,因为Lévy过程允许股票价格跳跃,能够有效的描述股票价格具有"尖峰厚尾"及对数收益非对称等特点.

跳-扩散模型^[2]是指数Lévy模型中一类非常重要的模型,该模型最早由Merton提出.跳-扩散模型在Black-Scholes模型基础上增加非连续的跳跃,连续部分表示的是每天价格过程的演变,跳跃代表稀少事件.与Black-Scholes模型相比,跳-扩散模型更具有灵活性^[3],从而能更好的模拟金融市场数据.然而,跳-扩散模型下市场是不完备的,此时存在许多等价鞅测度,从而对应着有许多无套利价格.为了得到一致的无套利价格,需要在所有等价鞅测度中选择一个最优(相对于某些标准)的鞅测度^[4],均值修正法是获取等价鞅测度的一种常用的方法^[5].另外,指数Lévy模型下封闭形式的概率密度函数往往是未知的,但其特征函数通常是可以得到的,因此利用特征函数对期权进行定价是非常必要的. Carr和Madan提出基于特征函数计算期权价格的快速傅里叶变换算法(FFT)^[6],在此基础上, Fang和Oosterlee提出了利用傅里叶余弦级数计算欧式期权价格的COS方法^[7],该方法利用了特征函数与密度函数傅里叶余弦展开系数之间的特殊关系,其收敛速度和计算难度都优于FFT方法.

在金融市场上,无论是对标准期权还是对奇异期权进行定价,都需要对模型中的未知参数进行估计,使得理论值与市场观测值尽可能的相匹配,这就是模型的参数校准问题.参数校准问题是金融中的反问题,往往是不适定的:解不稳定、不唯一或者根本不存在.此时,需要增加罚函数对校准问题进行正则化,从而得到校准问题唯一解. Cont和Tankov提出了相对熵正则化方法,并对指数Lévy模型中的部分模型进行校准^[8-10].

本文对跳-扩散模型下期权定价及参数校准问题进行研究,推导出了均值修正等价鞅测度下跳-扩散模型的风险中性特征函数,利用COS方法对跳-扩散模型进行定价,并证明了算法的收敛性.然后,利用相对熵正则化方法对跳-扩散模型参数进行校准,数值实验表明该方法是准确、可靠的.最后利用S&P500期权市场数据对Merton跳-扩散模型进行参数校准,得到了较好的校准结果,利用校准参数重构的波动率曲面与市场实际波动率曲面相吻合.

假设市场上有两种资产,一种为无风险资产(如债券),其价格过程为$B_t={\rm e}^{rt}$;另一种资产为风险资产(如股票),其价格过程$S_t$为定义在概率空间$(\Omega, {\cal F}, {\Bbb P})$上的指数Lévy过程

(2.1)$\begin{equation} S_{t}=S_{0}\exp(X_{t}), \, 0\leq t\leq T, \end{equation}$

其中$S_0>0$是常数, $X_t$是一维Lévy过程,其特征三元组为$(\sigma^{2}, \nu, \gamma)$,其中$\sigma$扩散系数, $\nu$是$X_t$的Lévy测度, $\gamma$是漂移项.由Lévy-Khintchine公式可知, $X_{t}$的特征函数具有如下形式

$\Phi(u)={\bf E}^P[\exp({\rm i}u X_{t})]=\exp(t\varphi(u)), $

其中${\bf E}^P[\cdot]$为概率测度${\mathbb{P}}$下的数学期望算子, $\varphi(u)$为特征指数,且

$\varphi(u)={\rm i}\gamma u-\frac{1}{2}\sigma^{2}u^{2}+\int_{-\infty}^{+\infty}({\rm e}^{{\rm i}u x}-1-{\rm i}u xI_{|x|\leq1})\nu({\rm d}x).$

如果关于$X_t$的Laplace变换存在,有${\bf E}^P[\exp(u X_t)]=\exp{(t\phi(u))}$,其中

(2.2)$\begin{equation} \phi(u)=u\gamma+\frac{1}{2}\sigma^2u^2+\int_{-\infty}^{+\infty}({\rm e}^{u x}-1-u xI_{|x|\leq1})\nu({\rm d}x). \end{equation} $

由Lévy-Itô分解公式, $X_t$可表示为

(2.3)$ \begin{equation} X_t=\gamma t+\sigma W_t+\int_0^t\int_{|x|>1}xN({\rm d}s, {\rm d}x)+\int_0^t\int_{|x|\leq1}x[N({\rm d}s, {\rm d}x)-{\rm d}s\nu({\rm d}x)], \end{equation}$

其中$W_t$为一维布朗运动, $N({\rm d}s, {\rm d}x)$为泊松跳跃的随机测度. (2.3)式中第一项为漂移项,第二项为布朗运动,第三项为复合泊松过程,第四项为纯跳跃鞅.

Black-Scholes模型的驱动过程可由$X_{t}=\gamma t+\sigma W_{t}$给出,其特征函数为

(2.4)$\begin{equation}\Phi(u)={\bf E}^P[{\rm e}^{{\rm i}uX_{t}}]=\exp\bigg\{t({\rm i}\gamma u-\frac{1}{2}\sigma^{2}u^{2})\bigg\}.\end{equation}$

基于Lévy过程的跳-扩散模型是同时包含布朗运动和复合泊松过程的随机过程,可由如下方程给出

(2.5)$\begin{equation} X_{t}=\gamma t+\sigma W_{t}+\sum\limits_{i=1}^{N_{t}}Y_{i}, \end{equation} $

其中$N=(N_{t})_{t\geq0}$是参数为$\lambda$的泊松过程, $Y_i$是跳跃幅度.如果$Y_{i}\sim N(\mu, \delta^{2})$,即为著名的Merton跳-扩散模型,其特征函数为

(2.6)$ \begin{equation} \Phi(u)={\bf E}^P[{\rm e}^{{\rm i}uX_{t}}]=\exp\bigg\{t({\rm i}u\gamma-\frac{1}{2}\sigma^{2}u^{2}+\lambda({\rm e}^{{\rm i}\mu u-\frac{1}{2}\delta^{2}u^{2}}-1))\bigg\}. \end{equation} $

对于指数Lévy模型,均值修正法是获取等价鞅测度的常用方法,其基本思想是通过修正Lévy过程的均值,获取某个与市场概率测度${\mathbb{P}}$等价的概率测度${\mathbb{Q}}$,使得股票的贴现价格过程$\{{\rm e}^{-rt}S_t\}$为${\mathbb{Q}}$鞅,相应的概率测度${\mathbb{Q}}$称为均值修正鞅测度.

定理2.1  假设$X_t$为定义在概率空间$(\Omega, {\cal F}, {\mathbb{P}})$上Lévy过程,其特征三元组为$( \sigma^2, \nu, \gamma)$.对$m\in{\Bbb R}$,记$Y_t=X_t-mt$,令${\mathbb{Q}}$为$Y_t$的有限维分布族在可测空间$(\Omega, {\cal F})$生成的概率测度.那么,当$m=r-\phi(1)$时,股票的折现价格$\{{\rm e}^{-rt}S_t\}$为${\mathbb{Q}}$鞅.

证  由于这种变换只改变了$X_t$特征三元组的漂移项,扩散项$\sigma$和Lévy测度$\nu$没变,因此$X_t$在概率测度${\mathbb{Q}}$下是特征三元组为$(\sigma^2, \nu, \gamma+m)$的Lévy过程.要想使$\{{\rm e}^{-rt}S_t\}$为${\mathbb{Q}}$鞅,则对于$0\leq \tau\leq t\leq T$,应有${\rm e}^{-r\tau}S_{\tau}={\bf E}^Q[{\rm e}^{-rt}S_t|{\cal F}_{\tau}]$,其中${\bf E}^Q[\cdot]$为概率测度${\mathbb{Q}}$下的期望算子.事实上,

(2.7)$ \begin{eqnarray}{\bf E}^Q[{\rm e}^{-rt}S_t|{\cal F}_{\tau}]&=&{\rm e}^{-rt}S_{\tau}{\bf E}^Q[{\rm e}^{X_t-X_{\tau}}]\\&=&{\rm e}^{-rt}S_{\tau}{\rm e}^{(t-\tau)[m+\phi(1)]}\\&=&{\rm e}^{-r\tau}S_{\tau}{\rm e}^{(t-\tau)[m+\phi(1)-r]} \end{eqnarray} $

有唯一解$m=r-\phi(1)$,此时$\{{\rm e}^{-rt}S_t\}$为${\mathbb{Q}}$鞅.

引理2.1^[12]  均值修正鞅测度${\mathbb{Q}}$存在的条件为$\sigma>0$,且$\int_{|x|>1}|{\rm e}^x-1|\nu({\rm d}x)<+\infty$.

定理2.2  设均值修正鞅测度${\mathbb{Q}}$存在, $X_t$在${\mathbb{P}}$和${\mathbb{Q}}$的特征三元组分别为$( \sigma^2, \nu, \gamma)$和$({\sigma^Q}^2, \nu^Q, \gamma^Q)$,则${\mathbb{Q}}$是${\mathbb{P}}$等价鞅测度的充要条件是存在$\beta\in{\bf R}$,使得

(2.8) $ \begin{equation} \gamma^Q=\gamma+\beta\sigma^2. \end{equation} $

如果${\mathbb{Q}}\sim{\mathbb{P}}$,则${\mathbb{Q}}$关于${\mathbb{P}}$导数满足

(2.9)$ \begin{equation} \frac{{\rm d}{\mathbb{Q}}}{{\rm d}{\mathbb{P}}}\bigg|{\cal F}_t=\frac{\exp(\beta W_t)}{{\bf E}^P[\exp(\beta W_t)]}. \end{equation}$

证  在均值修正鞅测度${\mathbb{Q}}$下, $\gamma^Q=\gamma+m$, ${\sigma^Q}^2=\sigma^2$, $\nu^Q=\nu$.由文献[5]中引理2.1可知(2.8)和(2.9)式成立,且参数$\beta$满足方程

$\beta\sigma^2=\gamma^Q-\gamma= m=r-\phi(1).$

由$\sigma>0$,则$\beta$有唯一解$\beta=\frac{r-\phi(1)}{\sigma^2}$.

由引理2.1和定理2.2可知均值修正鞅测度${\mathbb{Q}}$是${\mathbb{P}}$的等价鞅测度,下面给出Black-Scholes模型及Merton跳-扩散模型在测度${\mathbb{Q}}$下的风险中性特征函数.

命题2.1   Black-Scholes模型中, $\phi(1)=\gamma+\frac{1}{2}\sigma^2$,风险中性特征函数为

(2.10)$\begin{equation}\Phi(u)={\bf E}^Q[{\rm e}^{{\rm i}uX_{t}}]=\exp\left\{t\left[{\rm i}(r-\frac{1}{2}\sigma^2) u-\frac{1}{2}\sigma^{2}u^{2}\right]\right\}.\end{equation}$

Merton跳-扩散模型中, $\phi(1)=\gamma+\frac{1}{2}\sigma^2+\lambda({\rm e}^{\mu+\frac{1}{2}\delta^2}-1)$,风险中性特征函数为

(2.11)$\begin{eqnarray}\Phi(u)&=&{\bf E}^Q[{\rm e}^{{\rm i}uX_{t}}]\\&=&\exp\left\{t\left[{\rm i}\left(r-\frac{1}{2}\sigma^2-\lambda({\rm e}^{\mu+\frac{1}{2}\delta^2}-1)\right)u-\frac{1}{2}\sigma^{2}u^{2}+\lambda({\rm e}^{{\rm i}\mu u-\frac{1}{2}\delta^{2}u^{2}}-1)\right]\right\}.\end{eqnarray}$

如果密度函数$f_Q(s, T)$已知,在风险中性鞅测度${\mathbb{Q}}$下,敲定价格为$K$,到期日为$T$的欧式看涨期权在$t\, (0\leq t\leq T)$时刻的无套利价格可表示为到期日收益$(S_T-K)^+$折现价格的期望

$\begin{eqnarray*} C_t(S_t, t, K, T, \sigma)&=&{\rm e}^{-r(T-t)}{\bf E}^Q[(S_T-K)^+\mid S_t=s]\\ &=&{\rm e}^{-r(T-t)}\int_0^\infty f_Q(s, T)(s-K)^+{\rm d}s\\ &=&{\rm e}^{-r(T-t)}\int_K^\infty f_Q(s, T)(s-K){\rm d}s. \end{eqnarray*} $

引理3.1^[1]  在Black-Scholes模型下,欧式看涨期权在$t\, (0\leq t\leq T)$时刻的无套利价格可表示为

(3.1)$\begin{equation} C_t^{BS}(S_t, t, K, T, \sigma)=S_tN(d_1)-K{\rm e}^{-r(T-t)}N(d_2), \end{equation} $

其中$ N(x)=\frac{1}{\sqrt{2\pi}}\int_{-\infty}^x {\rm e}^{-\frac{u^2}{2}}{\rm d}u, \, d_1=\frac{\ln{\frac{S_t}{K}}+\left(r+\frac{1}{2}\sigma^2\right)(T-t)}{\sigma \sqrt{T-t}}, \, d_2=d_1-\sigma \sqrt{T-t}.$

引理3.2^[2]  在Merton跳-扩散模型下,欧式看涨期权在$t\, (0\leq t\leq T)$时刻的无套利价格可表示为

(3.2)$\begin{equation} C_t^M(S_t, t, K, T, \sigma)=\sum\limits_{k=0}^\infty \frac{(\lambda(T-t))^k}{k!}{\rm e}^{-\lambda(T-t)}C_t^{BS}(S_{k}, t, K, T;\sigma_k), \end{equation} $

其中$ S_k=S_t\exp\left[k\mu+\frac{k\delta^2}{2}-\lambda\left( \exp(\mu+\frac{\delta^2}{2})-1\right)(T-t)\right], \, \sigma_k=\sqrt{\sigma^2+\delta^2\frac{k}{T-t}}.$

在等价鞅测度${\mathbb{Q}}$下,欧式期权在$t\, (0\leq t\leq T)$时刻的风险中性价格可表示为

(3.3)$\begin{equation} v(x, t)={\rm e}^{-r(T-t)}{\bf E}^Q[v(y, T)|x]={\rm e}^{-r(T-t)}\int_{\bf R}v(y, T)f(y|x){\rm d}y, \end{equation}$

其中$x, y$为状态变量,表示股票价格$S_{t}$与敲定价格$K$商的对数

$x=\ln\frac{S_{t}}{K}, \ \ \y=\ln\frac{S_T}{K}.$

$v(y, T)=[\alpha K({\rm e}^y-1)]^+, \alpha=\left\{ \begin{array}{ll} 1, & \mbox{欧式看涨}, \\ -1, &\mbox{欧式看跌}. \end{array}\right. $

$f(y|x)$是$x$条件下$y$的概率密度函数, $r$为无风险利率, $T$为期权的到期日.

首先,把(3.3)式的积分区间截断为某个满足精度要求的有限区间$[a, b]$.由于密度函数在$y\rightarrow\pm\infty$时快速的衰减到0,因此截断为有限区间上的积分对定价结果不会有显著影响,这里假设$[a, b]\subset{\bf R}$,且

(3.4)$\begin{equation} \int_{{\bf R}\setminus[a, b]}f(y|x){\rm d}y<{\rm TOL}. \end{equation} $

对于给定的误差限TOL, (3.3)式可近似表示为

(3.5)$ \begin{equation} v_1(x, t)={\rm e}^{-r(T-t)}\int_a^bv(y, T)f(y|x){\rm d}y. \end{equation}$

第二步,对密度函数$f(y|x)$在$[a, b]$进行傅里叶-余弦展开,可得

(3.6)$\begin{equation}f(y|x)=\sum\nolimits_{k=0}^{'\infty}A_k(x)\cos\left(k\pi\frac{y-a}{b-a}\right), \end{equation}$

其中$\sum'$表示和式中第一项系数为$\displaystyle\frac{1}{2}$,系数$A_k(x)$可表示为

(3.7)$\begin{equation} A_k(x)=\frac{2}{b-a}\int_a^bf(y|x)\cos\left(k\pi\frac{y-a}{b-a}\right){\rm d}y. \end{equation}$

把(3.6)式代入到(3.5)式,并交换求和与积分运算次序,得

(3.8)$\begin{equation} v_1(x, t)=\frac{1}{2}(b-a){\rm e}^{-r(T-t)}\sum\nolimits_{k=0}^{'\infty} A_k(x)V_k, \end{equation} $

其中$V_k$为$\upsilon(y, T)$在$[a, b]$的傅里叶-余弦级数的系数,即

$V_k:=\frac{2}{b-a}\int_a^bv(y, T)\cos\bigg(k\pi\frac{y-a}{b-a}\bigg){\rm d}y.$

由于当$k\rightarrow\infty$时,系数快速减小,因此可以进一步用有限项(前N项)的和近似代替无穷级数的和,得

(3.9)$\begin{equation} v_2(x, t)=\frac{1}{2}(b-a){\rm e}^{-r(T-t)}\sum\nolimits_{k=0}^{'N-1} A_k(x)V_k. \end{equation}$

第三步,定义条件特征函数

$\Psi(u;x):=\int_{{\bf R}}{\rm e}^{{\rm i}ux}f(y|x){\rm d}y.$

那么,系数$A_k(x)$可表示为

(3.10)$\begin{equation}A_k(x)=\frac{2}{b-a}{\cal R}\left\{{\rm e}^{-{\rm i}k\pi\frac{a}{b-a}}\int_a^b {\rm e}^{{\rm i}\frac{k\pi }{b-a}y}f(y|x){\rm d}y\right\}, \end{equation}$

其中${\cal R}[\cdot]$表示求函数的实部算子,由(3.4)式可知(3.10)式中的有限积分可由下式来近似

(3.11)$\begin{equation}\int_a^b {\rm e}^{{\rm i}\frac{k\pi }{b-a}y}f(y|x){\rm d}y\approx\int_{{\bf R}} {\rm e}^{{\rm i}\frac{k\pi }{b-a}y}f(y|x){\rm d}y=\Psi\left(\frac{k\pi}{b-a};x\right).\end{equation}$

因此, $A_k(x)\approx F_k(x)$,其中

$F_k(x):=\frac{2}{b-a}{\cal R}\left\{\Psi\left(\frac{k\pi}{b-a};x\right){\rm e}^{-{\rm i}k\pi\frac{a}{b-a}}\right\}.$

由$F_k(x)$替换(3.8)式中的$A_k(x)$,可得欧式期权的COS定价公式

(3.12)$\begin{equation}\tilde{v}(x, t)={\rm e}^{-r(T-t)}\sum\nolimits_{k=0}^{'N-1}{\cal R}\left[\Psi\left(\frac{k\pi}{b-a};x\right){\rm e}^{-{\rm i}\frac{k\pi a}{b-a}}\right] V_k.\end{equation}$

在指数Lévy模型下,条件特征函数$\Psi(u;x)=\Phi(x){\rm e}^{{\rm i}ux}$,其中$\Phi(u):=\Psi(u;0)$,因此欧式期权的COS定价公式可表示为

(3.13)$\begin{equation}v(x, t)\approx {\rm e}^{-r(T-t)}\sum\nolimits_{k=0}^{'N-1}{\cal R}\left[\Phi\left(\frac{k\pi}{b-a}\right){\rm e}^{{\rm i} k\pi \frac{x-a}{b-a}}\right] V_k.\end{equation}$

COS方法积分截断范围

(3.14)$\begin{equation} [a, b]=\left[(c_1+x)-L\sqrt{c_2+\sqrt{c_4}}, (c_1+x)+L\sqrt{c_2+\sqrt{c_4}}\right], \, L=10, \end{equation} $

其中$x=\ln S_t/K$, $c_n$是$\ln(S_t/K)$的$n$阶累积量

$c_n=\frac{1}{{\rm i}^n}\frac{\partial^n(t\varphi(u))}{\partial u^n}\bigg|_{u=0}, t\varphi(u)=\ln\Phi(u).$

命题3.1   Merton跳-扩散模型中, $\ln({S_t/K})$在测度${\mathbb{Q}}$下的$n$阶累积量$c_n$及漂移修正项$\omega$分别为

$c_1=(r-\frac{1}{2}\sigma^2+\omega+\lambda\mu)t, \ \ \ c_2=(\sigma^2+\lambda\mu^2+\lambda\delta^2)t, $

$c_4=\lambda(\mu^4+6\delta^2\mu^2+3\mu^4\lambda)t, \ \ \ \omega=\lambda (\exp(\mu+\frac{1}{2}\delta^2)-1).$

其中$\omega$满足$\exp(-\omega t)=\Phi(-{\rm i})$.

记

(3.15)$\begin{equation}\epsilon(x, t;N, [a, b]):=v(x, t)-\tilde{v}(x, t;N, [a, b]), \end{equation}$

由$\tilde{v}(x, t)$近似代替$v(x, t)$误差来自三部分:

(1)由$v_1(x, t)$代替$v(x, t)$的积分范围截断误差

(3.16)$\begin{equation} \epsilon_1(x, t;[a, b]):=v(x, t)-v_1(x, t;[a, b])=\int_{{\bf R}\setminus[a, b]}\upsilon(y, T)f(y|x){\rm d}y. \end{equation} $

由于期权价格$v(y, T)$在$[a, b]$上有界,因此

$\epsilon_1(x, t;[a, b])\sim O\bigg(\int_{{\bf R}\setminus[a, b]}f(y|x)\bigg)\sim O({\rm TOL}).$

(2)由$ v_2(x, t)$代替$v_1(x, t)$时,由前$N$项和代替无穷级数的截断误差

(3.17) $ \begin{equation} \renewcommand\arraystretch{1.5} \begin{array}{ll}\displaystyle\epsilon_2(x, t;N, [a, b])&:=v_1(x, t;[a, b])-v_2(x, t;N, [a, b])\\&\displaystyle=\frac{1}{2}(b-a){\rm e}^{-r(T-t)}\sum\nolimits_{k=N}^{'\infty}A_k(x)\cdot V_k.\end{array}\end{equation} $

对于给定$x$,密度函数$f(x)\in C^{\infty}([a, b])$时, $\epsilon_2$呈指数式收敛

$|\epsilon_2|<M\cdot\exp(-(N-1)\upsilon), $

其中$\upsilon>0$为常数, $M$为常数,且小于以$N$为变量的指数函数;当概率密度函数具有不连续导数时,傅里叶余弦展开以代数形式收敛,即

$|\epsilon_2|<\frac{\bar{M}}{(N-1)^{\xi-1}}, $

其中$\bar{M}$为常数, $\xi\geq n\geq 1$($n$为$V_k$的代数收敛指数)^[7].

(3)由$\tilde{v}(x, t)$代替$v_2(x, t)$时, $A_k(x)$替换为$F_k(x)$的近似误差

(3.18) $ \begin{eqnarray}\epsilon_3(x, t;N, [a, b])&:=&v_2(x, t;N, [a, b])-\tilde{v}(x, t;N, [a, b])\\&\displaystyle=&{\rm e}^{-r(T-t)}\sum\nolimits_{k=0}^{'N-1}{\cal R}\left\{\int_{{\bf R}\setminus[a, b]}{\rm e}^{{\rm i}k\pi \frac{y-a}{b-a}}f(y|x){\rm d}y\right\}V_k.\end{eqnarray} $

$\epsilon_3$由积分被截断的误差,可由下式界定

$|\epsilon_3|<|\epsilon_1|+D\bigg|\int_{{\bf R}\setminus[a, b]}f(y|x){\rm d}y\bigg|, $

其中$D$为不依赖于$N$的某个常数.

由三部分误差限及三角不等式,可得

(3.19) $ \begin{eqnarray}\displaystyle\epsilon(x;N, [a, b])&\displaystyle\leq&|v-v_1|+|v_1-v_2|+|v_2-\tilde{v}|\\&\displaystyle\leq&\bar{D}\cdot\bigg|\int_{{\bf R}\setminus[a, b]}f(y|x){\rm d}y\bigg|+|\epsilon_2(x;N)|, \end{eqnarray} $

其中$\bar{D}$是不依赖于$[a, b]$和$N$的常数.当积分区间$[a, b]$足够大时,级数的截断误差$\epsilon_2(N)$是误差的主要部分,此时,对于平滑的密度函数而言, $\epsilon$的收敛速度呈指数形式(更多详细的证明可参看文献[7]).

下面通过数值实验将COS方法的数值解与Merton模型下欧式看涨期权的精确解进行比较,参数选取如下:

$\begin{eqnarray*} S_0=100, r=0.1, q=0, T=1, \, \sigma=0.2, \mu=0.5, \delta=0.1, \lambda=1. \end{eqnarray*} $

$K=90, 100, 110$,由公式(3.2)计算期权的精确价格分别为$31.4863104397\cdots, $$ 27.6088525266$$\cdots, $$ 24.3602541845\cdots . $

由表 1可以看出,随着截取数列项数$N$的增大COS方法定价结果与精确解的误差迅速减小,当$N=256$时, COS方法的数值结果与精确解的误差下降为约3.5e-9. COS方法收敛速度非常快,而且程序运行时间也比较短,因此COS方法是一种有效的定价方法.

参数校准问题与给定模型参数进行期权定价过程是两个互逆的问题,参数校准是通过观测到的期权市场价格求解模型中的未知参数,即寻找合适的参数$\theta^Q$,使得模型价格与市场观测的期权价格一致.

问题1  对于给定的一组到期日为$T_{i}\, (i=1, 2, \cdot\cdot\cdot, n)$,执行价格为$K_j\, (j=1, 2, \cdot\cdot\cdot m)$的欧式看涨期权的市场观测价格$\{C_{ij}\}$,寻找一个风险中性模型参数向量$\theta^Q$,使得

$C_{ij}= C^{\theta^Q}(T_i, K_j)={\rm e}^{-rt}{\bf E}^Q[(S_{T_{i}}-K_{j})^+|S_{t}=S], $

其中$C^{\theta^Q}(T_i, K_j)$为含有参数$\theta^Q$的模型价格, Merton跳-扩散模型中$\theta^Q=\theta(\sigma^Q, \lambda^Q, \mu^Q, \delta^Q)$.

由于市场价格存在一定的噪声,市场观测价格与模型价格不可能精确匹配,因此校准问题转化为两种价格的最佳近似问题,通常可采用非线性最小二乘方法求解.

利用非线性最小二乘(NLS)方法求解问题1,即对下列函数最小化

(4.1)$\begin{equation} \theta^{Q*}={\rm arg}\min\limits_{\theta^Q}\sum\limits_{i=1}^{n}\sum\limits_{j=1}^{m}\omega_{ij}\left|C^{\theta^Q}(T_{i}, K_{j})-C_{ij}\right|^{2}, \end{equation} $

其中$C_{ij}^{{\rm bid}}\leq C_{ij}\leq C_{ij}^{{\rm ask}} $, $C^{\theta^Q}(T_i, K_j)$为Merton模型下计算得到的风险中性的看涨期权价格, $\omega_{ij}$为权重,通常情况下权重$\omega_{ij}$有两种选法

(4.2)$\begin{equation} \omega_{ij}=\frac{1}{\left|C_{ij}^{{\rm bid}}-C_{ij}^{{\rm ask}}\right|^{2}}, i=1, 2, \cdots , n;\, j=1, 2, \cdots , m. \end{equation}$

(4.3)$ \begin{equation} \omega_{ij}=\frac{1}{({\rm vega}_{ij})^{2}}=\frac{1}{\left|K_{j}{\rm e}^{-rT_{i}}N(d_{2})\sqrt{T_{i}}\right|^{2}}, i=1, 2, \cdots , n;\, j=1, 2, \cdots , m. \end{equation}$

图 1给出了当固定参数$\sigma$和$\delta$时, Merton跳-扩散模型期权价格的定价误差曲面,期权价格由COS方法计算得出.从图中可以看出,非线性最小二乘函数目标函数是非凸的,这种情况下解对初值比较敏感,因此无法满足解的唯一性和稳定性,即非线性最小二乘方法是不适定的.因此需要增加罚函数项来增加目标函数的凸性,从而确定目标函数的唯一最优解,即对校准问题进行正则化.

本文选择相对熵(Kullback-Leibler距离) $\varepsilon({\mathbb{Q}}|{\Bbb P})$作为惩罚项^[8].因为相对熵具有良好的性质:

(1)相对熵函数是凸函数;

(2)相对熵函数保证了校准测度${\Bbb Q}$关于先验测度${\mathbb{P}}$的绝对连续性;

(3)如果先验测度是风险中性测度,相对熵函数保证了校准测度${\Bbb Q}$是${\mathbb{P}}$的等价鞅测度;

(4)相对熵函数易于计算.

定义4.1^[10]  令${\mathbb{P}}$和${\mathbb{Q}}$是$(\Omega, {\cal F})$上的两个等价概率测度. ${\mathbb{Q}}$关于${\mathbb{P}}$的相对熵定义为

(4.4)$\begin{equation}\varepsilon({\mathbb{Q}}|{\mathbb{P}})={\bf E}^Q\bigg[\ln\frac{{\rm d}{\mathbb{Q}}}{{\rm d}{\mathbb{P}}}\bigg]\={\bf E}^P\bigg[\frac{{\rm d}{\mathbb{Q}}}{{\rm d}{\mathbb{P}}}\ln\frac{{\rm d}{\mathbb{Q}}}{{\rm d}{\mathbb{P}}}\bigg].\end{equation}$

令$f(x)=x\ln x$,则

$\varepsilon({\mathbb{Q}}|{\mathbb{P}})={\bf E}^P\bigg[f(\frac{{\rm d}{\mathbb{Q}}}{{\rm d}{\mathbb{P}}})\bigg].$

由$f(x)=x\ln x$为严格凸函数,可知相对熵$\varepsilon({\mathbb{Q}}|{\mathbb{P}})$是关于${\mathbb{Q}}$的凸函数,并且$\varepsilon({\mathbb{Q}}|{\mathbb{P}})\geq 0$,当且仅当${\mathbb{Q}}={\mathbb{P}}$时等号成立.下面给出Merton跳-扩散模型的相对熵.

引理4.1^[10]  假设${\mathbb{P}}$和${\mathbb{Q}}$是指数Lévy模型下的等价概率测度,特征三元组分别为$(\sigma^{2}, \nu^P, \gamma^P)$和$(\sigma^{2}, \nu^Q, \gamma^Q)$,且$\sigma>0$.那么,在风险中性条件下, ${\mathbb{P}}$和${\mathbb{Q}}$相对熵函数可表示为

(4.5)$\begin{equation} \varepsilon({\mathbb{Q}}|{\mathbb{P}})=\frac{T}{2\sigma^{2}}\left[\int_{-\infty}^{+\infty}({\rm e}^{x}-1)(\nu^Q-\nu^P)({\rm d}x)\right]^2+ T\int_{-\infty}^{+\infty}(\frac{{\rm d}\nu^Q}{{\rm d}\nu^P}\ln(\frac{{\rm d}\nu^Q}{{\rm d}\nu^P}) +1-\frac{{\rm d}\nu^Q}{{\rm d}\nu^P})\nu^P({\rm d}x). \end{equation} $

定理4.1  假设$X_t$是由(2.5)式给出的Merton跳-扩散模型,特征三元组为$(\sigma^2, \nu^P, 0)$. ${\mathbb{Q}}$是与${\mathbb{P}}$等价的风险中性测度,其特征三元组$(\sigma^2, \nu^Q, 0)$,其中

$\nu^P(x)=\frac{\lambda^P}{\sqrt{2\pi}\delta^P}\exp\left(-\frac{(x-\mu^P)^2}{2\delta^{P^2}}\right){\rm d}x, \quad \nu^Q(x)=\frac{\lambda^Q}{\sqrt{2\pi}\delta^Q}\exp\left(-\frac{(x-\mu^Q)^2}{2\delta^{Q^2}}\right){\rm d}x.$

那么, Merton跳-扩散模型的相对熵函数为

(4.6) $ \begin{eqnarray}\varepsilon({\mathbb{Q}}|{\mathbb{P}})&=&\displaystyle\frac{T}{2\sigma^{2}}\left[\lambda^Q\left({\rm e}^{\mu^Q+\frac{1}{2}{\delta^Q}^2}-1\right)-\lambda^P\left({\rm e}^{\mu^P+\frac{1}{2}{\delta^P}^2}-1\right)\right]^2 +\displaystyle T\lambda^Q\ln\left(\frac{\lambda^Q\delta^P}{\lambda^P\delta^Q}\right)\\ &&+T\lambda^P+T\lambda^Q\left(-\displaystyle\frac{3}{2}+\frac{{\delta^Q}^2+\left(\mu^Q-\mu^P\right)^2}{2{\delta^P}^2}\right). \end{eqnarray} $

证  由引理4.1可知

(4.7)$\begin{eqnarray}\displaystyle\frac{1}{T}\varepsilon({\mathbb{Q}}|{\mathbb{P}})&=&\displaystyle\frac{1}{2\sigma^2}\left[\int_{-\infty}^{+\infty}({\rm e}^x-1)(v^Q-v^P)\right]^2 +\displaystyle\int_{-\infty}^{+\infty}v^P{\rm d}x-\displaystyle\int_{-\infty}^{+\infty}v^Q{\rm d}x+\int_{-\infty}^{+\infty}\ln(\frac{{\rm d}v^Q}{{\rm d}v^P}){\rm d}v^Q\\ &=&\displaystyle\frac{1}{2\sigma^2}[\lambda^Q{\bf E}^{Q}({\rm e}^x-1)-\displaystyle\lambda^P{\bf E}^{P}({\rm e}^x-1)]^2+\lambda^P-\lambda^Q\\&&+\displaystyle\int_{-\infty}^{+\infty}\left(\displaystyle\ln(\frac{\lambda^Q\delta^P}{\lambda^P\delta^Q})-\frac{(x-\mu^Q)^2}{2{\delta^Q}^2} +\frac{(x-\mu^P)^2}{2{\delta^P}^2}\right){\rm d}v^Q\\ &=&\displaystyle\frac{T}{2\sigma^{2}}\left[\lambda^Q\left({\rm e}^{\mu^Q+\frac{1}{2}{\delta^Q}^2}-1\right)-\lambda^P\left({\rm e}^{\mu^P+\frac{1}{2}{\delta^P}^2}-1\right)\right]^2 +\lambda^P-\lambda^Q\\&&+\displaystyle\lambda^Q\ln\left(\frac{\lambda^Q\delta^P}{\lambda^P\delta^Q}\right)-\displaystyle\int_{-\infty}^{+\infty}\frac{(x-\mu^Q)^2}{2{\delta^Q}^2}{\rm d}v^Q +\displaystyle\int_{-\infty}^{+\infty}\frac{(x-\mu^P)^2}{2{\delta^P}^2}{\rm d}v^Q.\end{eqnarray}$

而

(4.8)$\begin{eqnarray}\displaystyle\int_{-\infty}^{+\infty}(x-\mu^Q)^2{\rm d}v^Q&=&\displaystyle\int_{-\infty}^{+\infty}(x^2-2\mu^Qx+{\mu^Q}^2)^2{\rm d}v^Q\\&=&\displaystyle\lambda\left({\bf E}^{Q}(x^2)-2\mu^Q{\bf E}^{Q}(x)+{\mu^Q}^2\right)\\&=&\displaystyle\lambda\left({\delta^Q}^2+{\mu^Q}^2-2{\mu^Q}^2+{\mu^Q}^2\right)\\&=&\displaystyle\lambda^Q{\delta^Q}^2, \end{eqnarray}$

(4.9)$\begin{eqnarray}\displaystyle\int_{-\infty}^{+\infty}(x-\mu^P)^2{\rm d}v^Q&=&\displaystyle\int_{-\infty}^{+\infty}(x^2-2\mu^Px+{\mu^P}^2)^2{\rm d}v^Q\\&=&\displaystyle\lambda\left({\bf E}^{Q}(x^2)-2\mu^P{\bf E}^{Q}(x)+{\mu^P}^2\right)\\&=&\displaystyle\lambda\left({\delta^Q}^2+{\mu^Q}^2-2\mu^Q\mu^P+{\mu^P}^2\right)\\&=&\displaystyle\lambda^Q\left({\delta^Q}^2+(\mu^Q-\mu^P)^2\right).\end{eqnarray}$

把(4.8)和(4.9)式代入到(4.7)式可知定理4.1成立.

增加相对熵作为惩罚函数,问题1可转化对如下函数进行最小化

(4.10)$\begin{equation}J(\sigma^Q, \lambda^Q, \mu^Q, \delta^Q):=\sum\limits_{i=1}^{N}\sum\limits_{j=1}^{M_i}\omega_{ij}\left|C^{\theta^Q}(T_{i}, K_{j})-C_{ij}\right|^{2}+\alpha \varepsilon({\mathbb{Q}}|{\mathbb{P}}), \end{equation}$

其中$\varepsilon({\mathbb{Q}}|{\mathbb{P}})$是惩罚函数项, $\alpha$是正则化参数, $\omega_{ij}$是权重.

先验参数向量$\theta^P$和正则化参数$\alpha$对校准结果非常重要,因此必须谨慎选取.先验测度的选取可采用Cont和Tankov^[8-9]提出的通过求解非正则化的NLS校准问题得到风险中性测度$\theta^{Q*}$来充当先验参数.正则化参数对于校准结果的准确性和稳定性都非常重要,正则化参数$\alpha$与数据的扰动水平有关,因此不能由先验的一些固定数据来确定.本文应用Morozov偏差原则来决定$\alpha$,步骤如下:

首先,计算非线性方差问题(4.1)式得到扩散估计$\hat{\sigma}_{\alpha=0}$.

然后,固定$\sigma=\hat{\sigma}_{\alpha=0}$,利用梯度下降法计算非线性方差问题(4.1)得到模型内在的二次定价误差

$\hat{\varepsilon}_0^{2}=\inf\limits_{\nu^Q}\sum\limits_{i=1}^n\sum\limits_{j=1}^m\omega_i\left|C_{ij}^{\nu^Q}-C_{ij}\right|^2.$

接下来, $\varepsilon^{2}(\alpha)=\hat{\varepsilon}_0^{2}+\alpha \varepsilon({\mathbb{Q}}|{\mathbb{P}})$.由于相对熵函数项的增加, $\varepsilon^{2}(\alpha)>\hat{\varepsilon}_0^{2}$.在与模型误差相同数量级时, Morozov原则允许通过下降校准精确度来换取校准结果的唯一性和稳定性.因此可利用梯度下降算法求解下面的方程来确定$\alpha$的值

(4.11)$\begin{equation} \varepsilon^2(\alpha^*)=c\hat{\varepsilon}_0^{2} , (c>1). \end{equation}$

本节对前面给出的校准方法进行数值实验.首先,在给定参数条件下,利用COS方法生成Merton模型下期权价格的模拟数据.然后,利用Matlab R2012b,调用fmincon函数对目标函数进行优化求解.

假设市场服从Merton模型,特征三元组为$(\sigma, \nu, 0)$,其中$\nu=\nu(\mu, \delta, \lambda)$,取$\sigma=0.2, $$\mu=0.05, \delta=0.1, \lambda=1$,并假设股票在$t=0$时刻的价格为$S_0=100$,到期日$T=1$,敲定价格取$K\in[70,120]$,等间距取100个点.事实上,期权的真实价格与Merton模型计算价格$C_{ij}^*$是有差异的,因此在Merton期权价格的计算上加一项噪声来模拟真正的市场数据$C_{ij}$,即

$C_{ij}=C_{ij}^*(1+\eta z), $

其中$z$是一个标准正态分布的随机数,这里取$\eta=0.03$.

为了计算简单,固定参数$\sigma^Q=\sigma=0.2$和$\delta^Q=\delta=0.1$,对参数$\lambda^Q$和$\mu^Q$进行校准.

首先,利用NLS方法求解(4.1).运算结果表明:当初值取$\lambda_0=1.2, \, \mu_0=-1$时,得最优解为$\lambda^Q=1.7991, \, \mu^Q=-0.0198$;当初值取$\lambda_0=0.6, \, \mu_0=0.4$时,得最优解为$\lambda^Q=0.3084, \, \mu^Q=0.1574$.这种情况说明了由NLS方法得到的解是不唯一的、不稳定的,即非线性最小二乘校准方法是不适定的.

接下来,利用相对熵方法进行正则化校准.数值计算主要包括以下四步:

(1)选择权重: $\omega_{ij}=\frac{1}{({\rm vega}_{ij})^{2}}$;

(2)先验参数选取: $\delta^P=\delta=0.1, \, \sigma^P=\sigma=0.2, \mu^P=0.046, \, \lambda^P=0.98$;

(3)选择正则化参数: $\alpha=0.08$;

(4)在给定的$\alpha$和${\mathbb{P}}$下,对(4.10)进行最小化.通过运算得参数$\lambda^{{\rm Q}}=0.9917, \mu^{{\rm Q}}=0.0461$,且结果不依赖于给定的初值,因此校准结果是稳定的.

图 4给出了真实Lévy测度与校准测度,通过比较可以发现校准所得Lévy测度与真实的Lévy测度曲线是非常接近的,校准结果具有较高的精度.因此相对熵正则化方法是一种可行、有效的方法.

本节利用相对熵正则化方法对市场数据进行实证分析.数据选自2015年4月3日S&P500股指期权,当日指数收盘价格$S_0=2066.96$,分别选取了到期日为4月17日($T_1=0.038$)、5月15日($T_2=0.115$)、6月19日($T_3=0.211$)、9月18日($T_4=0.460$)的交易量非零的看涨期权数据,敲定价格$K$的选取范围为1700到2300.根据美国国债计算的无风险利率为$r=0.001$.

含有未知参数的Merton模型的期权价格$C^{\theta^Q}(T_i, K_j)$由COS方法得到.数值计算具体步骤如下:

(1)选择权重:选用市场买卖差价利用(4.2)计算;

(2)先验参数选取:利用非正则化的NLS方法对(4.1)最小化的解作为先验参数;

(3)选择正则化参数$\alpha$:令$c=1.2$,利用市场数据通过求解(4.11)得到;

(4)在给定的$\alpha$和${\mathbb{P}}$条件下,对函数$J$进行最小化求解.

我们通过下列根均值误差(RMSE)来检验校准的质量

(5.1)$\begin{equation}{\rm RMSE}=\sqrt{\sum\limits_{i=1}^n\sum\limits_{j=1}^{m}\frac{(C^{\theta}(T_i, K_j)-C_{ij})^2}{mn}}.\end{equation}$

表 2给出了利用不同到期日的S&P500欧式看涨期权市场数据分别对Black-Scholes模型和Merton模型的参数校准结果.可以看出Merton跳-扩散模型的RMSE值小于Black-Scholes模型的RMSE的值,这意味着Merton模型比Black-Scholes模型能更好的拟合市场数据.这主要反映在自由度上,具有更多自由度的模型比具有较小自由度的模型能更好的模拟市场数据.另一方面,可以看出对于Merton模型,参数值在不同的到期日得到的结果是不同的,并且较长的期限,平均跳跃幅度及跳跃强度都随持有期间长度而减少.

图 5给出了利用校准参数所得欧式期权价格与市场期权价格的比较,可以看出Merton跳-扩散模型下期权的价格与S&P500股指期权的市场数据是非常接近的,而由Black-Scholes模型计算的期权的价格与S&P500市场数据有较大的差异,这进一步说明了Merton跳-扩散模型比Black-Scholes模型能更好的模拟S&P500市场.

图 6给出了四个不同到期日的校准的隐含波动率与市场真实波动率的比较.由图可以看出Merton跳-扩散模型有效地捕捉了S&P500指数期权的隐含波动率曲线特征,校准结果是比较满意的.但当距离到期日时间较短时,校准结果相对较差,主要是因为临近到期日期权交易频繁,波动率变化较大,特征很难被有效的捕捉.而距到期日在一个月以上的,波动率变化较小,校准结果比较理想.

图 7给出了Merton跳-扩散模型下重构的隐含波动率曲面以及市场真实波动率曲面,可以看出Merton跳-扩散模型重现了隐含波动率的"微笑"和"偏斜"的特征,并且与市场真实波动率曲面拟合效果较好.

本文主要研究了跳-扩散模型下的期权定价问题及参数校准问题.首先,给出了Merton跳-扩散模型下期权定价的COS方法,分析了COS方法的定价误差,并通过数值实验验证了COS方法的有效性.其次,对Merton跳-扩散模型参数进行校准,结果表明非线性最小二乘方法得到的校准结果是不稳定的.采用相对熵方法进行正则化可以增加目标函数的凸性,从而可以得到校准问题的唯一性.然后,通过数值模拟实验验证了相对熵正则化方法的准确性和有效性.最后,利用相对熵正则化方法对S&P500股指期权数据进行参数校准,并利用校准参数重构了Merton跳-扩散模型下隐含波动率的曲面,得到了比较满意的校准结果.