碼上歡樂
相關內容    簡體    繁體

QuantLib 金融計算——高級話題之模擬跳擴散過程

本文轉載自   查看原文   2019-03-31 11:46   524    QuantLib 金融計算

目錄

如果未做特別說明，文中的程序都是 C++11 代碼。

QuantLib 金融計算——高級話題之模擬跳擴散過程

跳擴散過程

1976 年，Merton 最早在衍生品定價中引入並分析了跳擴散過程，正因為如此 QuantLib 中和跳擴散相關的隨機過程類稱之為 Merton76Process，一個一般的跳擴散過程可以由下面的 SDE 描述，

\[\frac{dS(t)}{S(t-)} = \mu dt + \sigma dW(t) + dJ(t)\\ J(t) = \sum_{j=1}^{N(t)}(Y_j-1) \]

其中 \(Y_j\) 是隨機變量，\(N(t)\) 是計數過程。\(dJ(t)\) 表示 \(J(t)\) 在 \(t\) 時刻發生的跳躍幅度，如果發生一次跳，則幅度為 \(Y_j-1\)；如果沒有發生跳，則幅度為 0。

在應用於衍生品定價時，需要對上述 SDE 中的項做出一些特殊約定，常見的約定有：

  1. \(N(t)\) 是參數等於 \(\lambda\) 的 Poisson 過程；
  2. \(\log Y_j\) 服從正態分布 \(N(\mu_{jump}, \sigma_{jump}^{2})\)；
  3. \(\mu = r - \lambda m\)，其中 \(m = E[Y_j] - 1\)，\(r\) 代表無風險利率

模擬算法

令 \(X(t) = \log S(t)\)，那么

\[\begin{aligned} X(t_{i+1}) = & X(t_i) + \left(\mu - \frac12 \sigma^2 \right)(t_{i+1} - t_i)\\ & +\sigma[W(t_{i+1}) - W(t_i)] + \sum_{j = N(t_i)+1}^{N(t_{i+1})}\log Y_j \end{aligned} \]

在 \(X(t_i)\) 的基礎上模擬 \(X(t_{i+1})\) 的算法如下：

  1. 生成 \(Z \sim N(0,1)\)；
  2. 生成 \(N \sim \text{Poisson}(\lambda(t_{i+1}-t_i))\)，若 \(N=0\)，則令 \(M=0\)，並轉到第 4 步；
  3. 生成 \(\log Y_1,\dots,\log Y_N\)，令 \(M = \log Y_1+\dots+\log Y_N\)；
  4. 令 \(X(t_{i+1}) = X(t_i) + \left(\mu - \frac12 \sigma^2 \right)(t_{i+1} - t_i) +\sigma \sqrt{t_{i+1} - t_i} Z + M\)

那么

\[\begin{aligned} S_{t_{i+1}} =& S_{t_i} e^{(r-\lambda m -\frac{1}{2}\sigma^{2})\Delta t+ \sigma Z \sqrt{\Delta t} + M}\\ =& S_{t_i} e^{(r-\frac{1}{2}\sigma^{2})\Delta t+ \sigma Z \sqrt{\Delta t}}e^{(-\lambda m)\Delta t+M} \end{aligned} \]

其中，\(\Delta t = t_{i+1} - t_i\)，而 \(e^{(-\lambda m)\Delta t+M}\) 是跳擴散相對於一般 Black-Scholes-Merton 過程的修正項。

面臨的問題

目前 QuantLib 中提供的 Merton76Process 類只能配合“解析定價引擎”使用，本身不具備模擬隨機過程路徑的功能。究其原因，問題出在 QuantLib 的編碼約定和 StochasticProcess1D 提供的接口兩方面：

  1. QuantLib 中約定 StochasticProcess 派生出的子類僅能描述 SDE 的結構信息，也就是 SDE 的參數、漂移和擴散項的函數形式，子類不攜帶有關隨機數生成的信息，所有隨機數生成的相關信息均由 Monte Carlo 框架中其他組件控制；
  2. 生成隨機過程路徑的核心函數是 evolve 方法，StochasticProcess1D 提供的接口是 evolve(Time t0, Real x0, Time dt, Real dw)。如果 Merton76Process 按約定實現 evolve 方法的話，形式必須是 evolve(Time t0, Real x0, Time dt, const Array &dw)，因為模擬跳需要額外的隨機性，所以 dw 必須是一個 Array。很明顯，不匹配。

“臟”的方法

在不改變當前接口的前提下，若要實現模擬跳擴散過程，需要用比較“臟”一點兒的手段，即打破約定，讓隨機過程類攜帶一個隨機數發生器，為模擬跳提供額外的隨機性。

具體來說，需要聲明一個 Merton76Process 的派生類，該類攜帶一個高斯隨機數發生器。因為從數學上來講跳擴散過程推廣自一般 Black-Scholes-Merton 過程，添加了一個修正項，所以遵循“適配器模式”（或“裝飾器模式”）的思想，絕大部分計算可以委托給一個 BlackScholesMertonProcess 對象，僅需要對 drift 和 evolve 方法作必要的修改。

