『Kernel Random Projection Depth for Outlier Detection』をPythonで実装した話

機械学習 python

はじめに

最近，下記の論文がアップロードされていた．

『Kernel Random Projection Depth for Outlier Detection』

arxiv.org

Pythonによる簡単なデモンストレーションを実装したので，公開する．

以降，Kernel Random Projection DepthをKRPDと呼ぶ．

KRPDの中身

ざっくり言えば，Kernel PCAによる非線形な変換をかけた後に，random projection depthを計算するアルゴリズムである． そのdepth値に-1を掛けた値を外れ値スコアとして利用する．

実験および実装

データセット

いくつかのtoy dataset上で，外れ値スコアの等高線を表示するデモンストレーションを用意した． 具体的には 以下の4つである．

  • Blobs: 3つのガウス分布から発生させたデータセット
  • Cross: 十字の形に分布するデータセット
  • Circles: scikit-learnに付属のmake_circles利用
  • Moons: scikit-learnに付属のmake_moons利用

それぞれのデータセットのサンプル点は総数が400個，うちinlierが300個，outlierが100個である．

実装

KRPDの実装およびデモンストレーションは以下のNotebookである．要PyODライブラリ．

その他の代表的な検知器（LOFやOCSVM，Isolation Forestなど; 下記参照）によってどのような等高線が得られるかを比較するNotebookたちも用意した． 検知器の実装はPyODライブラリによった．

ちなみに検知器たちは以下の8つである．

  • Angle-based Outlier Detector (ABOD)
  • Histogram-base Outlier Detection (HBOS)
  • $k$-Nearest Neighbor (KNN)
  • Local Outlier Factor (LOF)
  • Isolation Forest (IForest)
  • Cluster-based LOF (CBLOF)
  • Minimum Covariance Determinant (MCD)
  • One-Class Support Vector Machine (OCSVM)

実験結果

KRPDの等高線プロットを図1に示す．オレンジ色で示された領域は外れ値の比率（25%）から決められた，外れ値スコア上位の領域を示す．いずれのデータセットにおいても，この領域はinlierの分布をうまく囲うことができている．

図1：KRPDの等高線プロット

以下，図2から図5まで比較のための等高線プロットをそれぞれ示す．

図2：Blobs上での各検知器による等高線プロット

図3：Cross上での各検知器による等高線プロット

図4：Circles上での各検知器による等高線プロット

図5：Moons上での各検知器による等高線プロット

おわりに

Python （具体的にはPyOD）でKRPDを実装した．諸課題は論文に書かれている．

おまけ

InlierのみでKRPDをfitするColab notebookも作った．

Deep Divergence Learning (ICML 2020) の論文に掲載された実験結果を検証する試み −分布クラスタリング 前編−

深層学習 機械学習 距離学習

はじめに

距離学習および深層距離学習について述べた記事の最後に，

では距離関数自体を学習可能にする枠組みは存在するのか？→Yes.

とあった．その枠組みに関する論文が，"Deep Divergence Learning" と題して ICML2020にて発表された．

Cilingir, H.K., Manzelli, R. and Kulis, B., "Deep Divergence Learning," Proceedings of the 37th International Conference on Machine Learning, 2020.

proceedings.mlr.press

著者による実験のソースコードが長らく公開されておらず，空のリポジトリのみが存在していたが，2022年の2月にjupyter notebookが公開された．

github.com

残念ながら，このnotebookは論文に掲載された実験を再現するには余りにも不十分だった．それならば，最初から全て PyTorch & scikit-learn で書いて実験を再現してやろうと決めた（再現できるとは言ってない）．本記事はそれら実験結果のうち特に分布クラスタリングに関して再現実験を試みたときの報告である． 当初は分布クラスタリングで単独の記事のつもりとしていたが，実験単位で記事を2つに分けることにした．本記事はそのうちの最初の記事である．

論文の概要

