Non-Euclidean Manifold上での近似最近傍探索(論文紹介)

こんにちは。@xiangze750です。Machine Learning Advent Calendar 2012の13日目の投稿になります。
今回はコンピュータビジョンにおける最近傍探索と幾何学についての論文紹介です。
Fast ANN Methods for Non-Euclidean Manifolds with Applications to Human Activity Analysis in Videos(pdf)という(題名どおり)人間の動作の分類を近似再近傍探索を用いて行うという論文を読んでいたのですが、リーマン幾何などの比較的高度な数学的概念が前提とされ、先行研究も多数あげられていたので内容理解のための個人的な覚え書き,疑問点の整理をかねてその内容をまとめます。

課題、問題点

近年動画からの人体の姿勢や動作の分類、推定に注目が集まっています。
推定においてはまず動画から特徴点を抽出してその軌跡をとる、オプティカルフローの分布ヒストグラムをとる(Histogram of Oriented Optical Flow )などして得られた多次元のデータに対して最近傍探索を行い似ている既知の動作を同定する方法があります。

しかし特徴点を並べたベクトルやヒストグラムは一般に高次元となるため、近傍探索は困難な物になります(次元の呪い)。
そこでそれらの点をいかに効率よく、元の分類の性質を壊さないように低次元に射影、圧縮するかが重要になってきます。
また次元を削減するのみならずベクトルの各要素のとりうる値も少なくする、可能であれば２値(バイナリ)にすることなどがなされ、Hashingと呼ばれます。
動作、姿勢を特徴付ける多次元データの分布が平らな空間ではなく歪んだ形をした高次元空間(非ユークリッド空間)に分布していると仮定をすることでデータを局所的には低次元の空間に埋め込み、短いビット列でクエリーとなる多次元データベクトルと近い点を表すことができるという主張を実証したものが今回紹介する論文です。
論文内ではまず既存の手法としてLocally sensitive Hashing(LSH), semantic hashing, Spectral hashingなどについて触れられてるのでまず(非常に簡単ですが)それらの手法を順番に解説いたします。その後本論文の提案手法である非ユークリッド空間上に分布するデータのSpectral hashingを説明します。

Locally sensitive Hashing(LSH)

LSHとは近似的な最近傍探索の手法でクエリーベクトルqに対しデータベクトルpが近いか遠いかに応じてもとのベクトルより次元の低い値を返す関数をたくさん用意して、それらの値で近傍かどうかを判定する方法です。
より詳細には
p21に対してクエリーベクトルq、データベクトルpが
$||p-q|| <= R \Longrightarrow Pr_{h \in H} [h(q) = h(p) ] >= p1$
$||p-q|| >= cR \Longrightarrow Pr_{h \in H} [h(q) = h(p) ] <= p2$
を満たすようなhash関数hの集合Hを用意しようと言うことで、近い値の入力が同じハッシュ値を持つことを確率的に保証します。
これは $h \in H$ の性質についての要請であり、その具体的構成方法としてはp-stable hashという方法が知られています。
ガウス分布やコーシー分布などの分布(p-stable分布)に従ったランダムなベクトルaへのデータベクトルxの射影をハッシュ関数h(x;a)の値とします。
$h(x;a)=int((a \dot x +b)/r)$
rの値を小さくすることでh(x;a)のとりうる値は多くなり、精度が高くなることになります。
すると２つのベクトルのx,qに対してh(x;a)-h(q;a)<1となる確率がr,c=|q-x|に応じて厳密に求まります。
ハッシュ関数の値を要素とするK次元ベクトル $g=\{h(a_1), \cdots ,h(a_K) \}$ をL個用意し、あらかじめデータベクトルに対してh(x,a_i)を計算しておいてクエリーベクトルqに対するハッシュ値 $h(q;a_1), \cdots ,h(q;a_K)$ との差とその分布をみることでx,qが近似しているかどうかを判定します。gの次元Kは検索精度、要素数Lが検索時間に影響します。
LSH,あるいはその一部であるベクトル間のコサイン類似度をハッシュ値が同一である確率とするSimHashについては既に広く応用されており、また解説文献も多数あります。日本語の説明文献としては
http://www.slideshare.net/JavaDM/lsh-pstable-1565961
が非常にわかりやすく書かれています。
本論文で参照されている論文"Near-Optimal Hashing Algorithms for Approximate Nearest Neighbor in High
Dimensions(pdf)"では当初のLSHの検索時間と正確さの間のトレードオフ関係1/cを1/c^2とするような構成方法が具体的に示されているそうです。

