拓端tecdat|R语言广义二次跳跃、非线性跳跃扩散过程转移函数密度的估计及其应用

原文链接：http://tecdat.cn/?p=23551 

原文出处：拓端数据部落公众号

什么是跳跃扩散(Jump Diffusion)？

跳跃扩散模型是一种用来对期权合约进行估价或定价的模型，它混合了两种定价技术：一种是更传统的扩散模型，在这种模型中，因素以平稳和相对一致的方式发挥作用；另一种是跳跃过程模型，在这种模型中，一次性事件会引起重大变化。

理论上，跳跃扩散就是这样产生的跳转扩散是一种用来对期权合约进行估价或定价的模型期权定价是在期权合约上设定一个客观价值的技巧，这是一个交易者在未来某个日期以固定价格购买完成资产出售或购买的权利的一种金融协议，不同的模型试图计算出影响合约对持有人价值的不同因素。这些因素可能包括标的资产的当前价格、资产价格的波动性以及期权到期前的剩余时间。

许多交易者会使用定价模型决定他们可以为一个期权支付什么价格，并在他们从期权中赚的钱和不值得行使期权从而浪费购买价格的风险之间获得一个很好的价值平衡最常见的期权定价形式可以描述为基于扩散的定价，其原理是市场事件对资产价格的影响相对较小，总体趋势和模式仍将持续最著名的基于扩散的期权定价形式是Black-Scholes模型，它的主要优点是这样一个模型可以相对简单和直接地操作，一种对比的模型被称为跳跃过程基于这样一个基础，即市场不会一直朝着一个总体平稳的方向移动，且偏差较小，而是更容易受到一次性事件中方向和节奏的巨大变化的影响。

使用跳跃过程的模型，例如二项期权定价模型，尝试更多地考虑不可预测事件的可能性。这使得模型更加复杂，尽管期权到期前剩余的时间越短，布莱克-斯科尔斯（Black-Scholes）估值和二项式期权估值之间的差距就越小经济学家罗伯特C.默顿将这两种模型结合起来，具体称为默顿模型，通常称为跳跃扩散模型。它试图涵盖这样一种观点，即市场具有总体趋势、微小的日常变化和重大冲击的组合关于跳跃扩散的研究后来被纳入了一个改编的布莱克-斯科尔斯模型中，该模型在1997年获得了新的经济学。

跳跃扩散过程引入

扩散的意思是过程可以包含布朗运动的成份，甚至写成布朗运动的积分形式。同时，这样的过程可以包含跳跃。

最基础的跳过程是泊松过程。复合泊松过程也是一个泊松过程，但跳跃大小是随机的。

引入纯跳跃过程的概念。也就是一个初始为零，并且在时间区间上跳跃由有限多次，并且没有跳跃时保持常数的随机过程。

二次跳跃扩散的矩方程

二次跳跃扩散的矩方程是围绕着一种近似于Kolmogorov正向方程解的算法建立的（Hanson 2007）。

其中ϕ(.)给出了跳跃分布，元素∇(Xt+J(Xt,z˙t))和|δ(z˙t)|有特殊含义。这里，∇(Xt+J(Xt,z˙t))的作用是对跳跃J(Xt,z˙t)进行映射，达到作为跳跃的结果状态Xt。

例如，如果J(Xt,z˙t)=z˙t，那么∇(Xt+J(Xt,z˙t))=Xt-z˙t。

因此，函数∇(.)将过程的状态恢复到跳跃发生之前的状态。另一方面，元素|δ(z˙t)|只是对应于函数∇(.)的雅各布系数。例如，如果J(Xt,z˙t)=z˙t，∇(Xt+J(Xt,z˙t))=Xt-z˙t，那么|δ(z˙t)|=1。该近似方法是基于计算扩散过程随时间变化的条件矩轨迹，随后将这些矩带入一个合适的密度中。然后，代理密度被用来近似计算Kolmogorov方程的解。