本記事に関連する範囲で論文の概要（という名のメモ/抜粋そして私的な解釈）を述べる．

図1に提案手法と既存の距離学習手法との関係図を示す．図は上記 "Deep Divergence Learning"から引用する． 従来の距離学習における距離関数はユークリッド距離に留まり，また特徴量変換は線形なものに留まっていた（左上）．特徴量変換をニューラルネットによる非線形なものに拡張すると deep metric learning となる（左下）． 一方，距離関数自体をBregman divergenceに（学習可能な形で）拡張すると Bregman divergence learning となる（右上）． 提案手法は素のBregman divergenceではなく，functional Bregman divergenceを活用したdeep divergence learningとなる（右下）．

図1: 提案手法（deep divergence learning）と既存の距離学習手法との関係

Bregman divergenceは2点（ベクトル）の距離を測る道具であり，凸関数$\phi$で特徴づけられる．定義を書いておこう． $\phi$をある凸集合$\Omega$上で定義された連続的微分可能な狭義凸関数とするとき，2点$\mathbf{x}, \mathbf{y} \in \Omega$に関するBregman divergenceは次式で与えられる．

$$ D_\phi(\mathbf{x}, \mathbf{y}) := \phi(\mathbf{x})-\phi(\mathbf{y})-\langle \mathbf{x} - \mathbf{y}, \nabla \phi (\mathbf{y}) \rangle. $$

機械学習や信号処理の分野でもよく出会うユークリッド距離，KL-divergence，板倉−斎藤距離などは凸関数のとり方をそれぞれ変えて得られる*1 ．

functional Bregman divergenceは関数ペアの間のdivergenceを測る道具である．こちらも定義を書いておこう． 2つの関数ペア$p$および$q$，そして狭義凸汎関数$\phi$を考える．凸汎関数は関数の凸集合で定義されており，実数に値を取るものとする．functional Bregman divergenceは次式で定義される．

$$ D_\phi(p, q) := \phi(p) - \phi(q) - \int [p(x) - q(x)] \delta \phi (q) (x) dx $$

ここで$\delta \phi (q)$は$\phi$の汎関数微分（Fréchet 微分）である．本論文では特に分布関数のペアを考える．

論文ではfunctional Bregman divergenceが対称となる必要十分条件を与えている．ここでfunctional Bregman divergenceが対称とは，（その定義域に入る）全ての関数$p, q$に対して$D_{\phi}(p, q) = D_{\phi}(q, p)$が成り立つものとして定義する．

定理　functional Bregman divergence $D_{\phi}(p, q) $が対称となるためには，対称な半正定値関数$\psi (x, y)$を用いて以下のように書けることが必要十分である． $$ D_{\phi}(p, q) = \iint (p(x) - q(x))(p(y) - q(y)) \psi(x, y) dx dy $$

そこから従来の（深層）距離学習に現れる距離関数たちが特殊な具体例として導出される．Linear Metric Learning，Kernel Metric Learningの距離関数，そして深層距離学習との関連でmoment-matching 型の距離関数が挙げられている．

$$ D_{\phi}(p, q) = \left\lVert \mathbb{E}_{p} [f_{W}] - \mathbb{E}_{q} [f_{W}] \right\rVert ^{2} $$

ただし$p$と$q$はそれぞれ分布関数，$\phi$は凸汎関数である．$f_{W}$は$W$をパラメタに持つニューラルネットによる埋め込みを表す．これは上記定理で$\psi (x, y) = f_{W}(x) ^{T} f_{W}(y)$と置いたケースである．特に$p$と$q$がそれぞれ観測データ集合$P$，$Q$に関する経験分布だとすると，次式で書ける．

$$ D_{\phi}(p, q) = \left\lVert \frac{1}{|P|} \sum_{x \in P} f_{W} (x) - \frac{1}{|Q|} \sum_{y \in Q} f_{W} (y) \right\rVert ^{2} $$