Semantic Hashing(pdf)

Deep learningで有名なHinton先生らの研究であり、ニューラルネットの制限された形であるRestricted Boltzmann Machine(RBM)を多段に重ねることによって学習を行うのはDeep learningと同じです。
各階層の学習をpretrainingとして行い、pretrainingで得られた階層を結合しfine tuningを行う構成もDeep learningと同様ですが、pre trainingでの目的関数(エネルギーの式)は少々異なるようです。

これによってデータ転換の距離関係を保ったまま次元圧縮を行うような非線形関数 f(x|W)を学習するという手法でこれをNeighborhood Component Analysis (NCA)と称しています*1。
具体的にはRBMで構成される各層では隠れ変数 ${h_i}$ , 観測(下位層からの入力)変数 ${x_i}$ にたいして
$p(x)=\frac{exp(-E(x,h))}{ \sum_{y,z} exp(-E(y,z))}$
$E(x,h)= -\sum_i b_i x_i - \sum_j b_j h_j - \sum_{i,j} x_i h_j W_{ij}$
$f(x|W)=W \cdot x$
を最大化するようなパラメータ{W,b}を学習データをxへ代入することで学習させます*2。
gradient ascentによるパラメータWの更新式
[tex: \Delta W_{ij} =_{data} - < x_i h_j >_{model} ]
で出てくるmodelの計算においてContrastive Divergenceと呼ばれる手法で計算を行っています。
詳しい方法は論文や解説(Restricted Boltzmann Machine の導出に至る自分用まとめ)を参考にしていただくとして概要を述べますと、隠れ層の変数h,の変数xをそれぞれの条件付き確率
p(x|h), p(h|x) を用いて交互に変更してmodelの分布関数
$p_{model}(x,h)=\frac{1}{N}\sum_k \delta(x-x_i )*p(h|x_a)$
を求めるというMCMC的な手法です。
Learning a Nonlinear Embedding by Preserving Class Neighbourhood Structure
では手書き文字データセット MNISTでLDA、PCAとの比較を2次元に射影するという形でわかりやすく既存アルゴリズムに対する優位性を示しています。
近傍数ごとのerror,Semantic Hashing(上)と既存手法(下)の結果を２次元に射影したものの比較

Spectral Hashing

Spectral Hashing
Spectral clusteringは離散グラフG={V,E}の各ノードを行、列要素を連結の有無とする隣接行列Aと対角行列Dを用いて
L=D-A
と書かれるグラフラプラシアンに対する固有値問題を解き、その固有関数の値でノードを分類する手法で一般にはNP困難であるクラスタリングの近似を行うことができます。
これと同様に離散グラフの代わりに多次元空間に分布するデータ点間の類似度
$W(x_i,x_j)=exp(|| x_i - x_j ||/d)$
を要素とする行列を用いてグラフラプラシアンと対応する行列を作り、それをベクトルの制約
$\sum_{i=1}^N y_i =1$ , $\frac{1}{N}\sum_{i=1}^N y_i y_i^T=I$ (単位行列)
のもとに対角化して得られる固有関数の値をバイナリ化するというHashing手法も考えられます。
しかし”Spectral Hashing”を表題とした論文ではこの行列の固有値問題を解くことはせずにデータ点の分布をまず超直方体内の一様分布と仮定したうえでその固有関数の符号でhashを作ることを行っています。
そしてその超直方体の辺の方向はPCAによって決めています。
一方kernel PCA同様における射影軸と同様データ分布の固有関数はN個のデータ点を用いて作られるグラム行列M_Nとその固有ベクトルv_i、固有値l_i を用いて
$f_{k,n}(x)=\frac{\sqrt{N}}{l_k}\sum_{i=0}^N v_{ik} M(x,x_i)$
(Nystorm公式)
と近似でき、元となるデータがある分布p(x)に従っている場合データ点を無限個とってきた極限ではグラム行列(ここでは類似度行列)はLaplace-Beltrami operatorと呼ばれる関数に対する作用素、 $f_{k,n}$ はその固有関数に確率収束することが知られています*3。
無限個の極限における問題全体は以下のような最適化問題として書かれます。
$minimize \int || y(x_1) -y(x_2) ||^2 W(x_1,x_2)p(x_1)p(x_2) d_x1 d_x2$
$y(x) \in \{ -1,1 \} , \int y(x) p(x) dx =0 , \int y(x)y(x)^T p(x)dx = I$
Nystorm公式の結果を用いることでグラム行列の要素を用いて固有関数を近似することができるのですが、Spectral Hashing論文ではこの近似式の計算を嫌って一様分布という大胆な近似を行っているようです。