例如，给定一个标量的广义二次元扩散。

其中

为J(Xt,z˙t)对Xt的多项式，其中计数过程Nt的强度形式为

有可能推导出一个微分方程组，该方程组控制着该过程的时刻随时间的演变。

让mi(t)表示Xt的第i个累积值，考虑到该过程具有初始值Xs。然后，对于J(Xt,z˙t)=z˙t

受制于初始条件mi(s)=Xis for s<t 和 i=1,2,...，其中ui表示随机变量z˙t的第i个非中心时刻。通过标准Runge-Kutta方案对这些方程进行数值求解，并将得到的数值插入代用密度，如正态分布，就有可能精确地近似转移密度。

转移密度的过度因子化

由于跳跃扩散过程的二分法性质，在试图通过将其时刻插入到一个合适的代理密度中来近似转移密度时出现了一个问题。在短的转移期内，跳跃动力学有可能创造出与纯粹的扩散过程大不相同的密度特征。例如，在一个足够短的转移水平线上，一个纯粹的扩散过程是近似正常分布的，而跳跃扩散（取决于跳跃机制的性质）可能有厚尾，在过程的初始值周围有一个非常尖峰的分布。因此，简单地将跳跃扩散的时刻轨迹插入代用密度以近似转移密度并不总是足够的。为了这些目的，我们将转移密度计算为一个混合密度，将该过程在其纯粹的扩散动力学和超额分布方面进行因子化，这说明了跳跃机制的影响。也就是说，我们计算。

其中fJ(Xt|Xs)表示跳跃扩散的转移密度，fD(Xt|Xs)表示跳跃扩散的无跳跃对应物，fE(Xt|Xs)表示超额分布，P(Nt-Ns=0)被计算为方程的解。

为了近似转移密度，计算fJ(Xt|Xs)和fD(Xt|Xs)的矩，然后可以推导出fE(Xt|Xs)的矩。随后，fD(Xt|Xs)和fE(Xt|Xs)的矩与合适的代理密度一起使用，以近似转移密度。

---

非线性跳跃扩散

基本仿射跳扩散过程

金融学中经常使用的非线性跳跃扩散模型的一个例子，被称为 "基本仿射跳扩散过程"（BAJD）（Eckner 2009），它在流行的CIR过程的基础上进行了扩展，假设该过程的路径经历了指数分布的跳跃，其强度不变。因此，相应的SDE是由系统给出的。

其中Nt有一个恒定的强度，λ(Xt,t)=λ0，z˙t∼Exp(ν)。

按照跳跃扩散的GQD框架，我们可以用R代码在工作空间内定义BAJD。

1.  
    
2.  
    
3.  
   # 定义模型。
4.  
    
5.  
   # 扩散部分
6.  
   G0 <- function(t){a*b}.
7.  
   G1 <- function(t){-a}。
8.  
   Q1 <- function(t){sigma}.
9.  
   # 跳跃部分
10.  
    Lam0 <- function(t){lam_0}
11.  
    Jlam <- function(t){nu}.
12.  
     
13.  
     

随后，我们可以使用density()函数对转移密度进行近似。

density(Xs,Xt,s,t,dt, 'Exponential')

 persp(x=Xt,y=time,z=densit)

有趣的是，BAJD的转移密度在表面上似乎与CIR过程的转移密度非常相似。然而，正如我们在下面的例子中所看到的，跳跃机制的存在极大地影响了转移密度的性质。

具有随状态变化的跳跃强度的时间不均匀CIR

在用扩散过程对现实世界的现象进行建模时，一个常见的陷阱是构建模型的假设与数据中观察到的假设之间的差异，例如Black-Scholes方法中的恒定波动率假设。然而，为了使一个给定的理论在数学上具有可操作性，往往需要这样的假设。在跳跃扩散模型的背景下，更现实的情况是，数据中跳跃的强度与任何时候的过程水平有关。例如，考虑一个由SDE支配的跳跃扩散模型。

