扩散过程的无偏蒙特卡罗模拟

（22）3.2校正项利用方程（10）和（11），我们的选择（21）对应于一个具有漂移和波动性的过程tk+t、 Stk+f（k）(TW）= u（tk，Stk）+σ（tk，Stk）Su（tk，Stk）W^σtk+t、 Stk+f（k）(TW）= σ（tk，Stk）1 + Sσ（tk，Stk）W+Su（tk，Stk）T.由于贴现率r（t）是确定的，我们也取^r（t，St）=r（t）。将这三个函数用于k- 1与tk=tk-1=tk- tk-1和W=工作-1=周- 工作-1根据我们有1个+Htkλ=1+ukλStk+Ckλ斯特奎斯uk=u（tk，Stk）- ^utk-1+ tk-1、Stk-1+f(tk-1.工作-1)Ck=σ（tk，Stk）- ^σtk-1+ tk-1、Stk-1+f(tk-1.工作-1).通过函数f（k）的显式选择，我们得到uk=u（tk，Stk）-u（tk）-1、Stk-1） +σ（tk）-1，Stk）Su（tk-1、Stk-1)工作-1.Ck=σ（tk，Stk）- σ（tk）-1、Stk-1)1 + Sσ（tk）-1、Stk-1)工作-1(23)+Su（tk-1、Stk-1)tk-1..3.2.1中间时间如果tkis不是最后一个泊松时间，即k<p，香港演艺学院onuk（tk，Stk）=ZdStk+1块，tk+1（Stk，Stk+1）1 +Hk+1λutk+1（Stk+1）。（24）更一般地，我们将考虑二阶微分算子ak=1+ASk对该表达式的作用Stk+ASSkStk.为了计算表达式（24）相对于stkWe的导数，考虑变量从（t，S）到（t，W）的变化，如等式（22）所定义：S=S*+ f（k）（t）- tk，W- W*) .和S*= 斯特坎德W*= Wtk。我们介绍*而W*明确表示一旦选择函数f，S*而W*是恒定的，不应被区分。变量的这种变化导致了第一个导数W=(WS）S=Wf（k）（t）- tk，W- Wtk）S=^σS.（25）第二次区分我们的WW=^σS+(W^σ）S=^σS+(Wf（k））s

（26）对于函数f，我们的特殊选择是W=^σS+σ（tk，Stk）Sσ（tk，Stk）把这些方程倒过来S=^σWS=^σW-W^∑∑因此我们可以用变量WAk=1+ASk^σ（t+k，Wtk）来写W+ASSk^σ（t+k，Wtk）W-W^σ（t+k，Wtk）^σ（t+k，Wtk）W.再次使用^σ（t+k，Wtk）=σ（tk，Stk）W^σ（t+k，Wtk）=σ（tk，Stk）Sσ（tk，Stk）。这个becomesAk=1+询问σ（tk，Stk）-混蛋Sσ（tk，Stk）σ（tk，Stk）W+ASSkσ（tk，Stk）W.（27）变量从（t，S）变为（t，W）后，方程（24）becomesuk（tk，Wtk）=ZdWtk+1bU（W）tk，tk+1（Wtk，Wtk+1）1 +Hk+1λutk+1（Wtk+1）。（28）由于W是布朗运动，就新变量而言，演化算子bU（W）是贴现因子与高斯核的乘积：bU（W）tk，tk+1（Wtk，Wtk+1）=e-Rtk+1tkr（s）ds~n（tk+1- tk，Wtk+1- Wtk）带(TW）=√2πte-W2t、 （29）我们必须计算微分算子对ukgivenin公式（28）的影响，我们扩展asuk（tk，Wtk）=ZdWtk+1bU（W）tk，tk+1（Wtk，Wtk+1）utk+1（Wtk+1）+ZdWtk+1bU（W）tk，tk+1（Wtk，Wtk+1）Hk+1λutk+1（Wtk+1）。（30）当作用于本公式中的第二项时Hk+1λ，Wtkare上的导数替换为Malliavin权重：WtkbU（W）tk，tk+1（Wtk，Wtk+1）=WW(tk，Wk）bU（W）tk，tk+1（Wtk，Wtk+1）WtkbU（W）tk，tk+1（Wtk，Wtk+1）=WW W(tk，Wk）bU（W）tk，tk+1（Wtk，Wtk+1），其中高斯核（29）的显式形式给出了ww(tk，（工作）=工作将军澳(tk，（工作）=工作- tk蒂克。对于表达式（30）的第一项相对于Wtk的导数，我们使用高斯核的对称性，因此对于Wtk和Wtk+1，我们使用ofbU（W）的对称性来写入WtkbU（W）tk，tk+1（Wtk，Wtk+1）=-Wtk+1bU（W）tk，tk+1（Wtk，Wtk+1）。然后我们在（30）的第一项中进行部分积分，得到WtkZdWtk+1bU（W）tk，tk+1（Wtk，Wtk+1）utk+1（Wtk+1）=ZdWtk+1bU（W）tk，tk+1（Wtk，Wtk+1）Wtk+1utk+1（Wtk+1）和类似的二阶导数。

然后我们从W返回到S变量时间tk+1，使用（25）和（26），在我们的例子中是Wtk+1=σ（tk，Stk）1 + Sσ（tk，Stk）Wk+Su（tk，Stk）tkStk+1Wtk+1=σ（tk，Stk）1 + Sσ（tk，Stk）Wk+Su（tk，Stk）tkStk+1+σ（tk，Stk）Sσ（tk，Stk）Stk+1。综合所有术语，我们最终得出AKUK=ZdStk+11 +询问σ（tk，Stk）-混蛋Sσ（tk，Stk）σ（tk，Stk）工作tk+ASSkσ（tk，Stk）工作- tktk但是，tk+1（Stk，Stk+1）Hk+1λutk+1（Stk+1）+ZdStk+1bUtk，tk+1（Stk，Stk+1）1 +问-混蛋Sσ（tk，Stk）σ（tk，Stk）1 + Sσ（tk，Stk）Wk+Su（tk，Stk）tkStk+1+ASSk1 + Sσ（tk，Stk）Wk+Su（tk，Stk）tkStk+1+ASSkSσ（tk，Stk）σ（tk，Stk）Stk+1utk+1（Stk+1）。我们可以用Ak+1=1+ASk+1+1重写asAkuk=ZdStk+1，tk+1（Stk，Stk+1）Ak+1utk+1（Stk+1）Stk+1+ASSk+1Stk+1ASk+1=1+bk问-bkSσ（tk，Stk）σ（tk，Stk）ASSk+dk(tk，（工作）uk+1λ（31）ASSk+1=1+bkASSk+dk(tk，（工作）Ck+1λ和bk=Sσ（tk，Stk）Wk+Su（tk，Stk）tkdk(TW）=1+询问σ（tk，Stk）-混蛋Sσ（tk，Stk）σ（tk，Stk）Wt+ASSkσ（tk，Stk）W- Tt、 （32）使用公式（31），我们可以递归地累积修正项。我们从t=t时的A=1开始。然后从每个日期开始，通过方程（22）模拟Stk+1，并计算上文定义的Ak+1，直到最后一个泊松时间tp。3.2.2根据最终泊松时间tp，p应采取行动onutp（Stp）=ZdSTbUtp，T（Stp，ST）h（ST），其中h是支付函数。一种幼稚的做法是：。根据方程式（22）中给出的Stp模拟Stp:ST=Stp+f（p）(tp，可湿性粉剂）tp=T- tpandWp=WT- 高斯方差变量tp。然后计算支付（ST），并使用Malliavinweights计算导数。这意味着将贴现支付乘以dp(tp，Wp）给定不等式（32）：PT=dp(tp，Wp）h装货单.然而，这个术语的行为就像O总磷很小t而不是我们想要的O（1）。取而代之的是，我们将在最后一个时间步骤中使用对偶采样，如HenryLabord`ere等人（2015年）所述。

更准确地说，我们计算（+）T=Stp+f（p）(tp，Wp）S（0）T=^Etp装货单= Stp+u（tp，Stp）tp（33）S(-)T=Stp+f（p）(tp，-Wp），然后再乘以dp=dp(tp，Wp）hS（+）T+dp(tp，-Wp）hS(-)T+1.-dp(tp，（可湿性粉剂）-dp(tp，-（可湿性粉剂）HS（0）T哪一个简单的例子toPT=dp(tp，Wp）hS（+）T+dp(tp，-Wp）hS(-)T-助理σ（tp、Stp）可湿性粉剂- 总磷tphS（0）T. （34）h中的最后一个术语S（0）T预期值为0，第一项和第二项对期权价格的贡献相反。如果Payoff函数h是光滑的，并且具有泰勒展开式，则可以检查O（1）阶的PTI。对于罢工附近的看涨期权或看跌期权支付，情况并非如此。在这种情况下，我们有PT=O√总磷. 然而，两人之间发生罢工的可能性(-)最后一个时间步测量asO（ptp）。因此，差异仍然有限。我们最终会在所有时间步长上打折。在第一种情况下，当我们采用确定性贴现率时，这一因素被表示为乘以e-RTtr（t）dt。3.3蒙特卡罗方案总结我们现在用确定性贴现因子总结一维情况下的整个蒙特卡罗方案。我们考虑一种到期日为T且付息（ST）的欧式期权。对于每一条路径，我们都要执行以下操作：1。从S=Sat时间t=t开始。定义=0ASS=0.2。在每个日期tk（a）绘制下一个泊松时间tk+1=tk+强度为λ的TKW。例如，画一个0到1之间的随机均匀变量q，然后取tk=-log（q）λ。如果tk+1>T，则设置p=k并转至步骤3。（b） 画一个高斯变量有差异的工作蒂克。

通过方程（21）和（22）得到Stk+1:Stk+1=Stk+f（k）(tk，工作）。带f（k）(TW）=u（tk，Stk）t+σ（tk，Stk）W+σ（tk，Stk）Sσ（tk，Stk）(W- t） +σ（tk，Stk）Su（tk，Stk）TW（c） 计算uK根据等式（23）：uk=u（tk，Stk）-u（tk）-1、Stk-1） +σ（tk）-1，Stk）Su（tk-1、Stk-1)工作-1.Ck=σ（tk，Stk）- σ（tk）-1、Stk-1)1 + Sσ（tk）-1、Stk-1)工作-1+Su（tk-1、Stk-1)tk-1..（d） 从等式（31）中的ASk和ASSkas计算ASk+1和ASSk+1:ASk+1=1+bk问- bkSσ（tk，Stk）σ（tk，Stk）ASSk+dk(tk，（工作）uk+1λASSk+1=1+bkASSk+dk(tk，（工作）Ck+1λ，其中bk=Sσ（tk，Stk）Wk+Su（tk，Stk）tkdk(TW）=1+询问σ（tk，Stk）-混蛋Sσ（tk，Stk）σ（tk，Stk）Wt+ASSkσ（tk，Stk）W- Tt、 三,。准时tp（a）绘制一个高斯变量方差分析tp=T- tp。计算（+）T，S（0）和S(-)式（33）中的Tas:S（+）T=Stp+f（p）(tp，Wp）S（0）T=Stp+u（tp，Stp）tpS(-)T=Stp+f（p）(tp，-可湿性粉剂）。（b） 获取支付（ST）的未贴现路径贡献（34）：PT=dp(tp，Wp）hS（+）T+dp(tp，-Wp）hS(-)T-助理σ（tp、Stp）可湿性粉剂- 总磷tphS（0）T.4.乘以贴现系数：e-RTtr（t）dtPT。我们最终对所有路径进行平均，以获得期权价格的无偏蒙特卡罗估计。4数值结果4。1收敛为了在数值上检验我们的蒙特卡罗格式，我们采用了一个模型，我们可以用一个闭合公式作为参考值。因此，我们选择BlackScholes模型ds=uSdt+σsdwt，其对应于局部漂移波动率u（t，S）=uSσ（t，S）=σS。该模型已经存在一个无偏蒙特卡罗方案，使用对数作为变量。为了测试的目的，我们忽略了这一点，我们天真地应用了你的方案，并将其与欧拉方案进行了比较。我们考虑现货S=100且波动率σ=50%的标的证券。我们提高利率，漂移率r=u=5%。我们对到期日为T=1年的看跌期权进行定价，行使期为K=80。

布莱克-斯科尔斯公式给出的期权价格为7.8909。图1显示了λ=3的无偏方案的模拟值与路径数的函数关系。路径数98。587.5200000150000100000理论规则Runbiased500000期权价格图1：无偏方案和Euler方案的模拟期权价格相对于路径数的收敛性。为了进行比较，它还显示了Euler模式的相同收敛图，其平均时间步数n=4。我们看到后一种模式趋向于一个有偏差的值。4.2与Euler和Milstein方案相比，使用Euler方案时，可以使用更多的时间步长来减少偏差。这反过来又增加了计算时间。另一方面，使用无偏方案的代价是增加λ小值的方差。增加λ值会使方差变小，但也会增加计算时间，因为平均时间步长更高。为了评估可实现的性能增益，我们考虑了第4.1节中描述的80%看跌期权，并用两种方案对其定价。由于路径生成与Milstein方案有相似之处，我们还将其与Milstein方案进行了比较。这样可以更好地隔离无偏方案中纠正术语的影响。对于无偏格式，我们使用λ值从.01增加到29。对于Eulerand-Milstein方案，我们将时间步长从1增加到300。在所有情况下，我们绘制N=100万条路径，并用数值计算估计的期权价格、经验标准偏差和计算时间。结果如图2所示。

根据对数标度的计算时间，以99%置信区间绘制估计价格。775885995200 2000 20000Option price Computing time（对数标度）TheoreticalerMilsteinnbiased图2：使用无偏方案获得的欧洲看跌期权的模拟价格，与具有相同数量Monte Carlopath的Euler和Milstein方案进行比较。价格以对数比例以99%的置信区间与计算时间绘制。在非常低的泊松强度（λ=0.01）下，无偏格式的蒙特卡罗方差高于Euler和Milstein格式。然而，当λ增大时，它迅速减小。对于λ≥ 0.3甚至比theEuler和Milstein方案的噪声更低。超过λ=1时，蒙特卡罗噪声几乎是平稳的：与支付的基本方差相比，校正项增加了可忽略的方差。对于λ的所有值，我们检查最终估计值是否与理论值一致，直至蒙特卡罗统计误差。就计算时间而言，λ=1似乎是最佳值，对应于1年期间的平均时间步数n=2。相反，当时间步长较小时，欧拉格式表现出较大的偏差。这种偏差随时间线性减小T=T/n。进行加权最小二乘回归，我们估计偏差表现为渐近2。99n。为了使偏差等于欧拉蒙特卡罗标准偏差0。因此，我们需要n~ 230个时间步。这相当于计算时间比λ=1的无偏格式长43倍。Milstein方案具有较小的偏差，渐近-在本例中为0.648。因此，我们需要n=50个时间步才能获得与蒙特卡罗噪声相同大小的偏差。

这意味着计算速度比无偏方案慢10倍。5多维过程我们现在考虑一般的多维情况，其中抛物线PDE（1）中的所有参数都可以位于t和X上：tut（Xt）+μα（t，Xt）αut（Xt）+Cαβ（t，Xt）αβut（Xt）=r（t，Xt）ut（Xt）。（35）根据Feynman-Kac定理，它对应于具有随机贴现率r（t，Xt）的多维过程dxαt=μα（t，Xt）dt+σαa（t，Xt）dwat。σ（t，Xt）是一个波动矩阵，使得cαβ（t，Xt）=σαa（t，Xt）σβa（t，Xt），可以通过Cholesky分解得到。它们是独立的标准布朗过程。5.1蒙特卡罗路径我们绘制强度为λ的泊松乘以tkw。在这段时间内，我们模拟了ad维布朗过程Wt。这对应于抛物线偏微分方程tvt（Wt）+δabA.bvt（Wt）=0。（36）在两个泊松时间tk和tk+之间，我们将考虑两个变量的变化：o空间变量。我们使用函数fα（k）从（t，W）到（t，X）(TW）：Xαt=Xαtk+fα（k）（t）- tk，Wt- Wtk）。如果我们将由空间变量和价格构成的总空间视为一个纤维束，这对应于基础空间和纤维上变量的变化纳姆·埃莱尔。

我们使用函数g（k）来改变数值(TW）ut=vteg（k）（t）- tk，Wt- Wtk）。将股份转让给衍生工具的金额变更如下：tvt=e-G啧啧- (tg）e-肠胃avt=e-Gaut- (ag）e-肠胃A.bvt=e-GA.但是- (A.bg）e-内脏+(（德国）(bg）e-肠胃-(ag）e-G但是- (bg）e-G奥特。使用符号a=佤α= Xα我们定义αa=afα和cαab=A.bfα。我们还将使用逆矩阵eaα，定义性质为αaeaβ=Δαβeaαeαb=δab。空间变量的变化导致导数发生以下变换：a=eαaαA.b=eαaeβbαβ+cαabαt | W=t | X+tfαα.根据新变量t、X和ut，PDE（36）因此变为图坦卡蒙αut+^Cαβαβut=^rut（37）带^μα=tfα- (ag）eαa+cαaa^cαβ=eαaeβa（38）^r=甘油三酯+A.银-(（德国）(ag）。我们想找到函数fα（k）(TW）和g（k）(TW）使- ^u，C-^C和r- ^r表现为O(t） 很小t、 为了简化符号，我们在不需要（k）索引时删除它们。我们假设我们可以把fα（k）和g（k）展开为幂级数t和W:fα（k）(TW）=fα+（fα）aWa+fαt+2！（fα）ab佤Wb+（fα）aTWa+3！（fα）abc佤WbWc+·g（k）(TW）=g+（g）aWa+gt+2！（g） ab佤Wb+（g）aTWa+3！（g） abc佤Wb在这个表达式中，（fαij）abc··或（gij）abc··是在空间指数a、b、c中对称的张量。t=tkis fα（0，0）=0时的连续性约束。它转换为tofα=0。此外，我们选择数值，使得utand和Vt在周期开始时重合，t=tk。数学上这是g（0，W）=0，这表示所有j（g0j）=0。考虑到这些约束条件，方程（38）具有泰勒展开式^μμα=fα+（fα）aWa+（fα）bb+（fα）bbaWa+O(t） （39）^Cαβ=（fα）b（fβ）b+（fα）b（fβ）baWa+（fβ）b（fα）baWa+O(t） （40）^r=g+（g）aWa+O(t） 。

（41）另一方面，PDE（35）参数的泰勒展开式给出了μα=μα（tk，Xtk）+eγa（tk，Wtk）γ-α（tk，Xtk）Wa+O(t） （42）Cαβ=Cαβ（tk，Xtk）+eγa（tk，Wtk）γCαβ（tk，Xtk）Wa+O(t） （43）r=r（tk，Xtk）+eγa（tk，Wtk）γr（tk，Xtk）Wa+O(t） 。（44）从定义来看，eγaiseγa=afα=（fγ）a+（fγ）abWb+O(t） 。在期初t=0和W=0，这是iseγa（tk，Wtk）=（fγ）a。我们想要得到u- ^u，C-^C和r- ^r消失到O(t） 条款。让我们从方程（40）和（43）中的常数项开始。从分解αβ（tk，Xtk）=∑αb（tk，Xtk）σβb（tk，Xtk）得到一个解（fα）b=∑αb（tk，Xtk）。（45）这也给出了αa（tk，Wtk）=σαa（tk，Xtk）。它的逆iseaα（tk，Wtk）=σ（tk，Wtk）-1.aα。将方程（40）和（43）的一阶项相等，我们有eαb（tk，Wtk）（fβ）ba+eβb（tk，Wtk）（fα）ba=eγa（tk，Wtk）γCαβ（tk，Xtk），其中Cαβ在指数α和β中是对称的，而（fβ）在指数a和b中是对称的。我们将这个方程乘以ecα（tk，Wtk）和edβ（tk，Wtk），也使用eaαeαb=δab:fdca+fcda=Cacd（46），其中我们引入了符号sfabc=eaα（tk，Wtk）（fα）bcCabc=eαa（tk，Wtk）ebβ（tk，Wtk）ecγ（tk，Wtk）αCβγ（tk，Xtk）。最后两个指数中对称的fabcand-Cabcare张量。我们为指数的三个循环排列写出方程（46），getfbca+fcba=Cabc（47）fcab+facb=Cbca（48）fabc+fbac=Ccab。