異なるクラスタに属する離散グラフのノードを色分けして表示したものと　要素数無限大の極限での分布関数のイメージ
超直方体で(ここでは2次元)固有関数の値と固有関数の符号によるクラスタリングのイメージ

計算を高速化する別の手法としてAnchor Graphがあります。こちらはデータの中から代表となるAnchorを抽出してそれとの類似度の計算でクエリーデータに対する近傍探索を高速化する手法です(ビームさんの解説を参照してください)。
Spectral Hashingでは既存データの分類は可能ですが、クエリーベクトルを持ってきた場合にそれの再近傍データ点を探索することはできません。
データ点間の類似度をグラム行列の要素として、固有値問題の形として圧縮を行うところはkernel PCAに類似していますが、PCAではデータを分散が小さい軸にそって押しつぶして次元を削減するのに対してSpectral Hashingでは複数の固有関数の値を使ってデータを色分けするような感覚であると言えるかもしれません。
Spectral Hashing論文では人工的なデータ、LabelMeのデータに対して比較評価を行っておりその結果も示されています。
Laplace-Beltrami operatorといういかにも難しそうな数学的概念が出てきます。これは平坦でない空間(リーマン多様体)上でのラプラシアンの一般化なのですが、手法を実行する際にはデータが何らかの多様体上を分布しているということを意識することはあまりありませんし、実質的にはラプラシアンに相当するものといってよいと思います。

一方Diﬀusion Maps, Spectral Clustering and Reaction Coordinates of Dynamical Systemsでは離散グラフ上のマルコフ過程の遷移行列としてグラフラプラシアン、類似度行列をとらえており、対応する確率分布の発展方程式(Fokker-Plank方程式)、作用素との関係、また遷移行列の作り方を変えることで発展方程式を平坦な空間+ポテンシャル上のものか、歪んだ空間上のものかにmappingできる(4章)ということについて解説されています*4。

Spectral Hashing on Non-Euclidean Manifolds

(やっと)本題の論文です。上記spectral hashingを一般的な非ユークリッドな多様体上に拡張し、任意のクエリーデータに対しての近傍探索が行えるようなものとして

Riemannian Spectral Hashing
Distributed Kernel Spectral Hashing

の２つの手法のvariationが紹介されていますが、上記Spectral Hashingとの違いは表題の論文のFig. 1に集約されています。

Riemannian Spectral Hashing(RSH)

Spectral HashingではFig. 1.の(a)にあるようにデータ全体を(超)立方体に写し、その中での位置でHash値を決めていました。
Riemannian Spectral Hashingでは(b)のようにSpectral Hashingは部分ごとに行います。
データ全体がある多様体Mを形作るとして、ある点yに多様体の次元と同じ次元のユークリッド空間を対応づけたものを接空間 $T_y M$ と呼びます。
この多様体M上(のyに近い点の集合)から $T_y M$ への写像をlogarithm map (log-map)、その逆関数をexponential mapと呼んでいます。
そして接空間に移されたデータ点に対してSpectral Hashingを行います。
具体的手順としては以下のようになります。

学習フェーズ
- Riemannian k-meansを用いてデータをクラスターに分割する。
- クラスターごとに所属するデータ点のlog-mapを計算する。
- クラスターごとにSpectral hashingを行う。
検索フェーズ
- 各クラスターの重心のうち最もクエリーに近い物を探索し、そのクラスターをクエリーの所属クラスターとする。
- クエリーを探索の結果得られたクラスターの接空間にlog-mapで写す。
- Spectral hashingでクエリーに対応するcodeを計算する。
- (近似)最近傍探索を行う。

上記曲がった空間でのk-meansであるRiemannian k-meansは以下のような手順です。

データ点集合からクラスターの中心を選び出し、{c_j} とする(Riemannian k-means)。
各データx_iに対し d(c_j,x_i) =||log_{c_j} (x_i) ||を計算し、その値が最小となるようなlabel jをx_iの属するクラスターとする
各クラスターjに対しc_j=mean{x_l| w_l =j} を再計算する
繰り返し