其中，Nt具有依赖于状态的强度λ(Xt,t)=λ0+λ1Xt，z˙t∼Nu(μz,σ2z)。这里我们假设跳跃到达是线性地依赖于过程的状态的。

在GQD框架内，我们用其系数函数来定义该模型。

1.  
    
2.  
   # 一些参数值。
3.  
   kap <- 1; beta_0 <- 5; beta_1 <- 3; sigma <- 0.15;
4.  
    
5.  
   # 界定模型。
6.  
   remove()
7.  
    

1.  
    
2.  
   # 扩散部分
3.  
   G0 <- function(t){kap*(beta_0+beta_1*sin(2*pi*t)) }
4.  
   G1 <- function(t){-kap}。
5.  
   # 跳跃部分
6.  
   Lam0 <- function(t){lam_0}
7.  
   Jsig <- function(t) {sigma_z}
8.  
    

再次注意添加了跳跃分布Jmu(t)和Jsig(t)的系数函数，这里反映了正态跳跃密度的参数。

随后，通过向density()函数提供初始值，我们可以计算出转移密度近似值。

density(Xs,Xt,s,t,dt,  'Normal')

 persp(Xt,time,density)

请注意，由于转移期的范围很大，我们没有指定需要使用因式分解。如果我们希望调用因式分解，我们可以通过因式分解参数来指定。

1.  
   # 绘制t =0.2时的转移密度和零点的演变。
2.  
   # 跳跃概率
3.  
   plot(density~Xt)
4.  
    
5.  
   # 在概率轨迹上叠加短水平线。
6.  
   abline(v=t,lty='dotted')

注意由于跳跃机制的存在，转移密度的偏斜形状。在短的转移期内，跳跃机制对转移密度的性质有很大影响。更多的时候，转移密度将表现出重尾。因此，在使用矩截断方法对短转移期的转移密度进行近似时，必须使用因子化。例如，将跳跃扩散的转移密度与它的扩散进行比较。

1.  
   # 适用于纯扩散型GQD的软件包
2.  
   density(Xs,Xt,s,t,dt)

1.  
   plot(density[,20]~Xt)
2.  
   lines(density~Xt)

---

随机的跳跃强度

使用GQD框架，有可能构建复杂的时间不均匀结构，甚至在跳跃机制内。例如，有可能构建一个跳跃机制，其强度过程受一些外部随机过程的支配。例如，再次考虑一个具有加性跳跃的CIR过程，但在这种情况下，让强度参数成为一个双状态连续时间马尔科夫链（CTMC）。就是说。

强度参数y˙t的动态变化是由连续时间马尔科夫链给出的，Z˙∼N(θ5,θ26)和Nt∼PoiP(y˙t)。

与转移矩阵:

对于这个例子，我们让θ={2,5,2,1,1,0.25}，λ={λ1,λ2}={1,3}，β={β1,β2}={0.25,1}，初始值为X0=4，y˙0=λ1。值得注意的是，强度过程的初始状态将对转移密度的演变产生影响。例如，如果强度过程的初始值被设定为λ2，那么在转移期内发生跳跃的概率将大大增加。因此，考虑到我们从高频跳跃制度开始，转移密度与低频区制有明显的不同。此外，如果强度参数y˙t的初始状态是未知的，我们可以通过考虑转移密度的混合分布来解释这个延迟。在这种情况下，我们将通过计算两种强度区制的转移密度来构建转移密度，然后根据CTMC的平稳分布（或适当的初始分布）对每个分布加权。

由于强度过程现在取决于一个外部过程，我们需要稍微修改Kolmogorov公式。

由此可以看出，为了评估这个过程的矩方程，我们需要强度过程的期望值随时间的演变。在本例中，强度参数的期望值的分析表达式为：。 

因此，在R:

JGQD.remove()

