最近傍相関行列

nAGライブラリには最近傍相関行列の計算に関連した様々な機能があります。この記事では、最近傍相関行列の問題を見ていき、背景を説明し、それらを解くルーチンを紹介します。

イントロダクション

相関行列は、次のような実数の正方行列であることを特徴としています。

• 対称
• 対角線上がすべて１
• 非負の固有値を持つ

非負の固有値を持つ行列は正半定値と呼ばれます。 行列Cが相関行列である場合、C,c_ijの要素は、エンティティiとエンティティjのペアワイズ相関、すなわち、2つの間の線形関係の強さと方向を表します。

文献（参考文献）には相関行列の使用例が数多くありますが、私たちが最もよく遭遇したのは、金融の分野で、様々な銘柄間の相関関係が賢明なポートフォリオを構築するために使用されている場合です。残念ながら、様々な理由により、相関行列であるはずの入力行列が半定値にならないことがあります。例えば、相関関係は、ある期間に測定された銘柄間のものであり、データが欠落している場合があります。個々の相関が、それらが共通に持っている観測データを使用して計算され、これがすべての変数にわたって変化する場合、それは不定行列を生じさせます。さらに金融からの引用ですが、実務家は、過去の値から計算されたものとは異なる特定の資産間の相関を割り当てることによるポートフォリオへの影響を調査したいと思うかもしれません。これにより、計算された行列に負の固有値が生じる可能性もあります。

このような状況下では、結果は近似相関行列となります。しかし、真の相関行列であることに依存する後続の分析のために、これを修正する必要があります。 理想的には、「近い」の定義に基づき、近似相関行列に「最も近い」真の相関行列を見つけたいと考えます。これがここで言う最近傍相関行列問題です。

基本的な最近傍相関行列問題

nAG ライブラリルーチンnagf_correg_corrmat_nearest (g02aa) はイントロダクションで概説した基本的な問題を解くためのニュートンアルゴリズムを実装しています。このルーチンはフロベニウスノルムの近似入力行列Gに最も近い真の相関行列Xを求めます。つまり以下の最小化を行います。

$$\|G-X\|_F .$$

このアルゴリズムは、Qi と Sun の論文 [8] に記載されており、以前に提案されたアプロー チよりも優れた収束特性を持っています。マンチェスター大学の Borsdorf と Higham [2] は、これをより詳細に調べ、さらなる改善を提案しました。それには、異なる反復ソルバー（共役勾配よりもMINRESが好まれました）と、線形方程式の事前条件付けの手段が含まれています。この改良されたアルゴリズムが、我々のライブラリに組み込まれています。

重み付き問題と正定値相関行列の強制

nAG ライブラリルーチン nagf_correg_correg_cormat_nearest_bounded (g02ab) は g02ab で提供されている機能を拡張しています。近似相関行列がある場合、行列の全てが実際には近似ではなく、一部だけが近似であると仮定するのが妥当です。同様に、いくつかの相関を他のものよりも信頼し、最終的な行列の入力値に近い値になるようにしたいと考えるかもしれません。

このアルゴリズムでは、Qi と Sun のオリジナルの研究を応用し、重み付きノルムを使用します。つまり以下の最小化を行います。

$$\|W^{1/2}(G-X)W^{1/2}\|_F.$$

ここで Wは重みの対角行列です。これは、我々が要素 $\surd w_{ij}(g_{ij}-x_{ij})\surd w_{jj}$ を最小化しようとしていることを意味します。このように、Wの要素を適切に選択することによって、我々は、Gのいくつかの要素を有利にし、Xの対応する要素をそれらに近づけることを強制できます。

この方法では、Gの行と列全体が重み付けされることになります。ただし、nagf_correg_corrmat_h_weight (g02aj) では要素ごとの重み付けが可能なので、このルーチンでは以下の最小を求めます。

$$\|H \circ (G-X)\|_F,$$

ここで $𝐶=𝐴\circ 𝐵$ は以下の要素を持つ行列Cを示します。 $c_{i,i}=a_{i,j}{b}_{i,j}$

従って、Hに適切な値を選択することで、Gに含まれる個々の要素を強調し、他の要素は加重されずに残すことができます。ここで採用されているアルゴリズムは、Jiang, Sun and Toh [7]によるもので、ニュートンアルゴリズムを核にしています。

g02ab と g02aj の両方で、計算された相関行列が正定値であること、つまり固有値がゼロより大きいことを指定することができます。これは、行列の条件を改善し、安定性を高めるために、いくつかのアプリケーションで必要とされます。

相関行列のランクの制約

低ランクの相関行列が必要な場合、例えば独立ランダム変数の数を制約するためにnagf_correg_corrmat_nearest_rank (g02ak) を使用することができます。これは、フロベニウスノルムの中から、指定されたランクの最大値を持つ最寄りの相関行列を求めます。このルーチンは、GaoとSun[4]によって提案されたMajorized Penalty Approachに基づいています。