ただし$|P|$はその集合$P$の濃度を表す． $|P|=|Q|=1$の場合，経験分布はデルタ関数となり，従来のユークリッド距離に基づく深層距離学習に帰着する．上式は一次モーメントの意味でマッチする分布の距離を測っている．

さらに論文は学習可能な凸汎関数を導入したdeep Bregman divergenceの定式化へと繋がっていくが，本記事は上記moment-matchingまでの範囲で実験を行う． deep Bregman divergenceの実装と再現実験は次回以降の記事に回す．というのも，deep Bregman divergenceはクラスタリングの実験には使われていないからである．

応用可能性 −分布クラスタリング−

通常のクラスタリングはデータ点の集合に対して行われるが，ここでは経験分布（関数）の集合に対するクラスタリングを考える．経験分布はまたデータ点集合を伴う．例えばオンラインショップにおける各商品の評価値（レーティング）の分布が挙げられる; Amazonの各商品はユーザ投票による星1つ〜星5つの5段階評価値の集合と，そのヒストグラムすなわち経験分布を有している．2つの経験分布の距離を測る道具を我々は既に手にしているわけなので，その距離に基づくクラスタリングアルゴリズム，例えば$\mathit{k}$-meansクラスタリングを駆動できるだろう，と考えるのは自然な発想である．

Davis and Dhillon (2006) *2 は教師なし分布クラスタリングアルゴリズムを提案している．ただし各分布が多次元ガウス分布という仮定と，分布の平均ベクトル群が各クラスタごとに線形分離可能という暗黙の仮定を置いていると著者らは指摘する．後者の仮定は，分布群の配置があまり複雑に入り組んでおらず，「クラスタリングしやすい」状況を指摘している．その一方でdeep divergence learningの本論文において著者らは，分布クラスタリングに関して「半教師あり」（semi-supervised）な枠組みを考察している．同じクラスタに属するべき分布ペア，もしくは別のクラスタに分けるべき分布ペアが訓練データ上で与えられるという想定のもとに，contrastive loss あるいは triplet loss に基づく深層距離学習の適用を見据えている※． たとえ元の分布群の配置が線形分離不可能であっても，埋め込みを実現するニューラルネット$f_{W}$によって，埋め込み空間上で線形分離可能な配置へと分布群を移動できると期待する．あとはfunctional Bregman divergenceに基づくクラスタリングアルゴリズムを駆動すればよい．加えて各分布が多次元ガウスである仮定を置く必要もない．"Davis and Dhillon" の手法とは，これらの点で区別される．

※ 著者らの「半教師あり」はむしろ「弱教師あり」（weakly supervised）の意味に解釈すべきだろう．contrastive loss あるいは triplet lossを効率よく計算しようとするならば，ラベルの所与を前提に置かねばならない．件の分布ペアを訓練時のミニバッチ内で都度作り出す（on-the-fly）場合は特にそうである．ただしラベルとは分布相互の近隣・遠隔に関する何らかの情報を指す．本記事の実験ではクラスタラベルをそのまま利用している．いずれにせよ，各分布をあたかも個々の観測データと見なすという点を除けば，従来の深層距離学習と同じ枠組みに帰着できる．

実験： Clustering Multivariate Gaussian Distributions

実験条件

データセット生成

分布クラスタリングのデモンストレーション実験のために，著者らは多変量ガウス分布の群を取り上げた．以下，図2を参照しながらデータの作り方を説明しよう．

図2: 分布の平均ベクトルおよび散布図

半径を変えた3つの同心円の円周付近にガウス分布の平均ベクトルを500個発生させる（訓練用; テスト用はこれとは別に200個）．同心円の半径はそれぞれ0.2, 0.6, 1.0であり，3つの同心円は3つのクラスタにそれぞれ対応する．1つの平均ベクトルは円周上の点をランダムにサンプリングし，2次元のガウスノイズの重畳により得られる．このときのノイズの大きさは論文に明記されていなかったが，本実験では各次元で標準偏差0.05とした．以上が図2の一番左に対応する．色付けにより，対応するクラスタを見やすくしている．