1.  
    
2.  
   G0=function(t){2*5+2*sin(1*pi*t)}
3.  
   # 强度假设为强度过程的期望轨迹。
4.  
   plot(Lam0(t)~t,type='l')
5.  
    

出于比较的目的，我们可以比较过程中各时刻的模拟轨迹，并将转移密度的模拟轨迹与近似值的模拟轨迹进行比较。这可以通过模拟算法来实现。 

1.  
   #' 现在模拟跳跃扩散
2.  
    
3.  
   mu <- function(x,t){G0(t)+G1(t)*x} # 漂移
4.  
   sigma <- function(x,t){sqrt(Q1(t)*x)} # 扩散性
5.  
   j <- function(x,z){z} # 跳跃
6.  
    
7.  
   res2 <- simulate()

随后，我们可以看看这样一个过程的典型轨迹是什么样子的，并比较模拟过程的轨迹和矩截断法的轨迹。

1.  
   par(mfrow=c(3,1))
2.  
   plot(res2$xtrak~res2$tt,type='l',col='blue',main='Trajectory',xlab = 'time',ylab ='X_t')
3.  
   plot(res2$jtrak~res2$tt,type='h',col='black',ylim=c(0,3),lwd=2,main='Jumps',xlab ='time',ylab ='Z_t')
4.  
   plot(res2$etrak~res2$tt,type='s',col='black',ylim=c(0,6),lwd=1,main='Intensities',xlab ='time',ylab ='Z_t')
5.  
   abline(h =c(l1,l2),lty='dotted',col='lightgrey')

1.  
   par(mfrow=c(2,2))
2.  
   for(i in 1:4)
3.  
   {
4.  
   plot(res2$MM[i,]~res2$tt,type='l',main='Moment trajectory',xlab='Time (t)',
5.  
   ylab=paste0('m_',i,'(t)'))
6.  
   lines(res$moments[i,]~res$time,lty='dashed',col='blue',lwd=2)
7.  
   }

最后我们可以将模拟的转移密度与模拟过程的密度进行比较。 

 persp(x=Xt,y=time,zdensit)

1.  
   par(mfrow=c(2,2))
2.  
   for(i in 1:4)plot(hists[[i]]

现在，我们可以在一个短的转移期做一个类似的实验。为了便于比较，我们将把扩散参数缩小一些，以夸大跳跃机制的影响。 

1.  
   # 缩减扩散系数。
2.  
   Q1=function(t){0.2}
3.  
    
4.  
   # 近似，但使用因子化。
5.  
   TT <- 0.5
6.  
   density(4,factorize = TRUE)

1.  
   # 重新模拟并记录新点的直方图。
2.  
   res2 <- simulate())
3.  
    
4.  
   persp(xXt,time,zdensity)

1.  
   par(mfrow=c(2,2))
2.  
   for(i in 1:4)
3.  
   {
4.  
   plot(res2$hists[[i]]$density~c(res2$hists[[i]]$mids-diff(res2$hists[[i]]$mids)[1] / 2),
5.  
   type = 's',lty = 'solid', lwd = 1, xlab = 'time',
6.  
   ylab = 'Density', main = paste('Density at time t =',i*0.05))
7.  
   lines(res$density[,i*10]~res$Xt,col='darkblue')
8.  
   }

---

最受欢迎的见解

1.用R语言模拟混合制排队随机服务排队系统

2.R语言中使用排队论预测等待时间

3.R语言中实现马尔可夫链蒙特卡罗MCMC模型

4.R语言中的马尔科夫机制转换(Markov regime switching)模型

5.matlab贝叶斯隐马尔可夫hmm模型

6.用R语言模拟混合制排队随机服务排队系统

7.Python基于粒子群优化的投资组合优化

8.R语言马尔可夫转换模型研究交通伤亡人数事故预测

9.用机器学习识别不断变化的股市状况——隐马尔可夫模型的应用

---