縮小投影と交互投影による相関関係の修正

ここでは、真の相関であることが知られている要素のいくつかを修正することに注意を向けます。前の4つのアルゴリズムのようにニュートン法を使う代わりに、ここでは縮小法を使います。

これが必要とされる一般的な例としては、変数のサブセット間の相関が信頼されており、それ自体が有効な相関行列を形成している場合があります。このような場合、これらを入力行列の先頭ブロックに配置し、その間に、残りの部分を修正することができます。これを固定ブロック問題と呼びます。 ルーチンnagf_correg_corrmat_shrinking (g02an) は近似行列のこのような先頭ブロックの相関を保持します。Higham、Strabic、及び Sego [6]の縮小法を用いて以下に示される相関行列を見つけます。

$$\mathrm{\alpha <!--\; \alpha \; -->}\left(\begin{array}{cc}A& 0\\ 0& I\end{array}\right)+(1-<!--\; -->\mathrm{\alpha <!--\; -->})G.$$

ここでもG は入力行列であり、区間 [0,1] で正の半定値結果を与える最小のαを求めます。α が小さければ小さいほど、元の行列に近づき、どのようなαであっても、正定値である必要がある先行部分行列 A を保持します。このアルゴリズムは，有限のステップ数で素早く収束する二等分法を用いています。

ルーチン nagf_correg_correg_corrmat_target (g02ap) は、縮小の考え方を一般化し、我々自身のターゲット行列を供給することを可能にします。ターゲット行列Tは $𝑇=𝐻\circ 𝐺$ とする重みの行列 H を指定することで定義されます。そして、以下の解を求めます。
$$\mathrm{\alpha <!--\; -->}T+(1-<!--\; -->\mathrm{\alpha <!--\; -->})G,$$

最小固有値の境界を指定することも可能です。Hに1の値を指定することによりGの要素が固定され、従ってXでは変更されません。

例えば，2つの対角線ブロックを固定する必要がある場合があります。そのような場合以下のようなHを選択できます。

$$H=\left(\begin{array}{cc}\left[\begin{array}{ccccc}1& & \dots <!--\; -->& & 1\\ ⋮<!--\; -->& & \ddots <!--\; -->& & ⋮<!--\; -->\\ 1& & \dots <!--\; -->& & 1\end{array}\right]& 0\\ 0& \left[\begin{array}{ccccc}1& & \dots <!--\; -->& & 1\\ ⋮<!--\; -->& & \ddots <!--\; -->& & ⋮<!--\; -->\\ 1& & \dots <!--\; -->& & 1\end{array}\right]\end{array}\right).$$

次に、アルゴリズムは、次の形式の正の半定値行列を与える最小のαを見つけます。

$$\mathrm{\alpha <!--\; -->}\left(\begin{array}{cc}{G}_{11}& 0\\ 0& G_{22}\end{array}\right)+(1-<!--\; -->\mathrm{\alpha <!--\; -->})G.$$

そして、入力の2つの対角線外のブロックのみを摂動させます。

縮小アルゴリズムは、その速度と入力と出力の間の距離が大きくなる可能性があることが特徴です。アンダーソン加速度を用いた交互投影は、固定ブロック問題を計算するために採用している別のアルゴリズムであり、その速度と近さの特徴は、縮小の場合とは逆です。

入力行列は、単位対角線を含む、保存したいエントリを持つ半定値行列と行列の集合の中で、最も近い行列に繰り返し、交互に投影されます。理論的には収束を保証するものではありませんが、実際には、このアルゴリズムは、これら2つの集合の交点に最も近い相関行列を見つけます。

ルーチン nagf_correg_corrmat_fixed (g02as) は、Higham and Strabic [5] の方法を採用しています。これは、任意の要素を固定し、オプションで最小固有値を設定しながら、フロベニウスノルムの最近傍相関行列を計算します。

最近傍相関行列ルーチンの選択

ルーチンを選択するとき、計算時間と元の行列から解の行列への距離のトレードオフがあります。ニュートンアルゴリズム（g02aa、g02ab、g02aj、およびg02ak）と交互射影アルゴリズム（g02as）は、重み付きアルゴリズムが入力にのみ影響を与えることを考慮しつつ、常に解いている問題に最も近い解を見つけます。縮小ルーチン(g02anとg02ap)は、より遠くの結果を見つけますが、より高速です。

基本的な問題では、g02aa は常に最も近い行列を見つけます。g02apを使用して、ターゲットとして同一性行列を指定すると、これよりも遠くの行列を求めることになりますが、これは解の形式を考えると理解できます。