続いて各ガウス分布から複数個の点をランダムサンプリングする; 各ガウス分布の共分散の大きさについては幸いにも論文に記述があり，対角共分散行列を仮定し，対角成分はそれぞれ0.1とした．$\sqrt{0.1}$でスケーリングした2次元ガウスノイズを独立に複数個発生させ，各平均ベクトルに加えているわけである．論文には具体的にいくつ点をサンプリングするのかについて明確な言及がなかったが，論文のFigure 1には各ガウス分布から50個の点をサンプリングしたときの図が掲載されていた点を鑑みて，本実験でも1つのガウス分布あたり50個の点をサンプリングした． 図2の真ん中3つの散布図はそれぞれ，各クラスタのガウス分布からランダムサンプリングした点群を示している．一番右はこれら3つの散布図を全て重ねた図である．平均ベクトルのサンプリングと，対応するガウス分布からのサンプリングという2段階の構造になっているので気をつけて欲しい．

最終的に，手元には500 × 50 = 25,000個からなる訓練用の点群と200 × 50 = 10,000個からなるテスト用の点群があり，各点は3つのクラスタのいずれかに属する（＝ラベルが付与されている），という設定である．図2を見て分かる通り，分散が大きいため各クラスタの点群は互いに大きく重なっている．

比較手法

の3つを比較する．従来の深層距離学習では，個々のサンプル点は独立な観測データとして扱われ，クラスタリングもまた観測データ単位で行われる．"Davis and Dhillon"および提案手法の枠組みでは経験分布をクラスタリングの対象にするために，1つのガウス分布から発生させた50個のサンプル点が1つの経験分布を構成するよう，各点と経験分布とを一意に紐づけて管理する（下記補足参照）．各経験分布のクラスタラベルは所与とする．

"Davis and Dhillon"では経験分布から計算される標本平均ベクトルと標本共分散行列に基づいてアルゴリズムを駆動する．教師なしのクラスタリングアルゴリズムなので，訓練時にクラスタラベルは使用されない．評価用のreference（正解ラベル）として利用する．つまり分布クラスタリング後に，経験分布に対してクラスタラベルが予測されるので，referenceとの比較により性能を評価する．

深層距離学習に関する2手法では所与のクラスタラベルを活用してサンプル点ペアないし分布ペアを作成する．訓練が完了した後に訓練済みのニューラルネットによって各データの埋め込みを計算する（非線形変換）．そして変換後の埋め込み空間上で$\mathit{k}$-meansアルゴリズムを駆動する．学習済みの「距離関数」に基づいた教師なしクラスタリングというわけである．予測されるクラスタラベルは，従来の深層距離学習では各点（観測データ）単位で算出されるが，提案手法では経験分布単位で算出されるので，後述する評価指標の計算時には注意が必要である．

ネットワークアーキテクチャ および 損失関数

埋め込みネットワークは入力層に続いて隠れ層が2層続く構成である．ユニット数は入力層から最終層に向かって順に(2)-(1000)-(500)-(2)とした．最終層のユニット数が2となっているのは埋め込み空間の可視化を容易にするためである．層の間の活性化関数はReLUを用いた．

損失関数はcontrastive lossおよびtriplet lossである．

評価指標

クラスタリングの性能評価のために以下の2指標を用いた．参考まで，scikit-learnへのページをリンクしておく．

  • Rand Index (RI) URL

  • Adjusted Rand Index (ARI) URL

これらはいずれも，同じデータセットに対する2つのクラスタリング結果の間の類似度を計算するものである．ARIは特にクラスタラベルのランダム割り当てにかかるバイアスを補正した（adjusted）RIである．実験では予測クラスタラベル（predicted）と真のクラスタラベル（ground truth, reference）の間の類似度を計算する．

