２. "Implementing QuantLib"の和訳
Appendix A Odds and Ends 周辺の話題
A.4 Math-related Classes 数値計算関連のクラス群 （続き）
A.4.4 Statistics 統計値
Statisticsクラス群は、Monte Carloシミュレーションで生成されたサンプルを集計する目的で作られました（Chapter VIでの説明を思い出して下さい）。このクラスが提供するすべての機能は、（Decoratorパターンを使って）一連の decorator で実装されており、各 decorator が、階層化されてメソッドを追加するようになっています。階層化が多重になっている理由は、私の推察では、decoratorの基本階層になるベースクラスが２つ用意されている為であり、その内のいずれかを選択できるようにする為かと思われます。階層化された設計方法により、ユーザーは、基本階層の共通インターフェースに基づいて、機能を追加していく事が出来ます。
ユーザーが選択できる基本階層のひとつは IncrementalStatistics クラスで、そのインターフェースを下記 Listing A-16 に示します。このクラスがどの様なインターフェースを提供しているかは簡単に推察がつくと思いますが、例えば、サンプル数を返したり、加重されたサンプルの場合はその加重合計であったり、その他様々な統計値を返したりします。
Listing A-16：IncrementalStatisticsクラスのインターフェース
class IncrementalStatistics {
public:
typedef Real value_type;
IncrementalStatistics();
Size samples() const;
Real weightSum() const;
Real mean() const;
Real variance() const;
Real standardDeviation() const;
Real errorEstimate() const;
// skewness, kurtosis, min, max...
void add(Real value, Real weight = 1.0);
template <class DataIterator>
void addSequence(DataIterator begin, DataIterator end);
};サンプルを追加する場合、1つずつ追加する事もできますし、iteratorで指定された配列として追加する事もできます。このクラスの特徴は、引数として渡されたデータを保持せず、そのまま統計値を計算する事です。これは、メモリーを節約するのが目的でした。（Libraryのプロジェクトが始まった）2000年頃は PC の RAMのサイズは 128MB から 256MBであった為、メモリー確保は今よりもずっと大きな関心事でした。このクラスで実際の計算を行うプログラムコードは、かつては独自に開発されたものでしたが、今では Boost の accumulator libraryを使っています。
ふたつ目の基本階層は、GeneralStatisticsクラスで、IncrementalStatisticsクラスと同じ名前を持つインターフェースに加え、いくつかの独自のメソッドが加わっています。このクラスは、渡されたデータを一旦保持し、その値を返す事が出来るので、追加されたメソッドはそうした機能を提供するものです。例えば、データのパーセンタイルを返す事や、データの並び替えが可能です。このクラスは、expectationValue( )というテンプレートメソッドを提供しており、ユーザーが独自に開発した計算アルゴリズムを使う事が出来ます。興味のある方は、この Sectionの最後にある“aside”の欄をご覧下さい。
Listing A17：GeneralStatisticsクラス
class GeneralStatistics {
public:
// ... same as IncrementalStatistics ...
const std::vector<std::pair<Real,Real> >& data() const;
template <class Func, class Predicate> std::pair<Real,Size>
expectationValue(const Func& f, const Predicate& inRange) const;
Real percentile(Real y) const;
// ... other inspectors ...
void sort() const;
void reserve(Size n) const;
};次に、外側 （訳注：DecoratorパターンにおけるDecorateする側）の階層について説明します。最初のクラスは、Expected ShortfallやValue at Riskといった、リスク量計算における統計値を計算できるクラスです。（訳注： Value at Riskは、2σ、99%信頼区間、99%パーセンタイル、といった、ある一定の閾値における損失額。Expected Shortfallは閾値を超える損失について、その発生確率をかけた条件付き期待値。） これらの値を計算する際の問題点は、サンプルデータのすべてを必要とすることです。IncrementalStatisticsクラスの場合、サンプルデータを保持していないので、正確な計算はあきらめて、近似値を計算する事になります。その方法のひとつは、標本集団がガウス分布と同じモーメントを持つと仮定して、サンプルの平均と分散を使い、解析的に計算する方法です。下記の GenericGaussianStatisticsクラスは、まさにそれを行っているクラスです。
Listing A18：GenericGaussianStatisticsクラスのインターフェース
template<class S> class GenericGaussianStatistics : public S {
public:
typedef typename S::value_type value_type;
GenericGaussianStatistics() {}
GenericGaussianStatistics(const S& s) : S(s) {}
Real gaussianDownsideVariance() const;
Real gaussianDownsideDeviation() const;
Real gaussianRegret(Real target) const;
Real gaussianPercentile(Real percentile) const;
Real gaussianValueAtRisk(Real percentile) const;
Real gaussianExpectedShortfall(Real percentile) const;
// ... other measures ...
};
typedef GenericGaussianStatistics GaussianStatistics; すでに述べた通り、またコードを見ても判りますが、このクラスはベースクラスに対する Decorator として実装されています。Decorateされる対象（ベースクラス）をテンプレート引数として取り、そこから派生するクラスを定義する事によって、ベースクラスのすべてのメソッドが使えた上で、さらに新しいメソッドが追加できます。QuantLibでは、テンプレート引数を特定して具体化されたクラス型として、GaussianStatisticsクラスを提供していますが、テンプレート引数として GeneralStatisticsクラスを指定しています。おかしいですね。テンプレート引数として IncrementalStatisticsクラスの方が指定されるべきですが、Libraryの中に見当たりませんでした。これには理由があります。少しお待ちください。
ベースクラスがサンプルデータをすべて保持しているのであれば、（分布から解析的に計算されるリスク量では無く）データから実際のリスク量を計算する Decoratorクラスを作る事ができるはずです。それが、下記 Listing A-19 にある GenericRiskStatisticsクラスになります。ここではガウス分布に従った統計量計算アルゴリズムの実装内容については説明しません（各自Libraryのコードで確認して下さい）。
Listing A-19： GenericRiskStatisticsクラスのインターフェース
template <class S> class GenericRiskStatistics : public S {
public:
typedef typename S::value_type value_type;
Real downsideVariance() const;
Real downsideDeviation() const;
Real regret(Real target) const;
Real valueAtRisk(Real percentile) const;
Real expectedShortfall(Real percentile) const;
// ... other measures ...
};
typedef GenericRiskStatistics<GaussianStatistics> RiskStatistics; コードを見れば判る通り、Decoratorの階層が重ねられています （訳注：コードの最後にある typedef を使ったテンプレートクラスの具体化の部分）。すなわち、ライブラリーでは、GenericRiskStatistics のテンプレート引数に GaussianStatisticsクラスを指定して、RiskStatisticsという具体的クラスを定義しています （訳注：それ自体が Decoratorである GaussianStatisticsクラスを、さらに Decorateするクラスを定義している）。従って、このクラスは、ガウス分布を前提にして解析的に導出した統計値に加え、実際のデータを使って導出した統計値も算出できます。これが、GaussianStatistics の Decorate対象として GeneralStatisticsクラスが選択された理由のひとつでした。
以上に加えて、さらに別の Decorator を重ねる事も可能です。QuantLibライプラリーの中でも、追加の Decorator がいくつか提供されていますが、ここではプログラムコードは示しません。そのひとつに、SequenceStatisticsクラスがあります。このクラスは、標本集団がひとつでは無く、（標本集団の）配列として存在している場合に使われます。このクラスは、内部的にはスカラーの Statisticsインスタンスの配列を使っており、各標本集団間の相関や共分散の計算アルゴリズムが追加されています。このクラスは、例えば、LIBOR Market Modelのような、異なる LIBOR期間のキャッシュフローについて、それぞれに標本集団が必要な場合に使われます。その他には、ConvergenceStatisticsクラスや DiscrepancyStatisticsクラスがあり、これらは、標本集団の配列に関する情報を計算して提供しています。このクラスは Libraryの他の所では使われていませんが、少なくとも単体テストはやりました。
< Asid : Extreme Expectations (極端な期待)
GeneralStatisticsクラスを見直してみて、このコードを自慢すべきなのか、恥ずべきなのか、よく判りません。なぜなら、当時私は、このクラスで取り入れた汎用的な設計について、うまくいったと浮かれていました。
それは演算方法の数学的な抽象化でした。当初は、（様々な統計値の計算アルゴリズムを）それぞれシンプルに実装するだけでした。しかし、（様々な Decoratorの階層で実装された） 平均、分散、さらにもっと複雑な計算アルゴリズム （訳注；おそらく３次や４次のモーメントである Kurtosis や Skew の計算の事）を見てみると、すべて同じような形の次のような式で表せると気付きました。
\[ \frac{\sum_{x_i∈R} f(x_i)ω_i}{\sum_{x_i∈R}ω_i} \]すなわち、この式は領域Rにおける関数 \(f(x)\) の期待値を計算しています。この式を使えば、平均の計算は、\(f(x)\) に identity関数 （単に x を返す関数） で領域 R はサンプル集合全体を指定すればいいでしょう。分散の計算は、\(f(x)\) に \((x-\tilde x)^2\) （但し \(\tilde x\) は標本の平均）を使い、同じ領域で計算すれば求まるでしょう。同様に、\(f(x)\) の内容を変えるだけでより高次のモーメントも求まります。その結果、expectationValue( )メソッドをテンプレートメソッドにし、引数として関数 \(f(x)\) と領域 R を渡し、計算結果とサンプル数を返すようにしました。（平均や分散値を使う）他のメソッドは、適当な引数を渡してこのメソッドを呼び出すように実装されています。例えば、平均値の計算は、次のように実装されており、一見しただけでは何を計算しようとしているのか判りにくいかも知れません。
return expectationValue(identity(), everywhere()).first; ところで、当時、私は functional programming について、もっと学んでおくべきでした。サンプルの有効領域は関数を使って渡されており、（訳注；第2引数のeverywhere( )がそれ）、各サンプルが領域の有効範囲内かどうかで true または false を返します。上記コード中の everywhere( )は、その機能を持った Helper関数のひとつです。このように（統計値の計算アルゴリズムを抽象化・一般化できた事は）うまくいったと感じていました。誰かが \((x-\tilde x)^2\) を次のように書くまでは。
compose (square(), std::bind2nd(std::minus(), mean())) なるほど、このような書き方は優れて汎用的で、ユーザーが新しい計算方法をプログラムしやすくなります。しかし、やや数学的になりすぎて、意味がわかりにくいプログラムコードでもあります。私の方法の方が、もしかしたら汎用性と解りやすさのバランスがうまく取れているかも知れません。C++11に準拠して上記のbind関数をラムダ関数に置き換えれば、もう少し解りやすくなるかも知れません。将来、そちら （C++11） に移行していく際に、どのようにコードが変わっていくか、見てみたいと思います。
一方で、計算速度については、大きな問題になっていません。expectationValue( )や compose( )、everwhere( )メソッドはテンプレートメソッドになっており、コンパイラーが、適切に in-line化するなど、効率的にコンパイルされるはずです。そうすると、コンパイル後の計算アルゴリズムは、平均や分散の計算アルゴリズムを直接プログラムで書いた場合と同じ様な、より単純なループを使った計算方法になっているはずです。
A.4.5 Linear Algebra： 線形代数（ベクトルと行列）
下記に Arrayクラスと Matrixクラスの概要を示しますが、これらについて説明すべき事は、それ程ありません。
Listing A-20：ArrayクラスとMatrixクラスの概要
class Array {
public:
explicit Array(Size size = 0);
// ... other constructors ...
Array(const Array&);
Array(const Disposable<Array>&);
Array& operator=(const Array&);
Array& operator=(const Disposable<Array>&);
const Array& operator+=(const Array&);
const Array& operator+=(Real);
// ... other operators ...
Real operator[](Size) const;
Real& operator[](Size);
void swap(Array&);
// ... iterators and other utilities ...
private:
boost::scoped_array<Real> data_;
Size n_;
};
Disposable<Array> operator+(const Array&, const Array&);
Disposable<Array> operator+(const Array&, Real);
// ... other operators and functions ...
class Matrix {
public:
Matrix(Size rows, Size columns);
// ... other constructors, assignment operators etc. ...
const_row_iterator operator[](Size) const;
row_iterator operator[](Size);
Real& operator()(Size i, Size j) const;
// ... iterators and other utilities ...
};
Disposable<Matrix> operator+(const Matrix&, const Matrix&);
// ... other operators and functions ... このクラスが提供しているインターフェースは、皆さんの予想通りのものです。すなわち、コンストラクターと代入演算子(=)、配列要素へのアクセス演算子（注： Arrayクラスは a[i] という記法を使い、Matrixクラスでは、ユーザーにより便利なように、m[i][j] と m(i,j) という２種類の記法を用意しています）、さらに配列や行列の要素毎に行われる算術演算子とユーティリティーメソッドです （注： a * b という算術演算子が、ドット積ではなく、要素毎の積を返すように設計されている点は、どうも気になります。でも、それは私だけのようで、他のプログラマー達は気にならないようです） （訳注： ソースコードを見ると、* 演算子は、スカラー積と内積の両方で使えるようオーバーロードされており、内積については、別途 DotProduct( ) というメソッドが Arrayクラスにだけ用意されている。）。 このクラスには、要素サイズを変更するメソッドや、通常のコンテナクラスが持っているようなメソッドはありません。理由は、これらのクラスは数学的演算の為に使われるクラスであり、コンテナクラスとして使われる事を想定していないからです。データの保管は、単純な scoped_ptr を使って行っており （訳注：上記コードでは、boost::scoped_array
Arrayクラスでは、Abs( )や Log( )といったメソッドも提供しています。C++ 開発者コミュニティーでの良き開発者でいるために、これらのメソッドを std の namespace にある同じ名前のメソッドをオーバーロードしないようにしています。なぜなら、C++ の標準で、それが禁止されているからです。より複雑な計算機能 （例えば、行列の平方根や、様々な行列の分解方法など）は、別のモジュールで提供されています。
端的にまとめると、このクラスはベクトルや行列の基本的な機能を忠実に実装しているだけです。一点だけ、意味が明解でない内容の部分は、Disposableクラステンプレートを使っている所ですが、この Appendixの別のセクションで詳しく説明します。ここでは、このクラスを使用したのは、C++11 の標準が発表される前に、move semantics と同様の効果を試しただけ、と述べておきます。すなわち、抽象化によるメモリー負荷を軽減したかったという事です。演算子のオーバーロードは大変便利なものです。 c = a + b という記述の方が、add(a,b,c) という記述よりもはるかに解りやすいですが、デメリット無しでという訳にはいきません。例えば、
Array operator+ (const Array& a, const Array& b);という関数宣言では、この演算子は新しいArrayインスタンスを生成し、それを返すようにしなければなりません。すなわち、その為にメモリーを確保し、そのメモリーをコピーする動作が必要になります。そうすると、演算の回数が増えるほど、それにかかる計算負荷は大きくなっていきます。
ライブラリーの最初のバージョンでは、その問題に対し、expression templates（式テンプレート） を使って対応しようとしていました。（以下にそれを簡単に解説しますが、詳細はVeldhuizenの“Techniques for Scientific C++”を参照して下さい。） その考え方は、演算子は Arrayのインスタンスを返すのではなく、式中の各項を示す参照を持つ parse tree のようなものを返します （訳注：parse treeとは、式の構文を解析して演算の順序を示すフローチャートに変換されたもの）。従って、例えば、2 * a + b といった算術表現は、実際にその計算をするのではなく、計算に必要な情報を持ったアルゴリズムの構造を生成するようにします。式テンプレートに配列が代入された時にだけ、算術式が解析され、その時点でコンパイラーが内容を解析し、計算ループをひとつ生成し、計算を実行して結果を代入された配列にコピーします。
このプログラム技術は今でも有効（コンパイラーの技術が進歩した今では、より重要）ですが、我々は開発開始から数年後に、この方法を止める事にしました。理由は、コードを読むのが難しく、メンテナンスが大変な事です。コードの複雑さについては、今の実装コードと、かつての下記のようなコードを比べてみて下さい。
VectorialExpression <BinaryVectorialExpression<Array::const_iterator, Array::const_iterator, Add>>
operator+(const Array& v1, const Array& v2);それと、すべてのコンパイラーがこの処理を実行できず、expression templatesとより単純な実装コードの両方をメンテナンスせざるを得なかったからです。従って、C++ の開発者コミュニティーで move semantics の話が開始され、その実装方法のアイデアが登場した際に、それをヒントにして Disposableクラスを使う事にしました。
既に述べた通り、コンパイラーの性能が近年大幅に向上し、すべてのコンパイラーが expression-template の実装に対応するようになり、その技術も向上しました。しかし、仮に今プログラムコードを設計するとしたら、問題は Array や Matrixクラスをそもそも設計する必要があったかどうかです。最低限の機能は、std;;valarray を使って自分で設計し実装できると思います。しかし、boostライブラリーが提供している uBLAS といった既存のライブラリーを使うのがいいでしょう。これらは、数値計算のコードの専門家によって書かれたものであり、すでに QuantLibライブラリーの中で、いくつか使われています。
<ライセンス表示>
QuantLibのソースコードを使う場合は、ライセンス表示とDisclaimerの表示が義務付けられているので、添付します。 ライセンス