Online PCA - ややプログラム紀行

最近僕はcontinual learningの文脈でタスクの数に依存せずにリプレイバッファを構築する方法について考えているんだけど、ここら辺を研究しているとなんとなくonline pca、つまりオンライン主成分分析と関連性があるような気がしてくる

リプレイバッファというのは言わば訓練データの要約のようなもので、それをタスクの数に依存しない定数個に抑えるためには常にデータの取捨選択をする必要があり、タスク単位での分布に対して何かしらの仮定を置かないと効率的なアルゴリズムは構築できないような気がしてくる

しかし、この状況設定はまさにオンライン最適化、特にオンライン主成分分析と同じなので、オンライン主成分分析についての知見はあるに越したことないんじゃ...と思って少し調べてみた

オンライン主成分分析の手法はいくつか存在していて、サーベイ論文とかがあまり見つからないから体系的な理解がまだ出来ていないんだけど、手始めに読んだ論文がなかなか面白いアルゴリズムだった

(わかりにくいけど上の画像はpdfへの直リンク)

主成分分析

そもそも主成分分析とはn次元のデータを複数個受け取った時、それらのデータを最もよく表現できるような $k\ (<< n)$ 次元部分空間を見つけ出す*1ことで、オンライン主成分分析では複数個のデータがオンラインな形で、すなわちデータを1つずつ観測し、その度に対応する部分空間を更新していく

数式的にいうなら、 $t-1$ ステップまでに観測したデータをもとに $k$ 次元部分空間への射影行列 $P^{t-1}$ と平均値 $m^{t-1}$ を構築したとし、 $t$ ステップ目に新しいデータ $x^t$ を受け取るとき、"圧縮損失"と呼ばれる損失 $\|(x^t - m^{t-1}) - P^{t-1}(x^t - m^{t-1})\|_2^2$ を被る

以上のような段取りが各ステップごとに行われ、 $T$ ステップまでの圧縮損失の合計値を最小化するような射影行列と平均値を計算するアルゴリズムを考えることになる

ヘッジアルゴリズム

上の論文で提案されている論文は、オンライン最適化におけるエキスパート統合問題およびヘッジアルゴリズムと密接に関連しているので、まずはその説明をする

エキスパート統合問題においてプレイヤーは $n$ 人のエキスパートの答えを統合して次々出題される問題に答える必要があり*2、ナイーブなやり方としては各エキスパートの今までの問題に対する正答率を記録しておき、その正答率が最も高いエキスパートに従うなどが考えられる

ヘッジアルゴリズムはこの考えをもう少し洗練したもので、今までの問題に対するエキスパートの挙動に基づいて計算した確率ベクトル $w^{t-1}$ から確率的にエキスパートを選び出す

より具体的には、 $t$ ステップ目において
1. エキスパート $i$ を確率 $w_i^{t-1}$ で選び出し
2. それに対応する損失 $\ell^t_i$ を被り*3
3. 確率ベクトルを各 $i$ について $w^t_i = \frac{w_i^{t-1}\exp(-\eta \ell^t_i)}{\sum_{j=1}^n w_j^{t-1}\exp(-\eta \ell^t_j)}$ と更新する

確率ベクトルの更新式の分母はベクトルが確率になるように正規化を行なっており、損失を被らない答えを導き出しているエキスパートほど選ばれやすいように確率ベクトルが更新される

Randomized Online PCA Algorithms

もう一度オンライン主成分分析の圧縮損失を見返してみると、今までのデータの平均値を $\overline{x}$ とすれば射影行列 $P$ の圧縮損失は

$\displaystyle \begin{align*} \sum_{t=1}^T \|(x^t - \overline{x}) - P(x^t - \overline{x})\|_2^2 &= \sum_{t=1}^T (x^t - \overline{x})(I - P)^2 (x^t - \overline{x}) \\ &= \mathrm{tr} \left( (I-P)^2\sum_{t=1}^T (x^t - \overline{x})(x^t - \overline{x})^\top \right) \\ &= \mathrm{tr} \left( (I-P)\sum_{t=1}^T (x^t - \overline{x})(x^t - \overline{x})^\top \right) \end{align*}$

とかける( $I-P$ も射影行列なので $(I-P)^2 = I-P$ となる)

先程のオンライン最適化の文脈に合わせると、各部分空間 $I-P$ をエキスパートとした時に、上式の右辺をなるべく小さくするようなヘッジアルゴリズムを考えることができるが、当然そのままでは部分空間が無限個存在するため計算を行うことができず、この論文の賢い点はまさにこれを可能にした点にある

証明の詳細は省くが、 $n-k$ 次元への射影行列*4を $\frac{1}{n-k}$ 倍したものを考えるとこれは密度行列になっており、そのような無限個の密度行列の凸包を $\mathcal{A}^n_{n-k}$ とおくと、最大固有値が $\frac{1}{n-k}$ 以下である密度行列全ての集合 $\mathcal{B}^n_{n-k}$ と一致することが示され、特に $\mathcal{B}^n_{n-k}$ の要素は高々 $n$ 個の $\mathcal{A}^n_{n-k}$ の元の線型結合で表せることも示すことができ、この線型結合の係数は確率ベクトルになっているため、 $\mathcal{B}^n_{n-k}$ の要素を確率ベクトルと同一視することでヘッジアルゴリズムを適用することができる

具体的に手続きを述べると、初めに適当な $W^0 \in \mathcal{B}^n_{n-k}$ をとり、 $t-1$ ステップまでで $W^{t-1} \in \mathcal{B}^n_{n-k}$ が得られているとすると、
1. $W^{t-1} = \mathcal{W}\omega \mathcal{W}^\top$ と対角化
2. $(n-k)$ -cornerと呼ばれる、 $n-k$ 個の要素が $\frac{1}{n-k}$ で他が $0$ であるベクトル $r_j$ を用いて $\omega$ の対角ベクトルを $\sum_{j} p_j r_j$ と分解することができ、 $p_j$ は確率ベクトルになっている*5
3. 確率 $p_j$ でベクトル $r=r_j$ を選ぶ
4. $R = \mathcal{W}\mathrm{diag}(r)\mathcal{W}^\top \in \mathcal{A}^n_{n-k}$
5. $P^{t-1} = I - (n-k)R$ によって射影行列を得、損失 $\mathrm{tr}( (I-P^{t-1})x^t(x^t)^\top$ を被る*6
6. $\hat{W}^t = \frac{\exp(\log W^{t-1} - \eta x^t(x^t)^\top)}{\mathrm{tr} (\exp(\log W^{t-1} - \eta x^t(x^t)^\top) )}$
7. cappingにより $W^t = \mathrm{cap}_{n-k}(\hat{W}^t)$ を得る
という感じになる