その他の実験条件
項目 値
エポック数 50
最適化アルゴリズム Adam
学習率 1e-4
バッチサイズ 32
triplet loss のマージン 0.2
contrastive loss のマージン 0.2
補足：経験分布の構成

経験分布の構成を簡単に式に書いておこう（本記事オリジナル）．経験分布の集合$E$を

$$ E = \left\lbrace e_{1}, e_{2}, \ldots, e_{N} \right\rbrace $$

と書いてみる．$N$は経験分布の総数であり，訓練用は$N=500$である．テスト用には$N=200$として別に作成しておく．各経験分布が50個のサンプル点を持っている． つまり$e_{i}\; (1 \leq i \leq N)$ には

$$ \left\lbrace \mathbf{x}^{(i)}_{1}, \mathbf{x}^{(i)}_{2}, \ldots, \mathbf{x}^{(i)}_{50} \right\rbrace , \;\;\; \mathbf{x}^{(i)}_{j} \sim \mathcal{N} \left( \mathbf{m}^{(i)}; \begin{bmatrix} 0.1 & 0 \\ 0 & 0.1 \end{bmatrix} \right) \;\;\; (1\leq j \leq 50) $$

という点集合が付随している．ただし$\mathbf{m}^{(i)} \in \mathbb{R}^{2}$は平均ベクトルである．これらを用いて

$$ e_{i} (\mathbf{x}) = \frac{1}{50} \sum_{j=1}^{50} \delta \left(\mathbf{x} - \mathbf{x}^{(i)}_{j} \right),\; \; \; \mathbf{x} \in \mathbb{R}^{2} $$

書けるものとする．$\delta$はDiracのデルタ関数である．各経験分布はクラスタのラベルを有している．便宜上，3クラスタのラベルを$0, 1, 2$として表せば，$E$からラベル集合$C = \left\lbrace 0, 1, 2 \right\rbrace$への写像$g$が定義されている．

$$ g: E \to C; \;\; e_{i} \mapsto c \in C $$

それゆえにまた，各サンプル点にも，$e_{i}$を介してラベルが付与される．

実験結果および考察

実験結果

評価指標ごとに表を分けて示す．baselineは従来の深層距離学習，proposedは提案手法である． 平均および標準偏差は乱数の種（seed）を10回変えて得られたスコアから計算した． なおseedの値が関わるのはデータセット生成の部分のみである（経験分布を生成する際のサンプリング）． もしseedを固定すればたとえGPUを使用しても常に同じ結果が得られる．

  • RI
手法（損失関数） 平均 標準偏差
Davis & Dhillon 0.649826 0.008840
baseline (constrative loss) 0.650168 0.013206
proposed (constrative loss) 0.990864 0.007971
baseline (triplet loss) 0.657529 0.012002
proposed (triplet loss) 0.992834 0.006162


  • ARI
手法（損失関数） 平均 標準偏差
Davis & Dhillon 0.241574 0.009874
baseline (constrative loss) 0.238745 0.030552
proposed (constrative loss) 0.979434 0.017977
baseline (triplet loss) 0.238126 0.028162
proposed (triplet loss) 0.983881 0.013887

どちらの評価指標においても提案手法は"Davis & Dhillon"とbaselineを大きく上回った結果を確認できる．以上より，論文の実験結果をほぼ再現できたが，"Davis & Dhillon"の手法は今回の実験でより高いスコアを出した．

考察

すでに述べた通り，経験分布に属する各サンプル点の埋め込み平均によって，moment-matching型のfunctional Bregman divergenceを計算しながら深層距離学習を遂行するというのが提案手法の枠組みだった．今回のデータセットに関していえば，あたかも図2の一番左にプロットした各分布の平均ベクトルたちを対象に深層距離学習を行っているようなものである．平均ベクトルたちはクラスタ間で重なるようには生成されていないわけで，このデータセット上では提案手法が上手く動いたと説明できる．また経験分布単位のクラスタリングという観点で"Davis and Dhillon"と比較した場合，今回のような線形分離不可能な経験分布の初期配置であっても，クラスタリングがうまく働くように非線形変換が学習できており，提案手法は有利であるといえる．

