本文主要利用MATLAB实现Kim和Wee(2014)得出的对Geometric Fixed-strike亚式看涨期权估值模型。
Kim和Wee的模型发展于Heston的随机波动率模型之上，同样对股价和波动率进行建模：
![]()
![]()
![]()
其中，W(t)和Z(t)是非独立的维纳过程，r是利率，kappa、theta及sigma为正常数，满足 ![]() 。
具体公式推导本文不做详细展示，详见Kim和Wee在2014年的文章“Pricing of geometric Asian options under Heston's stochastic volatility model"
公式的显示解为文章内定理4.1,其中xi ![]() 为变量，其他符号为参数
![]()
其中，RE为Real函数，对计算结果取实数部分。其他展开函数同具体代码一同展示,
显示解中，变量w均为0，可以依据此省去公式中w乘积部分。
显示解函数具体代码如下：
function call=KimCallQuad(nu,theta,s0,r,T,t,sigma,rho,kappa,n,K)
warning off;
call=exp(-r.*(T-t)).*KimB(nu,theta,s0,r,T,t,sigma,rho,kappa,n,K);
end
function ret=KimB(nu,theta,s0,r,T,t,sigma,rho,kappa,n,K)
ret=(Kimpsi(1,nu,theta,s0,r,T,t,sigma,rho,kappa,n)-K)./2+...
(1/pi).*integral(@(x)arrayfun(@(xi)KimIntegral(xi,nu,theta,s0,r,T,t,sigma,rho,kappa,n,K),x),0,10^5);
%利用arrafun将函数
%变为支持向量输入，再用Integral代码解积分
end
function ret=KimIntegral(xi,nu,theta,s0,r,T,t,sigma,rho,kappa,n,K)
ret=real((Kimpsi(1+1i*xi,nu,theta,s0,r,T,t,sigma,rho,kappa,n)-...
Kimpsi(1i*xi,nu,theta,s0,r,T,t,sigma,rho,kappa,n).*K).*...
(exp(-1i*xi.*log(K)))./(1i*xi));
endKim psi函数具体代码如下：
![]()
Kim psi
function psi=Kimpsi(s,nu,theta,s0,r,T,t,sigma,rho,kappa,n)
a1=2*nu./(sigma.^2);
a2=2*kappa.*theta./(sigma.^2);
a3=((T-t)./T).*log(s0)+((r.*sigma-kappa.*theta.*rho).*(T-t).^2)/(2*sigma.*T)-(rho.*(T-t)/(sigma.*T)).*nu;
a5=(kappa.*nu+(kappa.^2)*theta.*(T-t))./(sigma.^2);
psi=exp(-a1.*(KimH(s,T,t,sigma,rho,kappa,n,2)./KimH(s,T,t,sigma,rho,kappa,n,1))-...
a2.*log(KimH(s,T,t,sigma,rho,kappa,n,1))+...
s.*a3+a5); % 1 2指定H函数类型
endKim H 函数具体代码如下：
![]()
Kim H
function H=KimH(s,T,t,sigma,rho,kappa,n,type)
h=ones(1,n+4);
h(1,1)=0; % 对应h-2
h(1,2)=0; % 对应h-1
h(1,3)=1; % 对应h0
h(1,4)=((T-t).*kappa)./2; % 对应h1
for j=5:n+4 % 大于等于2的h函数
h(1,j)=((T-t).^2/ (4*(j-3)*(j-4).*T.^2))*(-1.*(s.^2)*(sigma.^2).*(1-(rho.^2)).*((T-t).^2).*h(1,j-4)+...
s.*sigma*T.*(sigma-2*rho*kappa).*(T-t).*h(1,j-3)+...
T.*((kappa.^2).*T-2*s.*rho.*sigma).*h(1,j-2));
h1(1,j)=((j-3)./(T-t))*h(1,j);
end
if type==1
H=sum(h(1,3:n+4)) %H函数
else
H=sum(h1(1,5:n+4))+(1/(T-t))*h(1,4) %H波浪线函数
end给定参数后，利用MATLAB实现公式：S0=100， v0=0.09, t=0, r=0.05, theta=0.348, sigma=0.39, kappa=1.15, rho=-0.64, n=(10,20,30), T=(0.5,1,2), K=(90,100,110), 且无穷积分上限取10^5。
T=[0.5,1,2];
t=0;
sigma=0.39;
kappa=1.15;
rho=-0.64;
n=(10:1:38);
w=0;
nu=0.09;
r=0.05;
theta=0.348;
s0=100;
K=[90,100,110];
for i=1:29
for j=1:29
for k=1:29
t1=clock;
C(i,j,k)=KimCallQuad(nu,theta,s0,r,T(i),t,sigma,rho,kappa,n(j),K(k))
t2=clock;
te(i,j,k)=etime(t2,t1) %计算运算时间
end
end
end计算结果如下，基本复现论文计算结果：
![]()
计算结果 |
|
转播