データ全体からの代表的な点を選ぶところはAnchor graph hashingに共通していますが、Anchor graph hashingではハッシュ値の計算に必要な点を減らしますが、Spectral hashingの適用範囲はデータ全体になります。RSHではSpectral hashingは代表点(クラスタの重心)ごとに異なります。

Distributed Kernel Spectral Hashing (DKSH)

log-mapの式が明示的に書けない場合にはデータ点間の類似度行列をカーネルとしたMDS(多次元尺度法)を用いてユークリッド空間に埋め込んだ後にk-meansを行い、重心に最も近いデータ点をクラスターの代表点(pivot)とする手法(Non-Linear clusteringと呼んでいます)。
以後はRSHとほぼ同じです。

学習フェーズ
- Non-Linear clusteringでクラスタリング
- pivotを選ぶ
- クラスターごとにSpectral hashingを行う。
検索フェーズ
- 各クラスターのpivotのうちクエリーと類似度が最も高い物を選ぶ
- Spectral hashingでクエリーに対応するcodeを計算する。
- (近似)最近傍探索を行う。

実験と評価

まず人工データとして

超球面 $S^{99}$
Grassmann manifold - G10;3(10次元ユークリッド空間の3次元部分空間すべてをあわせたもの)

という形状(多様体)からデータ点をとってきたものを使っています。
図形の定義式があることからlogarithtic mapは閉じた形で求めることができ、超球面 $S^{99}$ の場合にはarccos関数を用いて
$log_x(y) =\frac{y-(x^T y))x}{||y-(x^T y))x||}arccos(x\T y)$
$exp_x(\Delta)= cos(||\Delta||)x + sin(||\Delta||)\frac{\Delta}{||\Delta||}$
と書けます。
Grassmann manifoldなどは純粋数学以外ではほとんどお目にかからないと思っていたのですが、
拘束のある人体の一部分の運動はこの多様体に沿っていると考えられるそうです(参考スライド(pdf) page 46以降 )。
リアルな人間の動作データでの性能評価としてはKTH human action datasetを用いています。このデータセットから時空間３次元のなかのキー領域として抽出したものとHistogram of Oriented Optical Flow(HOOF)*5に対して
Nearest Neighbor法(NN)をbaselineとしてKernel LSH (KLSH),Spectral hashing (SH), Kernel Spectral hashing (KSH)などの他の近似最近傍探索法とRSH,DKSHの比較を行っています。

HOOFのヒストグラム間の距離としてBhattacharyya distanceやEarth Mover’s Distanceなどが考えられ、それをカーネルk(x,y)としてみることができます。すると
$k(x,y)=\Phi(x)^T \Phi(y)$
と書かれる高次元の特徴空間への写像Φで状態空間モデルが
$x_{t+1}=Ax_t +B v_t$
$\Phi (y_{t})=Cx_t + w_t$
(v,wはノイズ)
と書けます。
行列 $A \in R^{n \times n} ,C \in R^{p \times n}$ を用いて時系列を記述した大きな行列
$O = [C^T , (AC)^T, (A^2C)^T, \cdots , (A^{n-1}C)^T ]^T \in R^{p*n*n}$
の列ベクトルの部分が $R^{p \times n \times n}$ の $R^{n \times p}$ 次元部分空間であることからグラスマン多様体上にあることになり、RSH,DKSHを使うことができることがわかります。*6

正確さ、学習にかかる時間、検索にかかる時間が Table 5,Table 6に示されています。*7
RSHは２つの結果の双方で近似最近傍探索の中では最良の一致率を示しています。Table 6 observability matricesに関しては一致率はKLSHとさほど変わりませんが、検索にかかる時間が非常に短いのが利点となっています。その分学習に要する時間は多いのですが、10回程度検索すれば元がとれることになります。
クラスタ(接空間)ごとにhashingを行えばよいためかコード長が短いことも利点としてあげられています。しかしその分学習には時間がかかってしまうのだと考えられます。一方DKSHの方は残念なことに他の手法と比べて一致率がよいという訳ではありません。

感想、個人的課題

test data setがあげられているのでそれを使った追試、比較評価を行いたいです。
データの読み込みや疎行列の取り扱い、同一条件での比較がスキル的に難しそうです。
オープンソースであり、手書き文字データMNISTの扱い、多様体学習の比較もグラフを使ってわかりやすく紹介されているscikit.learnを参考にnumpy,scipyを用いるのが良さそうに思えます。