baselineは25,000個全ての点を対象にクラスタリングしていたわけだが，もし，baselineについても埋め込み後に経験分布単位で平均ベクトルを計算し，得られた500個の点を対象にクラスタリングしたらどうなるだろうか？ 独自に実装したからこそ，このような実験も新たにできるわけである． その結果は次の表に示す通りとなった．

評価指標（損失関数） 平均 標準偏差
RI (constrative loss) 0.972352 0.017127
RI (triplet loss) 0.986995 0.005775
ARI (constrative loss) 0.937757 0.038647
ARI (triplet loss) 0.970741 0.013020

baselineのスコアは大幅に向上しproposedとほぼ同等になったが，「平均」ではわずかながらproposedが勝る結果となった．以上より，baselineにおいても深層距離学習は問題なく動いていることが確認できた．

おまけ　埋め込みの可視化

まずbaselineに関して，図3と図4に埋め込みの可視化をそれぞれ示す．ただし訓練データに対する埋め込みの可視化であり，埋め込みネットワークはtriplet lossで学習したときのものを使っている． 図3は25000点すべてをプロットした場合，図4は埋め込み後に経験分布単位で平均ベクトルを計算した場合にそれぞれ対応している．

図3ではクラスタ間の重なりが大きいように見える．したがってクラスタリングのスコアも大幅に低下したと考えられる．図4の経験分布単位ではクラスタごとに分かれた埋め込みが得られており，やや広がりを持って散布している．このような埋め込みが得られた理由は，損失関数が分布単位で計算されておらず，あくまでも単独のサンプル点単位で計算されているから，つまり「分布単位で点がまとまる」ように変換は学習されないから，と一応の説明がつけられる．

図3：埋め込みの可視化（baseline，すべてのサンプル点）

図4：埋め込みの可視化（baseline，経験分布単位の平均ベクトル）

最後にproposedに関して図5に埋め込みの可視化を示す．スコアの高さが示す通り，クラスタがきれいにまとまり，かつきれいに分かれた埋め込みが得られている．

図5：埋め込みの可視化（proposed）

実装

以下のリポジトリに置いた．深層距離学習の実装にはPyTorch Metric Learningが大いに役立った．動かし方の詳細はREADMEに書いてある．本記事では実装の詳細をいくつか抜粋して紹介する．

github.com

clustering/gaussianに置いてある．

  • moment-matching型 functional Bregman divergence のクラス MomentMatching:
    model.pyにある．PyTorch Metric LearningのBaseDistanceを継承し，compute_matとpairwise_distanceを独自に実装する．平均ベクトルを計算する処理が特別に入る以外はユークリッド距離の計算とほぼ同じであり，実装自体は容易である．

  • データセットのクラス DistributionDataset:
    dataset.pyにある． moment-matching型のdivergenceの計算には経験分布に関する点集合を伴う． それら点集合たちをリストで管理する．リストの各要素が2次元の行ベクトル × 50個つまり[50, 2]のarrayである．実験条件の箇所で述べた通り，リストの全長は訓練用が500, テスト用が200である．なお従来法（深層距離学習のbaseline）はこのような特別のデータセットは使わずとも，全ての点をフラットに扱う一つの巨大なarray（例えば訓練用は[25000, 2]のサイズ）を作り，TensorDatasetでラップするだけで済む． ただし推論時に経験分布単位のフェアな比較をする場合はbaselineと提案法ともにDistributionDatasetを使う．

  • 評価指標を計算する util.py
    元論文ではRIとARIを計算していたが，scikit-learnには幸い"Clustering metrics"のカテゴリにて種々の指標を計算するAPIが提供されているので，併せて計算できるようにしておいた．ただしPurity scoreは自前の実装である．

    • Adjusted Mutual Information URL
    • Normalized Mutual Information URL
    • Homogenity score URL
    • Purity score
    • クラスタリング後の分類性能という観点で，$\mathit{k}$-Nearest Neighbor (kNN) 分類器に基づくaccuracyとROC AUCも計算できるようにしてある．

  • "Davis & Dhillon"の分布クラスタリングを実行する entropic_clustering.py:
    以前の記事の実装をそのまま持ってきて，今回の実験で用いたデータセット上で訓練・評価する．この記事はまさに本記事の布石であった．Python実装が転がっているはずもないので，自前で実装する必要があったからだ．