固定ブロック問題を解きたい場合は、専門家向けのルーチン g02an が最も高速です。ニュートンアルゴリズムの中では、g02abは妥当な時間で解を見つけますが、行と列全体に重み付けを行うため、正しいブロックの外側では、いくつかの要素が強調されすぎてしまいます。g02ajを使えば、より正確な重み付けが可能で、近い解が得られます。しかし、このルーチンはかなり時間がかかります。最も近い解決策はg02asで見つかるでしょう。ニュートンアルゴリズムよりもはるかに遅いですが、速度の点ではg02ajよりも優れています。

2つの対角線ブロックを固定したり、任意の固定と重み付けを行う場合は、同じ速度と近さのトレードオフでg02aj、g02ap、g02asの間で選択します。g02asの交互投影は、要素を固定し、最も近い解を見つけます。縮小アルゴリズムは要素を固定し、高速ですが、ターゲット行列は正定値であり、有効な相関行列の一部を形成する必要があり、これは制限となります。g02aj は入力中の要素を重み付けするだけなので、保存したいブロックが正半定値に近いが、正半定値にならない場合、ここでは柔軟性があるかもしれません。

固有値の最小値を固定したい場合で、重み付けが不要な場合は、g02ab または g02ap を使用することができます。後者では、基本問題と同様に同一性目標を使用します。重み付けや固定も必要な場合は、上述の問題についても同様の結果が得られます。しかし、重み付けとの組み合わせでは、g02ap は大きなαの値を返すことがあります。 これは、入力行列の多くが失われ、それとはかけ離れた結果が返されることを意味します。

出力相関行列のランクを制約するには、g02akを使用します。

収束を定義するすべてのアルゴリズムで使用されている許容値は、明らかに実行される反復の数に影響を与え、その結果、速度と近さに影響を与えます。典型的な問題を表すデータを使用して、いくつかの実験を行うことをお勧めします。ルーチン g02aj は使用される重みに敏感になることがあるので、異なる値を試して近さと計算時間の両方を調整する必要があります。

因子構造を持つ最近傍相関行列

因子構造を持つ相関行列は、非対角要素がランクkの行列と一致する行列です。 相関行列Cは次のように書くことができます。

$$C=\text{diag}(I-<!--\; -->X{X}^T)+X{X}^T,$$

ここで、Xはn×kの行列であり、因子負荷行列と呼ばれることが多いものです。また、kは一般的にnよりもはるかに小さい値となります。 このような相関行列は、資産収益、債務担保証券、および多変量時系列の因子モデルで良く出現します。

ルーチンg02aeは最近傍因子負荷行列Xを計算します。これは、以下の最小値を見つけることにより、近似値Gに最も近い相関行列を与えます。

$$\Vert <!--\; -->G-<!--\; -->X{X}^T+\text{diag}(X{X}^T-<!--\; -->I){\Vert <!--\; -->}_F.$$

我々はBorsdorf, Higham and Raydan [3]が提案したBirgin, Martinez and Raydan [1]のスペクトル投影勾配法を実装しました。

機能表

以下の表は、最近傍相関行列ルーチンをすべて示し、近さの尺度と、それぞれで使用できる重み付けと固定方法を示しています。

ルーチン	フロベニウスノルムによる近さ	縮小アルゴリズム	因子構造を持つ最近傍行列	要素の重みづけ	要素の固定	最小固有値のリクエスト	最大ランクのリクエスト
g02aa	X
g02ab	X			X		X
g02ae	X		X
g02aj	X			X		X
g02ak	X						X
g02an		X			X
g02ap		X		X	X	X
g02as	X				X	X

参考文献

[1]   Birgin E G, Martinez J M and Raydan M $2001$ Algorithm 813: SPG-software for convex-constrained optimization ACM Trans. Math Software 27 340-349

[2]   Borsdorf R and Higham N J $2010$ A preconditioned $Newton$ algorithm for the nearest correlation matrix IMA Journal of Numerical Analysis 30$1$ 94-107

[3]   Borsdorf R, Higham N J and Raydan M $2010$ Computing a nearest correlation matrix with factor structure. SIAM J. Matrix Anal. Appl. 31$5$ 2603-2622

[4]   Gao Y and Sun D $2010$ A majorized penalty approach for calibrating rank constrained correlation matrix problems Technical report Department of Mathematics, National University of Singapore

[5]   Higham N J and Strabic N $2016$ Anderson acceleration of the alternating projections method for computing the nearest correlation matrix Numer. Algor. 72 1021-1042

[6]   Higham N J, Strabic N and Sego V $2014$ Restoring definiteness via shrinking, with an application to correlation matrices with a fixed block MIMS EPrint 2014.54 Manchester Institute for Mathematical Sciences, The University of Manchester, UK

[7]   Jiang K, Sun D and Toh K-C $2012$ An inexact accelerated proximal gradient method for large scale linearly constrained convex SDP SIAM J. Optim. 22$3$ 1042-1064

[8]   Qi H and Sun D $2006$ A quadratically convergent Newton method for computing the nearest correlation matrix SIAM J. Matrix Anal. Appl. 29$2$ 360-385