２. "Implementing QuantLib"の和訳
Chapter-VI The Monte Carlo Framework
6.1 Pathの生成 （つづき）
6.1.2 Stochastic Processes : 確率過程
RNGを使う目的は、当然ながら、生成された乱数を使って確率過程を生成していく事です。StochasticProcessクラスのインターフェースは （下記 Listing 6.4に示します）、その利用方法を念頭に置いて設計されています。
Listing 6.4: Partial implementation of the StochasticProcess class.
class StochasticProcess : public Observer,
public Observable {
public:
class discretization {
public:
virtual ~discretization() {}
virtual Disposable<Array> drift(const StochasticProcess&,
Time t0, const Array& x0,
Time dt) const = 0;
// same for diffusion etc.
};
virtual ~StochasticProcess() {}
virtual Size size() const = 0;
virtual Size factors() const {
return size();
}
virtual Disposable<Array> initialValues() const = 0;
virtual Disposable<Array> drift(Time t,
const Array& x) const = 0;
virtual Disposable<Matrix> diffusion(Time t,
const Array& x) const = 0;
virtual Disposable<Array> expectation(Time t0,
const Array& x0,
Time dt) const {
return apply(x0, discretization_->drift(*this,t0,x0,dt));
}
virtual Disposable<Matrix> stdDeviation(Time t0,
const Array& x0,
Time dt) const {
return discretization_->diffusion(*this,t0,x0,dt);
}
virtual Disposable<Matrix> covariance(Time t0,
const Array& x0,
Time dt) const {
return discretization_->covariance(*this,t0,x0,dt);
}
virtual Disposable<Array> evolve(Time t0,
const Array& x0,
Time dt,
const Array& dw) const {
return apply(expectation(t0,x0,dt),
stdDeviation(t0,x0,dt)*dw);
}
virtual Disposable<Array> apply(const Array& x0,
const Array& dx) const {
return x0 + dx;
}
}; このクラスでモデル化されている概念は、一般的な ｎ次元の確率過程で、次の式で表現されるものです。
\[ d{\bf x} ={\boldsymbol \mu }(t; {\bf x})・dt + {\boldsymbol \sigma }(t; {\bf x})・d{\bf w} \]（６Chapter目で初めて、この本の中に微分方程式が出てきました。私も、かつては物理学者でしたので、それらしく。）
しかし、このクラスのインターフェースは、時間を離散化した確率過程のサンプルPathを主に取り扱っています（訳注：上記式は連続時間の確率過程を示している）。
このクラスの最初の所で、polymorphic(多相的)な内部クラスとして discretizationクラスを定義しています。これは、Strategyパターンの一例で、確率過程の離散的なサンプルを取る方法に、様々なバリエーションを与える事ができるようにする為に、このパターンを使っています。discretizationクラス内のメソッド （訳注：drift( )や diffusion( )）は、StochasticProcessクラスのインスタンス、初期値 \( （t_0、x_0）\) 時間間隔 \( \Delta t\)、を引数として取り、離散化した確率過程のドリフト項と拡散項を返します。
（注： ここで、読者の方は、ジャンプ（訳注：jump diffusion process）を扱う為の仕組みを取り込んでないではないかと気づかれたと思います。わかっています。少し待って下さい。後で説明します。）
StochasticProcessのインターフェースに戻ります。このクラスではいくつかのインスペクター関数を定義しています。size( )メソッドは、その Processインスタンスがモデル化している確率変数の数を返します。例えば、通常の Black-Scholesモデルでは１ですし、ボラティリティー確率変動モデルでは２になります。あるいは、対象資産が N個あり、それぞれ一定の相関を持ってBlack-Scholes過程で変動するモデルであれば、Nを返します。factors( )メソッドは、確率変数を動かす確率変動成分の数を返します。このメソッドは、デフォールト値として size( )と同じ値を返しますが、後で考えてみるとそれで良かったか疑問です。デフォールト値を設定せずに、常に外生的に与えられる値を要求すべきであったかも知れません。最後のインスペクター関数は、initialValues( )メソッドで、t=0 時点の対象資産の値を返します。言い換えると、初期値であり、その値に確率的な性格はありません。StochasticProcessクラスの中の殆どのメソッドと同じように、これらのメソッドの返す値は Disposableクラステンプレートのインスタンスで受取ります。これは、配列や行列の型で値を返す場合に発生するメモリー上のコピーを避ける為です。興味のある方は Appendix Aを参照下さい。C++11に準拠したコンパイラーを使えるなら、これらの仕組みを “r-value参照” （下記[1]参照） に置き換える事ができます。しかし乍ら、QuantLibライブラリーの環境をC++11準拠の環境に移すには、あと数年かかるでしょう。（“環境を移す”の意味には、ユーザーを雨の中に放っておかずに古いバージョンのコンパイラーの環境でのサポートをストップする事も含まれます）
次の２つのメソッドは （ある意味インスペクター関数でもありますが）、Processインスタンスに関する、より詳細な情報を返します。drift( )メソッドは、上で示した方程式の中のドリフト項の係数 \( {\boldsymbol \mu }(t; {\bf x})\) を返す一方、diffusion( )メソッドは、拡散項の係数 \( {\boldsymbol \sigma }(t; {\bf x})\) を返します。 当然ながら、確率変動成分の数と確率変数の数を整合させる為、ドリフト項は N次元の配列で、拡散項は N × M 次元の行列になります 。M と N はそれぞれfactors( ) と size( )メソッドが返す値と同じになります。また、拡散項の行列には、相関の情報も含めなければなりません。
最後の方にあるメソッド群は、個別の Pathの離散化を扱います。expectation( )メソッドは、t 時の確率変数の値を \( \bf x\) （訳注： N個のベクトル）とした場合の、\( t+\Delta t\) 時における期待値 （訳注： N個の配列）を返します。デフォールトの実装は、内部クラスである discretizationのインスタンスを使って \( \Delta t\) の間の、ドリフト項\( \Delta \bf x \) を計算し、その値と初期値 \( \bf x\) を apply( )メソッドに渡しその計算結果を返します。（注： apply( )メソッドはもうひとつの仮想関数で、\( \bf ｘ\) と\( \Delta \bf x \)を引数で取り、単純に\( {\bf ｘ} + \Delta \bf x \) を返します。このような単純な計算の為に、なぜこのような多相的なメソッドが必要なのか疑問に思われるでしょう。これについては、悲しい話ですが、もう少し後で説明します。） stdDeviation( )メソッドは、\( (t, {\bf x})\) を初期値とする \( t + \Delta t \) 間の拡散係数 （訳注： N × M 行列） を返します。また covariance( )メソッドは、その値の2乗を返します。デフォールトの実装は、両方のメソッドとも、discretizationインスタンスにその計算を委託しています。最後に、evolve( )メソッドですが、やや抽象的なレベルのインターフェースを提供しています。このメソッドは、引数として確率変数の初期値 \( (t, {\bf x})\)、離散化した時間の間隔 \( \Delta t\)、そしてガウス分布する乱数（確率変動成分）の配列 \( \bf w\) を取り、終点時の確率変数の配列 \( \bf x^‘\) を返します。このメソッドのデフォールトの実装は、標準偏差の行列 （訳注：拡散項の係数行列 σ のこと）を取り込み、それを確率変動成分の行列に掛けてシミュレーションされた Pathの拡散項を導出し、それにドリフト項を加え、\( t + \Delta t \)時における確率変数の値を返します （その値には確率変動の中の確定項の部分も含みます）。
ここで見た通り、最後のメソッド群はすべてデフォールトの実装を備えていますが、時には discretizationストラテジーパターンを使って、その drift( )や diffusion( )メソッドを呼び出したり、別の時には同じグループの他のメソッドを呼び出したりしています。これは、多層レベルの Template Methodパターンの一種と言えます。この仕組みにより、StochasticProcessクラスの派生クラスにおいて、様々なレベルでの動作を実装する事が可能になっています。
そのような派生クラスでは最低限、ベースクラスで純粋仮想関数として宣言されている drift( )と diffusion( )メソッドを実装する必要があります。これが、確率過程の根幹部分を特定し、それだけで、稼働する StochasticProcessを生成するのに十分なパーツとなります。その他のメソッドはすべて、デフォールトの実装を使うか、内部クラスの discretizationインスタンスに計算を委託するようになっています。派生クラスでは個別の discretizationストラテジーを記述する事ができますが、何かデフォールトの実装をこしらえてもいいですし、ユーザー自身の選択にすべて任せてしまっても構いません。
何段階か外の階層の派生クラスのレベルで、StochasticProcessクラスの expectation( )と stdDeviation( )メソッドをオーバーライドする事も可能です。例えば、計算を解析解で導出できるような場合には、そのようなプログラムコードを書く事も出来ます。その場合、そのクラスの開発者にとっての選択肢は、コンストラクターに discretizationインスタンスの生成をさせず、解析解の式を必ず使うようにさせるか、あるいはコンストラクターにそのインスタンスを持たせるものの、オーバーライドするメソッドが呼ばれた時にそのインスタンスの有無をチェックし、仮に存在していれば自分自身のインスタンスを使わずベースクラスの実装に動作を投げる、という風に決める事ができます。確率変動成分から確率変数の変動を記述する動作や、ドリフト項と拡散項の変動の累計を合算する動作は、引き続き evolve( )メソッドに託されています。
クラス階層の最も外縁のレベルでは、StochasticProcessクラスが提供している大半のメソッドをオーバーライドした上で、evolve( )メソッドをオーバーライドする事も可能です。この方法は、今のクラス階層の中で実装されているインターフェースが全く使えないような確率過程を取り扱う際の、Hail Maryパス （訳注：マリア様お願いパス。アメリカンフットボールにおいて、終了間際の逆転を狙ったロングパス。最後の手段） のようなものです。例えば、確率過程の中に確率的なジャンプを加えたい場合に （訳注：jump diffusion modelを使う場合）、そのような確率変動成分を evolve( )メソッドに加える事も可能です。但しその場合、疑似乱数生成 RNG がガウス分布の乱数とポアソン分布の乱数を適切に組み合わせた乱数を生成して、そのメソッドに渡す必要があります。クラスのインターフェースの中に、そのような jumpする拡散項を扱えるものを加えた方がいいでしょうか？もし開発者がそのインターフェースの周りではなく、それそのものを使って作業をできるのであれば、その方がいいかも知れません。しかし、現時点で StochasticProcessクラスのインターフェースに jump する拡散項を扱うメソッドを加えるのは２つの問題があると思います。
ひとつは、適切なインターフェースがどのようなものか、私自身、確信がありません。例えば、jumpを取り込んだ一般的な計算式を1つ用意すればいいのか？後で正しくなかったと言われる可能性のあるインターフェースを今の段階でリリースする事はしたくありません。現段階では、それを汎用的に使えるようなプログラムコードを持っておらず、そうなる可能性が大です。
もうひとつは、次のバージョンと統合しようとする時、どうなるのでしょうか？インターフェースにそれを加えるのでしょうか？または、逆にスペックの詳細化ではなく汎用化をすべきなのでしょうか？すなわち、ベースクラスのインターフェースの数を増やすのではなく、減らす方法にすべきでしょうか？例えば、evolve( )メソッドの様なものはベースクラスに残すものの drift( )や diffusion( )の様なメソッドは派生クラスに移した方が良いのでしょうか？
（注： ところで、次のような方法が、Jump Diffusionのインターフェースを、コミットせずに取り込むにはいいかも知れません。すなわち、新しいメソッドとオーバーライドした evolve( )メソッドを持つクラスを実験的に作ってみて、動作するかどうか試してみれば、何等かの有用な情報が得られるかも知れません。）
少し、脇道に逸れ過ぎました。端的に言うと、インターフェースの中に Jump Diffusionを取り扱えるメソッドが無いのは、良くない事です。しかし、今それを加えようとすると、出来の悪いものしか出来ないでしょう。私としては、誰かが作動するプログラムコードを提供してくれるのを、窓の横に座って待っています。もしどなたか、助けて頂いてそのような実験的なクラスを提供して頂けるなら、大変感謝します。
いくつかの具体例に進む前に、最後の寄り道です。一次元の確率過程を実装する場合、Arrays（配列を扱うクラス）やそれに似たものが、若干邪魔になります。これらのデータ型を使うと、プログラムコードを書いたり読んだりするのが、より大変になります。理想的には、StochasticProcessクラスの中に、配列や行列の値ではなく、単なる実数値を入出力するようなメソッド群を加えた方がよかったかも知れません。しかし、そのようなメソッド群を既存の StochasticProcessクラスの関数をオーバーロードするような形で宣言する事はできません。そうしてしまうと、それらのメソッドは通常の Processに利用できなくなるからです。また、一次元の確率過程だけの為に別の独立した階層構造をつくるのもしたくありません。なぜなら、一次元の確率過程も StochasticProcessであるからです。
QuantLibライブラリーでは、次の Listing 6.5に示すような StochasticProcess1Dクラスを提供する事により、その問題を解決しています。
Listing 6.5: Partial implementation of the StochasticProcess1D class.
class StochasticProcess1D : public StochasticProcess {
public:
class discretization {
public:
virtual ~discretization() {}
virtual Real drift(const StochasticProcess1D&,
Time t0, Real x0, Time dt) const = 0;
// same for diffusion etc.
};
virtual Real x0() const = 0;
virtual Real drift(Time t, Real x) const = 0;
virtual Real diffusion(Time t, Real x) const = 0;
virtual Real expectation(Time t0, Real x0, Time dt) const {
return apply(x0, discretization_->drift(*this, t0,
x0, dt));
}
virtual Real stdDeviation(Time t0, Real x0, Time dt) const {
return discretization_->diffusion(*this, t0, x0, dt);
}
virtual Real variance(Time t0, Real x0, Time dt) const;
virtual Real evolve(Time t0, Real x0,
Time dt, Real dw) const {
return apply(expectation(t0,x0,dt),
stdDeviation(t0,x0,dt)*dw);
}
private:
Size size() const { return 1; }
Disposable<Array> initialValues() const {
return Array(1, x0());
}
Disposable<Array> drift(Time t, const Array& x) const {
return Array(1, drift(t, x[0]));
}
// ...the rest of the StochasticProcess interface...
Disposable<Array> evolve(Time t0, const Array& x0,
Time dt, const Array& dw) const {
return Array(1, evolve( t0, x0[0], dt, dw[0]));
}
}; 一方で、このクラスは、StochasticProcessクラスのインターフェースをそのまま写したインターフェースを宣言しています。デフォールトの実装までも、同じ内容ですが、（変数の取扱いをArrayではなく）実数の型を使っています。他方で、このクラスは StochasticProcessクラスから派生させています （１次元の確率過程であっても Stochastic Processの一種である”is a”という事を示したかったのです）。そして新しく定義されたインターフェースは、すべて Arrays（配列）を使っていた部分を Real（実数）に入れ替えており、その結果このクラスを使って開発を行う際に配列を気にする必要がなくなりました。端的に言えば、Adapterパターンをそっくりそのまま使ったという事です。スカラー値とベクトル値を使ったメソッド群は、完全に1対1対応するので （上記 Listingを見ての通りです）、1次元の確率過程は、いずれの（１次元あるいは多次元のクラスの）インターフェースを使っても同じように作動します （訳注：多次元のインターフェースの配列変数の要素の数を１とすれば同じ動作をする）。また、このクラスも多次元の StochasticProcessクラスと同じように、カスタマイズが可能になっています。
[1] H.E. Hinnant, B. Stroustrup and B. Kozicki, A Brief Introduction to Rvalue References, C++ Standards Committee Paper N2027, 2006.
***********************************************************************************
今から1つ、いえ３つの例を示します。
当然ながら、最初の例は、一次元の Black-Scholes過程です。読者の方は既に Chapter II で、ヨーロピアンオプションの PricingEngine に渡される引数の1つとして、これに対応するクラスに出会っていました （GeneralizedBlackScholesProcessという長ったらしい名前でした）。その時使ったプログラムコード（の問題点）については後程説明しますが、その時に、現在のデザインが持ついくつかの欠点と将来的な解決策についても議論します。ここではまず、現在の実装について見てみましょう。概要を下記の Listing 6.6に示します。
Listing 6.6: Partial implementation of the GeneralizedBlackScholesProcess class.
class GeneralizedBlackScholesProcess
: public StochasticProcess1D {
public:
GeneralizedBlackScholesProcess(
const Handle<Quote>& x0,
const Handle<YieldTermStructure>& dividendTS,
const Handle<YieldTermStructure>& riskFreeTS,
const Handle<BlackVolTermStructure>& blackVolTS,
const boost::shared_ptr<discretization>& d =
boost::shared_ptr<discretization>());
Real x0() const {
return x0_->value();
}
Real drift(Time t, Real x) const {
Real sigma = diffusion(t,x);
Time t1 = t + 0.0001;
return riskFreeRate_->forwardRate(t,t1,Continuous,...)
- dividendYield_->forwardRate(t,t1,Continuous,...)
- 0.5 * sigma * sigma;
}
Real diffusion(Time t, Real x) const;
Real apply(Real x0, Real dx) const {
return x0 * std::exp(dx);
}
Real expectation(Time t0, Real x0, Time dt) const {
QL_FAIL("not implemented");
}
Real evolve(Time t0, Real x0, Time dt, Real dw) const {
return apply(x0, discretization_->drift(*this,t0,x0,dt) +
stdDeviation(t0,x0,dt)*dw);
}
const Handle<Quote>& stateVariable() const;
const Handle<YieldTermStructure>& dividendYield() const;
const Handle<YieldTermStructure>& riskFreeRate() const;
const Handle<BlackVolTermStructure>& blackVolatility() const;
const Handle<LocalVolTermStructure>& localVolatility() const;
}; まず最初に、クラス名の長さについて説明すべきでしょう。”generalized”の部分は、このモデルの拡張性を示したものではありません。いくつかの同種の確率過程もカバーする事を意図していたという事です（その様な過程としては、本来の Black-Scholesプロセス、Black-Scholes-Mertonプロセス、Blackプロセスなどです）。短いクラス名は、それら細分類をカバーする為に残していました。
コンストラクターは、汎用的なオブジェクト構築に必要なすべての引数を Handle< >の型で取ります。即ち、リスクフリー金利の Term Structure、原資産の市場価格の Quote、配当利回りの Term Structure、そして Implied Black Volatilityの Term Structure、です。また、StochasticProcessクラスのインターフェースを実装する為、discretizationクラスのインスタンスも引数として取り込んでいます。金融商品の価格モデルに必要なすべての入力項目はコンストラクターによって生成されたインスタンスに保持され、インスペクター関数を使ってその内容を取り出す事が出来ます。
コンストラクターが、Black Volatility \( \sigma_B (t､k)\) を引数として取っているのに気が付かれたでしょうか。Black Volatilityは、拡散項の係数 \( \sigma (t､k)\) を返すために、直接その値が使われる訳ではありません（このクラスは、Observerパターンを使って、Black Volatilityから、局所的な Volatilityへの換算のタイミングを調整しています。再計算は Lazyに行われます。即ち、Observableのデータが更新されても直ちに再計算は行われず、局所的な Volatilityが必要となった場合にのみ再計算します）。Processクラスは、引数として渡された Black Volatilityの種類 （定数なのか、時間依存なのか、あるいは Volatility Smile があるのか）に合わせて、それぞれ異なるアルゴリズムを使って、その値を必要な形に変換しています。データ型は、dynamic cast演算子を使ってRun-timeでチェックされます。あまりスマートな方法ではありませんが、このケースでは Visitorパターンを使うよりは簡単です。今の所、変換するアルゴリズムは Processクラスの中にすでにハードコードされており、ユーザーによって内容を書き換える事はできません。あるいは事前計算された局所的な Volatility を提供する事も出来ません。
次のメソッド、drift( )と diffusion( )はこのクラスの根幹部分になります。このプロセスは、実際には原資産の価格そのものの確率過程ではなく、価格の対数の確率過程を表しています。従って、drift( )メソッドが返す値は、対数正規分布におけるドリフト項にあたり、式で表すと \( r(t) - q(t) – 0.5 \times \sigma^2(t､x)\) になります。この式の中にある２種類のレート （リスクフリー金利と配当利回り）及び局所 Volatilityの値は、それぞれ対応する Term Structureから引っ張ってきます （この点が、もうひとつの問題点ですが、後で説明します）。diffusion( )メソッドは局所 Volatilityの値を返します。
確率変数として価格の対数を選択した事により、いくつか問題となる影響が出ています。最初に、このクラスが取り扱う確率変数は、ベースクラスのインターフェースが外部とデータのやり取りを行う際に想定している値とは違ってきます。コンストラクターの引数で取る Quoteの値は、原資産の市場価格で、x0( )メソッドが返す値も原資産の価格を想定しています。例えば evolve( )などの他のメソッドでも、原資産の価格そのものを想定したデータのやり取りをしています。今振り返ると、確率変数に価格の対数を選択したのは、まずい考えでした （読者の方の“あ～あ”という声が聞こえます）。自分達ですでに開発したデザインに引きずられてこうなりました。しかし重要なピースが抜けているという事に気づくべきでした。
その影響は、これまで示したプログラムコードの Listingの中の至る所に出ています。価格の対数を確率変数として選択したことにより、（ベースクラスで） apply( )メソッドを仮想関数にする必要があり、かつそのメソッドが存在している理由そのものです。確率変数 x にドリフト項の変化 Δ を加える計算は、（価格の対数を使った場合は） x にそれを足すのではなく、\(x \times exp(\Delta )\) となります。expectation( )メソッドも同様に影響を受けます。このメソッドのデフォールトの実装は、apply()と同様 x に Δ を足すものですが、これも \( E[x]\) を返すのではなく \( exp(E[log(x)])\) を返す必要があります。このメソッドが間違った値を返さないようにする為、今の所このメソッドを使えない状態にしています（使おうとすると例外処理に飛びます）。その結果 evolve( )メソッドは、デフォールトの実装では expectation( )メソッドに依存していましたが、（それが使えない為）別の方法を使ってオーバーライドしなければなりませんでした。
（注： evolve( )メソッドは、計算速度の問題から、いずれにしてもオーバーライドする必要がありました。StochasticProcess1D でのデフォールトの実装は、\( x \times exp(\mu dt) \times exp(\sigma dt )\) を返すようになっていましたが、オーバーライドしたメソッドでは \( x \times exp(\mu dt + \sigma dt )\) を返すようになっています。（訳注：指数計算の演算は時間がかかるので、デフォールトでは２回行うところを、オーバーライドしたメソッドでは1回で済む））
上記のような、まずいデザインに加え、計算速度の問題も抱えていました。計算速度の問題については、QuantLibが提供しているモンテカルロ法を使った PricingEngineを試した事のあるユーザーや、さらには Wilmottフォーラム （訳注： Wilmottは英国にある金融工学に関する教育を提供している会社 https://www.wilmott.com/）の参加者すべて に知れ渡っていました （この人たちからは他にも様々なコメントを頂きました）。 evolve( )メソッドは、その中で、（指数計算を行う）apply( )メソッドを使っているので、離散化した時間の各ステップで、指数関数を計算する事になってしまいます。最悪なのは、drift( )と diffusion( )メソッドは、その時間ステップ毎に、繰り返し Term Structureインスタンスに値を求めています。Chapter III での説明を覚えておられるのなら、これは、少なくとも２階層以上の仮想関数の呼び出しのプロセスを経由しており、時には（非常に運が悪ければ）ゼロ金利の Term Structureから生成された Discount Factor の Term Structureを使って、そこからゼロ金利を取り出すような不要な手間をかけている事もあり得ます。
最後にもうひとつ、デザイン上の問題があります。前にも指摘していましたが、Chapter II での具体例で、AnalyticEuropeanEngineクラスの中で、この Black-Scholes過程を既に使っていました。Libraryにあるその PricingEngineのプログラムコードを見ればわかる通り、その Engineは StochasticProcessクラスを本来の目的通りに使っていません。そこでは StochasticProcessインスタンスを、市場価格のQuoteとTerm Structureのデータの入れ物として使っていました。この事は、StochasticProcessクラスをよろず屋のように使っている事を示しています。
さて、これをどのように修正して行きましょうか？おそらく、設計のやり直しが必要でしょう。このChapterを書く為にコードを見直していく中で、いくつかのアイデアが浮かびました。だいたい、次のようなものです。但し、現段階では仮説程度のものですし、実際のプログラムコードも出来ていませんが。
まず最初に（また再び、Overengineering（凝りすぎ）との批判を受けるのを覚悟の上で）、このクラスの‘確率過程を扱っている部分’と、Term Structureのような‘データを保存している部分’を分けます。後者については、おそらく BlackShcolesModelのようなクラス名を付けて、新しいクラスを作ります。このクラスのインスタンスが、Quoteや Term Structureなどのデータをすべて保持し、Pricing Engineに渡されてそこで使われる事を想定しています。それを目的としながらも、そのインターフェースは、ユーザーが Strategyパターンを使って、（デフォールトの方法で）局所Volatilityを提供したり、あるいはユーザーが独自に開発した方法で Black Volatilityを局所 Volatilityに変換し、それを提供できるようにしたりします。
その Model（クラス）から、StochasticProcessクラスを構築できるようにします。Modelインスタンスを引数で取り、StochasticProcessインスタンスを返すような Factoryクラスを作ってもいいかも知れません （クラス間の協働の度合いやクラスの数をどの程度にしたいかによりますが）。あるいは Modelクラスの中に、StochasticProcess返す Factoryメソッドをいくつか加えてもいいです。または、Modelクラスを引数で取る StochasticProcessクラスのコンストラクターを作ってもいいかもしれません。いずれにしても、こうする事によって、確率過程のシミュレーションをより効率的に行える StochasticProcessクラスを作ることが出来るでしょう。例えば、（どのような設計にしても）その Factoryは、シミュレーションの為の離散時間の Nodeを取り、その Node毎に事前に計算したレート ( \( r(t_i)とq(t_i) \) ) とともに StochasticProcessインスタンスを返します （レートは原資産の価値に依存しないので事前計算が可能です）。もし型チェックで Volatility Smileを持たないとわかっていれば、StochasticProcessインスタンスは \(\sigma (t_i)\) も事前計算できるようにしてもいいでしょう （訳注： Volatility Smileは、対象資産価格とStrikeの距離によってVolatilityの値が変わるので、導出するには対象資産価格が必要になり、事前計算できない）。
確率変数が原資産価格 x ではなく log(x) である点については、この StochasticProcessクラスを書き直して、明白に対象資産価格の対数を扱う StochasticProcessクラスとして設計し直した方が良いでしょう。x0( )メソッドは、価格の対数を返すようにし、expectation( )や evolve( )などのその他のメソッドについても、そのようにします。そうすれば、このクラスは整合性と明確さを取り戻せるでしょう。しかし、それによってこのクラスが原資産の価格の Pathを生成できなくなるので、一般的な Processクラスを使った Clientコードでは、確率変数の値（log(x)）を原資産の価格(x)に変換する関数あるいは関数オブジェクトが必要になるかもしれません。そのような関数は、StochasticProcessクラスのインターフェースに加えるか、Clientコードに任せるかでしょう。
***********************************************************************************
実際のプログラムコードに戻りましょう。行儀のいい StochasticProcessクラスの例として、OrnsteinUhlenbeckProcessクラスについて、下記 Listing 6.7を見てみましょう。
Listing 6.7: Implementation of the OrnsteinUhlenbeckProcess class.
class OrnsteinUhlenbeckProcess : public StochasticProcess1D {
public:
OrnsteinUhlenbeckProcess(Real speed,
Volatility vol,
Real x0 = 0.0,
Real level = 0.0);
Real x0() const;
... // other inspectors
Real drift(Time, Real x) const {
return speed_ * (level_ - x);
}
Real diffusion(Time, Real) const {
return volatility_;
}
Real expectation(Time t0, Real x0, Time dt) const {
return level_ + (x0 - level_) * std::exp(-speed_*dt);
}
Real stdDeviation(Time t0, Real x0, Time dt) const {
return std::sqrt(variance(t,x0,dt));
}
Real variance(Time t0, Real x0, Time dt) const {
if (speed_ < std::sqrt(QL_EPSILON)) {
return volatility_*volatility_*dt;
} else {
return 0.5*volatility_*volatility_/speed_*
(1.0 - std::exp(-2.0*speed_*dt));
}
}
private:
Real x0_, speed_, level_;
Volatility volatility_;
};
Ornstein-Uhlenbeck過程は、シンプルな確率過程であり、その特徴的な動きは、中心回帰するドリフト項 \( \theta (\mu -x)\) と定数となる拡散項 σ が、そのまま（時間軸に沿って）積分していく事が出来るという事です。従って、必ずオーバーライドする必要がある drift( ) と diffusion( )メソッドに加え、expectation( ) と stdDeviation( )メソッドについても、この過程が記述する確率変動の値をそのまま返すようにオーバーライドしています。variance( )メソッドは （そのメソッドを使ってstdDeviation( )が実装されています）、数値積分を行っていく際に起こり得る数値的不安定性を防ぐため、２パターンに分岐するようになっています。すなわち θ が非常に小さい場合は、varianceの計算は θ→0 となる極限の値を返すようになっています。
***********************************************************************************
最後に、複数の確率変動成分を持つ確率過程の例を示します。下記 Listing 6.8に示す、StochasticProcessArrayクラスを見てみましょう。
Listing 6.8: Partial implementation of the StochasticProcessArray class.
class StochasticProcessArray : public StochasticProcess {
public:
StochasticProcessArray(
const std::vector<shared_ptr<StochasticProcess1D> >& ps,
const Matrix& correlation)
: processes_(ps), sqrtCorrelation_(pseudoSqrt(correlation)) {
for (Size i=0; i<processes_.size(); i++)
registerWith(processes_[i]);
}
// ...
Disposable<Array> drift(Time t, const Array& x) const {
Array tmp(size());
for (Size i=0; i<size(); ++i)
tmp[i] = processes_[i]->drift(t, x[i]);
return tmp;
}
Disposable<Matrix> diffusion(Time t, const Array& x) const {
Matrix tmp = sqrtCorrelation_;
for (Size i=0; i<size(); ++i) {
Real sigma = processes_[i]->diffusion(t, x[i]);
std::transform(tmp.row_begin(i), tmp.row_end(i),
tmp.row_begin(i),
bind2nd(multiplies<Real>(),sigma));
}
return tmp;
}
Disposable<Array> expectation(Time t0, const Array& x0,
Time dt) const;
Disposable<Matrix> stdDeviation(Time t0, const Array& x0,
Time dt) const;
Disposable<Array> evolve(Time t0, const Array& x0,
Time dt, const Array& dw) const {
const Array dz = sqrtCorrelation_ * dw;
Array tmp(size());
for (Size i=0; i<size(); ++i)
tmp[i] = processes_[i]->evolve(t0, x0[i], dt, dz[i]);
return tmp;
}
private:
std::vector<shared_ptr<StochasticProcess1D> > processes_;
Matrix sqrtCorrelation_;
}; このクラスは、特定の確率過程をモデル化したものでは無く、互いに相関を持つN個の１次元確率過程を統合したものです。このクラスを使って、相関係数（行列）が、どのようにして確率過程の中に取り込まれていくのかを示します。
コンストラクターは、複数の1次元確率変数の配列と、それらの相関行列を引数として取ります。確率変数の配列は、メンバー変数に保存されますが、相関行列は保存されません。その代わり、その相関行列の平方根を取り、それをメンバー変数に保存します （行列Aの平方根とは、\( {\bf XX^T ＝ A}\) となるような行列 X ）。また、コンストラクターは、自らを Observerとして、引数で取った各StochasticProcess1Dインスタンスに登録し、データの変更があった場合に通知を受けられるようにします。
initialValues( )、 drift( )、expectation( )などのメソッドは、引数で取った各StochasticProcess1Dインスタンスの該当する同じメソッドをLoopの中で呼び出し、計算結果を配列にまとめるだけです。diffusion( )メソッドも同じく各ProcessインスタンスのメソッドをLoopの中で呼び出しますが、それらの計算結果を、相関行列を使って変換します。すなわち、相関行列の平方根に、対応する各確率過程の拡散項の値をかけます（それに、さらに拡散項の転置行列をかけると、馴染みのある共分散行列 \( σ_iρ_{ij}σ_j \) になります）。stdDeviation( )メソッドも同様の計算を行いますが、\( \Delta t\) を含む各確率過程の標準偏差の行列を使います。
残っているメソッドは evolve( )です。もし、すべての StochasticProcessインスタンスが想定通りに動くのなら （つまり、確率変数（の初期値）に計算された変動部分（ドリフトと拡散を足しているだけなら）、デフォールトの実装をそのまま継承すればいいだけです。すなわち expectation( )と stdDeviation( )の計算結果を取りだし、後者の値 （訳注：多変数なので行列の形で存在） に確率変動成分の配列を掛け、それにドリフト項の値を足せばいいだけでした。しかし、そのような保証はありませんし、evolve( )メソッドは別の方法で実装する必要があります。まず、ガウス分布する確率変動成分（の配列）に相関行列の平方根を掛け、相関する確率変動成分の配列を算出します。その値を使って(引数として)、各 StochasticProcess1Dインスタンスの evolve( )メソッドを呼び出して、その配列メンバー毎の計算結果を集計します。
<ライセンス表示>
QuantLibのソースコードを使う場合は、ライセンス表示とDisclaimerの表示が義務付けられているので、添付します。 ライセンス