tam5917.hatenablog.com

おわりに

今回の記事では分布クラスタリングに関して，実験の一部分について再現に成功した（と思う）．続く記事にて，分布クラスタリングに関して残された実験（"Human Activity Recognition Task"）についても報告する．

参考

著者らによるスライド＆プレゼン動画

slideslive.com

*1:ちなみに指数型分布族と（正則な）凸関数の族が対応することが知られており，またBregman divergenceを与えたとき，それをKL-divergenceに持つ指数型分布族が存在することも知られている．

*2:Davis, J. and Dhillon, I. S. Differential entropic clustering of multivariate Gaussians. In Advances in Neural Information Processing Systems (NIPS), 2006.

"Differential Entropic Clustering of Multivariate Gaussians"をNumbaを使って高速化してみた話

機械学習 python

はじめに

前回記事で実装した Differential Entropic Clustering をもう少し高速化したいなぁ，という話．

tam5917.hatenablog.com

実装

やり方は簡単で，numbaをインストールして，@jit デコレータをBurg matrix divergence およびMahalanobis距離を計算する関数につけるだけ．

@jit(nopython=True)
def comp_burg_div(mat_x, mat_y):
    """Compute Burg matrix divergence."""
    dim = mat_x.shape[0]
    mat_y_inv = LA.inv(mat_y)
    mat = mat_x @ mat_y_inv
    burg_div = np.trace(mat) - np.log(LA.det(mat)) - dim
    return burg_div

@jit(nopython=True)
def comp_maha_dist(vec_x, vec_y, cov):
    """Compute Mahalanobis distance."""
    delta = vec_x - vec_y
    m = np.dot(np.dot(delta, LA.inv(cov)), delta)
    return np.sqrt(m)

前回記事のdemoスクリプトを上記の関数たちに置き換えて，time moduleで簡単に計測してみたところ，

  • 高速化なし ... 58.50秒
  • Burg matrix divergence を計算する関数のみ高速化 ... 33.59秒
  • Mahalanobis 距離を計算する関数のみ高速化 ... 45.33秒
  • 両方の関数を高速化 ... 12.31秒

となった．

高速化に成功した．行列計算関連で繰り返し呼び出される関数には，jitデコレータをつけておくのが良さそうである．

以下，参考まで，今回の実装．main関数の前後で時間計測する（面倒なので）．

A demo script of "Differential Entropic Clustering of Multivariate Gaussians" in Python. This is a faster version thanks to numba and jit. · GitHub

おわりに

Numbaによる高速化の効果は大いにあったので，薦めたい．

Differential Entropic Clustering of Multivariate Gaussians (NIPS 2006) をPythonで実装した

機械学習 python

論文はこれ．

proceedings.neurips.cc

多変量ガウス分布の平均と共分散行列の集合が与えられたときに，分布を単位として（＝対応する平均と共分散をペアにして）クラスタリングするアルゴリズムが提案されている． 行列の距離を測るための"Burg matrix divergence"が利用されているのが興味深い．これは行列を対象にしたBregman divergenceに相当する．

実装してみたのがこれ．要progressbar2．