“干凈”的方法

當然，“干凈”的方法要改變當前接口：

  1. 聲明一個和 StochasticProcess1D 平行的新類 StochasticProcess1DJump，二者唯一的區別是 evolve 方法，在 StochasticProcess1DJump 中形式是 evolve(Time t0, Real x0, Time dt, const Array &dw)；
  2. 將 Merton76Process 改成繼承自 StochasticProcess1DJump。

實現

下面的代碼實現了前面提到的“臟”的方法，因為隨機數發生器的種類有很多，且沒有基類提供統一的接口，所以使用了模板技術讓類可以接受不同類型的隨機數發生器。同時，許多計算被委托給了一個 BlackScholesMertonProcess 對象。

#ifndef MERTON76JUMPDIFFUSIONPROCESS_HPP
#define MERTON76JUMPDIFFUSIONPROCESS_HPP

#include <ql/math/distributions/normaldistribution.hpp>
#include <ql/math/distributions/poissondistribution.hpp>
#include <ql/math/randomnumbers/boxmullergaussianrng.hpp>
#include <ql/math/randomnumbers/mt19937uniformrng.hpp>
#include <ql/processes/blackscholesprocess.hpp>
#include <ql/processes/merton76process.hpp>

namespace QuantLib
{
template<typename GAUSS_RNG>
class Merton76JumpDiffusionProcess : public Merton76Process
{
  public:
    Merton76JumpDiffusionProcess(const Handle<Quote>& stateVariable,
                                 const Handle<YieldTermStructure>& dividendTS,
                                 const Handle<YieldTermStructure>& riskFreeTS,
                                 const Handle<BlackVolTermStructure>& blackVolTS,
                                 const Handle<Quote>& jumpInt,
                                 const Handle<Quote>& logJMean,
                                 const Handle<Quote>& logJVol,
                                 const GAUSS_RNG& gauss_rng,
                                 const ext::shared_ptr<discretization>& disc =
                                     ext::shared_ptr<discretization>(new EulerDiscretization))
      : Merton76Process(
            stateVariable,
            dividendTS,
            riskFreeTS,
            blackVolTS,
            jumpInt,
            logJMean,
            logJVol,
            disc)
      , blackProcess_(
            new BlackScholesMertonProcess(
                stateVariable,
                dividendTS,
                riskFreeTS,
                blackVolTS,
                disc))
      , gauss_rng_(gauss_rng)
    {
    }
    virtual ~Merton76JumpDiffusionProcess() {}

    Real x0() const
    {
        return blackProcess_->x0();
    }
    Time time(const Date& d) const
    {
        return blackProcess_->time(d);
    }
    Real diffusion(Time t,
                   Real x) const
    {
        return blackProcess_->diffusion(t, x);
    }
    Real apply(Real x0,
               Real dx) const
    {
        return blackProcess_->apply(x0, dx);
    }
    Size factors() const
    {
        return 1;
    }
    Real drift(Time t,
               Real x) const
    {
        Real lambda_ = Merton76Process::jumpIntensity()->value();
        Real delta_ = Merton76Process::logJumpVolatility()->value();
        Real nu_ = Merton76Process::logMeanJump()->value();
        Real m_ = std::exp(nu_ + 0.5 * delta_ * delta_) - 1;

        return blackProcess_->drift(t, x) - lambda_ * m_;
    }
    Real evolve(Time t0,
                Real x0,
                Time dt,
                Real dw) const;

  private:
    const CumulativeNormalDistribution cumNormalDist_;
    ext::shared_ptr<GeneralizedBlackScholesProcess> blackProcess_;
    GAUSS_RNG gauss_rng_;
};

template<typename GAUSS_RNG>
Real Merton76JumpDiffusionProcess<GAUSS_RNG>::evolve(Time t0,
                                                     Real x0,
                                                     Time dt,
                                                     Real dw) const
{
    Real lambda_ = Merton76Process::jumpIntensity()->value();
    Real delta_ = Merton76Process::logJumpVolatility()->value();
    Real nu_ = Merton76Process::logMeanJump()->value();
    Real m_ = std::exp(nu_ + 0.5 * delta_ * delta_) - 1;

    Real p = cumNormalDist_(gauss_rng_.next().value);
    if (p < 0.0)
        p = 0.0;
    else if (p >= 1.0)
        p = 1.0 - QL_EPSILON;

    Real j = gauss_rng_.next().value;
    const Real n = InverseCumulativePoisson(lambda_ * dt)(p);
    Real retVal = blackProcess_->evolve(
        t0, x0, dt, dw);
    retVal *=
        std::exp(-lambda_ * m_ * dt + nu_ * n + delta_ * std::sqrt(n) * j);

    return retVal;
}
}
#endif // MERTON76JUMPDIFFUSIONPROCESS_HPP