Section 4.1 "Synthetic Data"を参考に，simplex上から平均ベクトルをクラスタ数だけランダムにサンプリングする．共分散行列は対角共分散にして，かつ対角成分の大きさも論文を参考にした．論文と設定が異なるのは，クラスタリングの対象とするのが標本平均と標本共分散という点である．つまり平均ベクトルと共分散行列で指定されるガウス分布から点を指定個数だけサンプルし，それらサンプル点に基づいて標本平均と標本共分散を計算する．論文よりも厳しい条件で分布クラスタリングを試していると言える．

クラスタリングの性能指標は以下を計算できるようにした．Purity Score以外は全てscikit-learnから利用できる． 参考までにAccuracyも計算している．

  • Rand Index
  • Adjusted Rand Index
  • Adjusted Mutual Information
  • Normalized Mutual Information
  • Purity Score
  • Homogenity Score

試しに動かしたらこんな感じ：

  • Rand Index = 0.776869 ± 0.071684
  • Adjusted Rand Index = 0.420437 ± 0.191200
  • Adjusted Mutual Information = 0.485865 ± 0.204090
  • Normalized Mutual Information = 0.494591 ± 0.200726
  • Purity Score = 0.630800 ± 0.135231
  • Homogenity Score = 0.485421 ± 0.196300

乱数の種を50回変えたときの結果（平均と標準偏差）である．データの次元数は4，クラスタ数も4である．今回はスコアのバラつきが大きくなった． スクリプト内のget_dataset関数をカスタマイズし，各クラスタの平均ベクトルの取り方や共分散行列の取り方を変えてみると，バラつきの傾向も変わるので色々と楽しめる．例えば，今回のように平均ベクトルをsimplex上に限定する必要はない．

クラスタリング結果の可視化はそのうち．

距離学習と深層距離学習の違い

距離学習 深層学習 機械学習

距離学習と深層距離学習の違いについて備忘録を残しておく．

共通していること

データセット中の2点 $\mathbf{x}, \mathbf{y}$間の距離関数$d(\mathbf{x}, \mathbf{y})$が与えられる． 距離関数はクラス分類やクラスタリング・異常検知などのタスクで活用される道具である． それらタスクに有利になる形で特徴量の変換を考えたい． その変換をパラメタ（行列など）で表現し，種々の損失関数 with 正則化項 に基づいてパラメタを最適化する． 距離に関する「望ましさ」は（不等式などの）制約によって表現される．例えば同一クラスのサンプルペアは近い距離，異なるクラスのサンプルペアは遠い距離にする制約である． したがって各サンプルには所属クラスなどのsupervisedなラベルが与えられることが典型である．

違い

代表的なユークリッド距離で考える．

  • 距離学習：$d(\mathbf{x}, \mathbf{y}) = ||G\mathbf{x} - G\mathbf{y}||_{2}$となる線形変換$G$（半正定値行列）を獲得する．

  • 深層距離学習：$d(\mathbf{x}, \mathbf{y}) = ||f_{\theta}(\mathbf{x}) - f_{\theta}(\mathbf{y})||_{2}$となる非線形変換$f_{\theta}$を獲得する． この非線形変換がニューラルネットにより実現され，$\theta$はそのパラメタである．

つまり距離学習は行列の最適化問題に帰着される．深層距離学習は当たり前だがニューラルネットのパラメタ最適化問題に帰着される． しばしば混同されるが，距離関数そのものはpre-definedであり最適化の対象ではない．この意味で距離関数の選択はハイパーパラメータといえる．

線形であれ非線形であれ，変換された特徴量はしばしば「埋め込み」（embedding）とも呼ばれる． 距離学習では特徴量エンジニアリングが依然として重要であるが，深層距離学習では埋め込みの獲得が特徴量抽出を兼ねている面がある．

では距離関数自体を学習可能にする枠組みは存在するのか？→Yes.