示例

下面模擬兩條曲線

#include <iostream>

#include <ql/math/randomnumbers/boxmullergaussianrng.hpp>
#include <ql/math/randomnumbers/mt19937uniformrng.hpp>
#include <ql/processes/blackscholesprocess.hpp>
#include <ql/quotes/simplequote.hpp>
#include <ql/termstructures/volatility/equityfx/blackconstantvol.hpp>
#include <ql/termstructures/yield/flatforward.hpp>
#include <ql/time/calendars/target.hpp>
#include <ql/time/date.hpp>
#include <ql/time/daycounters/actualactual.hpp>

#include "Merton76JumpDiffusionProcess.hpp"

int main() {
    using namespace std;
    using namespace QuantLib;

    Date refDate = Date(27, Mar, 2019);
    Rate riskFreeRate = 0.03;
    Rate dividendRate = 0.01;
    Real spot = 52.0;
    Rate vol = 0.2;
    Calendar cal = TARGET();
    DayCounter dc = ActualActual();

    ext::shared_ptr<YieldTermStructure> rdStruct(
        new FlatForward(refDate, riskFreeRate, dc));
    ext::shared_ptr<YieldTermStructure> rqStruct(
        new FlatForward(refDate, dividendRate, dc));
    Handle<YieldTermStructure> rdHandle(rdStruct);
    Handle<YieldTermStructure> rqHandle(rqStruct);

    ext::shared_ptr<SimpleQuote> spotQuote(new SimpleQuote(spot));
    Handle<Quote> spotHandle(spotQuote);

    ext::shared_ptr<BlackVolTermStructure> volQuote(
        new BlackConstantVol(refDate, cal, vol, dc));
    Handle<BlackVolTermStructure> volHandle(volQuote);

    // Specify the jump intensity, jump mean and
    // jump volatility objects

    Real jumpIntensity = 0.2;    // lambda
    Real jumpVolatility = 0.3;
    Real jumpMean = 0.0;

    ext::shared_ptr<SimpleQuote> jumpInt(new SimpleQuote(jumpIntensity));
    ext::shared_ptr<SimpleQuote> jumpVol(new SimpleQuote(jumpVolatility));
    ext::shared_ptr<SimpleQuote> jumpMn(new SimpleQuote(jumpMean));

    Handle<Quote> jumpI(jumpInt), jumpV(jumpVol), jumpM(jumpMn);

    ext::shared_ptr<BlackScholesMertonProcess> bsmProcess(
        new BlackScholesMertonProcess(
            spotHandle, rqHandle, rdHandle, volHandle));

    unsigned long seed = 12324u;
    MersenneTwisterUniformRng unifMt(seed);
    MersenneTwisterUniformRng unifMtJ(25u);

    typedef BoxMullerGaussianRng<MersenneTwisterUniformRng> GAUSS;

    GAUSS bmGauss(unifMt);
    GAUSS jGauss(unifMtJ);

    ext::shared_ptr<Merton76JumpDiffusionProcess<GAUSS>> mtProcess(
        new Merton76JumpDiffusionProcess<GAUSS>(
            spotHandle, rqHandle, rdHandle, volHandle,
            jumpI, jumpM, jumpV, jGauss));

    Time dt = 0.004, t = 0.0;
    Real x = spotQuote->value();
    Real y = spotQuote->value();
    Real dw;
    Size numVals = 250;

    std::cout << "Time, Jump, NoJump" << std::endl;
    std::cout << t << ", " << x << ", " << y << std::endl;

    for (Size j = 1; j <= numVals; ++j) {
        dw = bmGauss.next().value;
        x = mtProcess->evolve(t, x, dt, dw);
        y = bsmProcess->evolve(t, y, dt, dw);
        std::cout << t + dt << ", " << x << ", " << y << std::endl;
        t += dt;
    }

    return EXIT_SUCCESS;
}

參考文獻

  1. 《金融工程中的蒙特卡羅方法》


免責聲明！

本站轉載的文章為個人學習借鑒使用，本站對版權不負任何法律責任。如果侵犯了您的隱私權益，請聯系本站郵箱yoyou2525@163.com刪除。



猜您在找 levy過程、擴散過程、隨機過程帶跳 QuantLib 金融計算——隨機過程之概述 QuantLib 金融計算——QuantLib 入門 QuantLib 金融計算 QuantLib 金融計算——隨機過程之一般 Black Scholes 過程 QuantLib 金融計算——隨機過程之 Heston 過程 QuantLib 金融計算——原理之 Bootstrap QuantLib 金融計算——數學工具之插值 QuantLib 金融計算——基本組件之 Calendar 類 QuantLib 金融計算——基本組件之 Schedule 類
 
粵ICP備18138465號   © 2018-2025 CODEPRJ.COM