機械学習プロフェッショナルシリーズ-関係データ学習 　読後メモ

サマリー

関係データとは、N個のオブジェクトがあったときに、その中の任意のペアに対してどのような「関係」が存在するかを表すデータとなる。この関係データを表す形として行列を考えると、関係性は行列要素そのものとなり、関係データ学習とはこの行列データのパターンを抽出する学習とも言える。

これらが適用されるタスクには「予測(prediction)」と「知識抽出(knowledge xtraction)」の2つがある。

予測問題とは、観測データから学習、設計された統計モデルを用いて未観測のデータの値を推定する問題であり、知識抽出問題とは、関係データ自体の特性を解析することや、与えられた観測データを適切にモデリングすることによって、何らかの有用な知見やルール、知識につながる情報を抽出するタスクとなる。

ここではこの関係データ学習に関して機械学習プロフェッショナルシリーズ「関係データ学習」をベースに様々なアルゴリズムとその応用について述べる。

今回は読後メモについて述べる。

機械学習プロフェッショナルシリーズ-関係データ学習 　読後メモ

第1章　導入：関係データ解析とは

1.1 統計的機械学習
1.1.1 データからの学習
“ラベル”情報が重要
データおよびラベル情報からなる 観測データ{X,y}={(xi,yi),i=1,2…を用いて
未知のデータxからラベルyiを 計算する関数(写像)を推定する
1.1.2 教師有り学習と教師無し学習
画像データでは”ラベル”が付けられて、 それがなるべく再現できるような学習を行う
教師情報(データ、ラベル)
教師あり学習 (supervised learning)
教師なし学習 (unsupervised learning)
1.2 関係データとは
本種では教師なし学習を扱う
例 コミュニティ分析
ハイパーリンク、感染症の拡大等
同じ種類のオブジェクト間の関係だけではなく、 異なる種類のオブジェクト間の関係もある
例:購入履歴データ
「商品」と「顧客」
「文書」「キーワード」「参照文書」「著者」
オブジェクト間に非常に多数の「関係」が 定義可能な場合に対する技術は発展途上
30,48
関係データのデータセットの入手先
Stanford Network Analysis Project
1.3 関係データの表現
1.3.1 関係データのグラフ表現
グラフ理論 (graph theory)
1.3.2 関係データの行列(多次元配列)表現
隣接行列 (adjacency matrix) による表現
本書での扱い
接続行列 (connectivity matrix) による表現
関係データの行列表現
線形代数によるアプローチが可能となる
1.3.3 関係の値の表現
関係がある/なしでは2値の(離散)数値が与えられる
閾値以上か、閾値以下かの関係
各関係を離散的な記号(シンボル)で値づけることも可能
シンボルはインデックス付けされる
グラフの表示の場合は3つ組み(e=(v1,v2,v3))を使う
行列の場合は、xij=s(シンボル値), xij=0(整数)などの表記を行う
観測データと未観測データの区別が大事
観測できる(observed)
観測できない(missing)
1.4 関係データの種類
1.4.1 有向関係データと無向関係データ
関係性の方向
有向(directed)な関係データ
xij ≠ xji
無向(undirected)な関係データ
1.4.2 単一ドメインと複数ドメイン
関係データ内のオブジェクトのドメインについて
単一ドメイン(single domain)
全てのオブジェクトが同じ種類である
行と列が同じドメイン
正方行列
行と列のサイズが同じ
複数ドメイン(multiple domains)
異なる種類のドメイン
頂点ペア(i,j)と(j,i)は意味が異なる
データが2部グラフ(bipartite graph)の場合
ドメイン内のオブジェクトは関係が定義されない
ドメイン間のみ関係が定義される
1.4.3 対称関係データと非対称関係データ
対称関係データ (symmetric relational data)
行列表現した時に関係データ行列が必ず対称行列になるもの
単一ドメインかつ無向関係データ
非対称関係データ (asymmetric relational data)
複数ドメインあるいは有向関係データ
1.4.4 2項関係と多項関係
2項関係 (two-place relation, dyadicrelation, binary relation)
2つのオブジェクトの関係を表す
2項関係:2次元行列
多項関係 (military relation)
3項関係:3次元行列
N項関係(N-ary relation):N次元行列
多項関係の実用列
時系列データ(time series data)
多次元拡張としてのテンソル(tensor)
1.5 関係データ解析
関係データ解析のタスクカテゴリ
1.5.1 予測
予測(prediction)
観測データから学習、設計された統計モデルを用いて未観測のデータを推定する
典型的な予測問題
リンク予測(link prediction)問題
関係ネットワーク内に存在するリンクの有無の推定
購買データを用いたユーザごとのアイテム購入(sdoption)確率の推定
関係データ行列の欠損値の予測
ネットワーク内の情報伝搬(information dissemination)の推定
ネットワーク内の情報拡散(information diffusion)の推定
アイテム推薦の例
ユーザーがまだ試聴したことのない映画を推薦する問題
欠損値予測
1.5.2 知識抽出
知識抽出 (knowledge extraction)
知識発見 (knowledge mining, knowledge discovery)
グラフ特徴を計算することで関係データ自体の特性を解析する
与えられた観測データを適切にモデル化することにより何らかの有用な知見や知識につながる情報を抽出する
例
ネットワーク内のコミュニティ抽出(community extraction)
ネットワークのクラスタリング(clustering)
多くは教師なし学習
グラフの特徴量により関係データの性質を明らかにするという目的の場合
すでに定義されているグラフ特徴量を計算するだけなので教師データや学習は必要ない
本当にクラスタが存在するかどうかわからないので教師データが作成できない
「もしクラスタが存在するならばこのような性質に従うはずだ」という利用者の仮説(意図)を反映する
具体的な知識抽出例
関係データクラスタリング
特徴空間ないで近くにあるデータをまとめる(だけ)の操作
人間の目で簡単に区別がつくクラスタリング例
関係データのクラスタリング例
1.5.3 本書のアプローチ: 低次元構造をとらえる
実関係グラフは完全なランダムグラフではなく、 何らかの低次元構造が埋め込まれているととう仮説
モデルにマッチする部分(構造)と マッチしない部分(ノイズ)に分離
代表的な手法
関係データ行列のブロック構造
関係データをいくつかのブロックに分割
関係データの行インデックスと列インデックスを それぞれK<N個程度のグループに分割
同じグループのオブジェクトが 隣接するように行と列を並び替える
関係データ行列の低ランク性
NxNの関係データ行列をK<N個の基底ベクトルの線形和で表現
N個の関係パターンをK種類の製材パターンで表現
1.6 本書の目的と構成
行列の基礎
基本概念
行列とは
行列積
積(掛け算)が常に計算できるとは限らない
Aがaxb行列、Bがcxd行列(axb行列とは縦にa個、横にb個並んでいる)
積が計算可能であるとはb=c
積はaxd行列になる
零行列と単位行列
何を掛けても0になるもの
AO = OA = O
何を掛けても変化しないもの
AI = IA = A
正方行列しか存在しない
逆行列
行列には割算はない
AR = RA = I となるようなRをAの逆行列A-1
元の行列が正方行列でないと定義できない
直行行列
転置行列が逆行列と等しい
転置行列
行列の左上と右下を掴んで、息を吹きかけて180°回転させたもの
AtA = AAt = I
複素数の範囲で考えた直行行列をユニタリ行列と呼ぶ
対角行列
対角成分にだけ数字が並んでいるもの
例
ベクトル
縦か横かどちらかの長さが1であるような行列
例
操作
対角化
行列Aを対角化するということ
直行行列Pを持ってきて、対角行列Dについて、 A = PDPt とすること
D = PtAP
直行行列で挟み上げした結果対角行列になる
何が嬉しいのか
対角化することでAを単純なベクトルの組み合わせに変換可能
直行行列は空間上の回転・反転を表している
回転した先の世界で考えて後で戻す
対角行列は各方向だけが関与していて方向同士が相関を持っていない
行列の相似にも適用する
対角化、固有値分解、”特異値分解(Singular Value Decomposition, SVD)の概要とアルゴリズム及び実装例について“でも述べている特異値分解
固有値分解
特異値分解
固有値分解を一般の行列に拡張したもの
対角化、ジョルダン標準形

第2章　対称関係データのクラスタリング技術：スペクトラルクラスタリング

2.1 関係データのクラスタリングとは
知識発見タスクの一つとしてのクラスタリング(clustering)
多数のオブジェクトが与えられた時に、 何らかの意味で類似する特徴を持つオブジェクトを集めて、 いくつかのグループに分類するタスク
データセットx={x1,…xN}が与えられた時、データ群をK個のクラスタに分類
手法
K-means
データセットが通常のd次元連続特徴ベクトルの集合の時によく使われる
mean shiftクラスタリング
Latent Dirichlet allocation
文書データのような離散シンボルの集合で表される場合
関係データ解析でのクラスタリングの対象
オブジェクト(頂点、関係の主体となるもの)か エッジ(辺、リンク)
本書ではオブジェクトのクラスタリング
類似性の3つの定義
各オブジェクト自身i,jがもつ、関係によらない特徴量の類似性
各オブジェクト(頂点)が関係情報以外の 類似性を持つことが前提
K-means等で分析
各オブジェクトi,jが他のオブジェクトと どのような関係を根持っているかによる類似性
各オブジェクトが匿名化
それぞれのオブジェクトを特徴付けるような情報(属性等)も見当もつかない
利用できる情報は関係性のみ
辺で繋がった頂点集合の類似性等を利用
上記の2つの両方を考慮した類似性
2.2 対称関係データのオブジェクトクラスタリング法: スペクトラルクラスタリング
見つけたいクラスタの性質や 関係データの特徴に応じて利用可能な手法が異なる
2.2.1 コミュニティ検出と密結合グラフ
コミュニティ検出 (community detection)に着目
コミュニティ
そのクラスタ内部では密な結合を持つ一方で クラスタ外のオブジェクトとは非常に疎な結合を持つ頂点のクラスタ
密結合クラスタ (assortative cluster)
疎結合クラスタ (disassotative cluster)
疎結合クラスタでセレンディピティができる?
2.2.2 グラフカット
あるクラスタが密結合か否か判断したり、 密結合クラスタを発見したりするためには、 どのような特徴量を計測すれば良いか
特徴量の一つ
グラフのカット(cut)
グラフカット(graph cut)関数
関係データ行列あるいは類似度行列(affinity matrix)によるグラフデータの頂点の分割に対して定義される
例
N個のオブジェクトをK個のクラスタに排他的に分割
K番目のクラスタに所属する頂点のインデックス集合をPk
K番目以外のクラスタに所属する頂点のインデックス集合をṖk
分割P={P1, P2, …, Pk)に対するグラフカット関数
異なるクラスタに属する頂点間の辺の重みを表す
カットの最小化(mincut)
頂点集合を繋がりの弱いところで分割することで
グラフ内の密結合クラスタを抽出する
傾向としてネットワークから孤立したクラスタと その他のような分け方になる場合が多い
各クラスタがある程度の大きさに なるように正規化したカットを行う
2.2.3 スペクトラルクラスタリング (spectral clustering)
正規化を伴う最小化カットを求める問題はNP困難
現実的な時間で探索することが困難
離散的な正規化カット最小化問題を 連続緩和して解く
スペクトラルクラスタリング (spectral clustering)
グラフを隣接行列で表現
グラフ内の2つのノード間にリンクがあれば1, リンクがなければ0の値を持つ関係データ行列
隣接行列から導かれるグラフラプラシアン(graph Laplacian)と呼ばれる行列の固有値に注目
第2最小固有値およびそれ以降の固有値の値によって グラフ内のノード間の連結性が推定できる
固有値分解とk-meansクラスタリングによりノードのクラスタリングを実現
スペクトラルクラスタリングは対象関係データ行列にのみ適用可能
例
2.3 非正規化グラフラプラシアンによる スペクトラルクラスタリング
最も単純なアルゴリズムに基づく スペクトラルクラスタリングを説明
2.3.1 入力と出力
入力
N個のオブジェクトの対称関係データX=(xi,j)
ノード数(オブジェクト数) N∈ ℕ
対称データなのでインデックスi,jの取る範囲は一致
簡単のために、Xi, j ∈ {0, 1}として議論
想定されるクラスタ数K ∈ ℕ
解析のために事前に想定される潜在クラスタ数を与える必要がある
出力
各オブジェクトのクラスタリング結果
主産物
N個のオブジェクトを K個の密結合クラスタ(コミュニティ)に分類した結果
副産物
各オブジェクトiを、 K次元連続空間に埋め込んだ場合の特徴ベクトルvi
オブジェクトの「位置」にあたる
入出力例
2.3.2 次数行列とグラフラプラシアン
具体的なスペクトルクラスタリングの計算過程の説明
入力された関係データ行列Xから 次数行列(degree matrix)を計算する
次数行列D=(di,j)は対角行列
非対角要素は全て0
I番目の体格要素di,jはX内のオブジェクトiの持つ次数 (すなわちリンク数あるいはリンクの重みの総和)
次数は オブジェクトの受けるリンク総数である入次数(in-degree)と オブジェクトの伸ばしたリンク総数である出次数(outーdegree)と 区別する必要がある
スペクトルクラスタリングは対称関係データのみなので入次数と出次数は同じ
NxN対角行列である次数行列Dは
次数計算とラプラシアン計算を行なった例
次に、グラフラプラシアンと呼ばれる行列を計算する
グラフラプラシアンにはいくつか種類がある
単純なグラフラプラシアンL
L = D – X
非正規化グラフラプラシアン (unnormalized graph Lalacian)
2.3.3 固有値分解によるクラスタ抽出
グラフラプラシアンの結果に、 似たような値を持つ行(列)が現れる
グラフラプラシアンから 「似たような」パターンを計算できれば オブジェクトのクラスタリングができる
パターンの計算に 「固有値分解(eigenvalue decomposition)」を使う
固有値分解の復習
「固有値分解」は最も基本的な線形代数の一つ
対称行列M ∈ ℝNxNが与えられたとする
K < Nに対して
M =VΛVTが成立するような 列直行行列V∈ℝNxKと 対角行列Λ∈ℝKxKを求める
VTはVの転置(transpose)
対角行列Λの各対角要素を行列Mの固有値(eigenvalue)
Vの各列ベクトルを行列Mの固有ベクトル(eigenvector)
λkvk = Mvkが成り立つ
固有ベクトルはそれぞれが互いに線形独立
行列Mの各列ベクトルは固有ベクトルの線形結合で表現される
λ1≤λ2≤…≤λk
スペクトラルクラスタリングでは グラフラプラシアンLの固有値分解を行う
L = VΛVT
あるいはλkvk = Lvk
グラフラプラシアンのような半正定値行列では、 固有値は全て非負の実数となるので
Λ1≤λ2≤…≤λkが設立する
固有値にそれぞれ対応する固有ベクトル行列V=(v1,v2,…,vk)∈ℝNxKでは
各i番目の行が、i番目のオブジェクトのK次元圧縮表現となる
Kは想定しているクラス多数
なぜこれで良いのか?
直感的な解釈
K個のクラスタを持つグラフの場合
最低K個の固有値は0になる
それ以降の固有値は非ゼロの整数となる
最小K個の固有ベクトルがクラスタの指示ベクトルになる
2.3.4 クラスタ割り当てZの計算
上記は 「各クラスタ内は全結合している」と 「2つのクラスタ間には辺が一切存在しない」 という仮定での議論
実際の例ではあまりケースとしてない
固有値を計算しても最小値が0とはならない
できるだけ0に近いほうが理想的なクラスちに近づく
具体的な初期値を考えるためにk-meansを適用する
VをN個のK次元ベクトルと捉える
K次元空間でk=Kとするk-meansを実施
2.3.5 まとめ:アルゴリズム
固有値分解とk-meansだけで計算可能
入力パラメータはKだけ
2.4 正規化グラフラプラシアンによる スペクトラルクラスタリング
2.4.1 正規化カットとスペクトラルクラスタリングアルゴリズムとの関係
頂点の数が大きすぎず小さすぎず、 かつ密結合な特性を持ったクラスが少数あるような状態が理想
各クラスタの大きさに制約をかける
数多くある正規化の中で有名なもの
レシオカット (Ratio Cut)
各クラスタkに関するカットを クラスタのサイズ(クラスタ内の頂点数)で除して正規化
クラスタが大きいほどカットが小さくなる
単一の頂点からなるクラスタは排除される
ノーマライズカット (Normalized Cut, NCut)
クラスタ内のリンクの重みの総和で除して正規化
クラスタが大きくかつ相互に密に結合しているほどカットが小さくなる
典型的なコミュニティ検出ではこちらが適している
2.4.2 対称正規化グラフラプラシアンに基づくスペクトラルクラスタリング
ノーマライズドカットで、グラフラプラシアンの行列を用いた2次形式に書き換え
その2次形式を最小化する解を制約付きで解く
新しいスペクトラルクラスタリングアルゴリズムになる
Lsym = D-1/2LD-1/2
対称正規化グラフラプラシアン (symmetric normalized graph Laplacian) アルゴリズム
2.4.3 酔歩正規化グラフラプラシアンに基づくスペクトラルクラスタリング
シーとマリクによって提案された、 異なる正規化グラフラプラシアンに基づく スペクトラルクラスタリング
提案されたグラフラプラシアン
Lrw = D-1L
酔歩正規化グラフラプラシアン (Random-walk normalized graph Laplacian)
なぜ酔歩(random-walk)になるのか
ノーマライズドカット最小問題が、 グラフを複数のクラスタに分割したときに、 クラスタ間の遷移確率を最小化する、 互いに独立したコミュニティを求める問題に帰結する
アルゴリズ
2.5 実データへの適用例
Zachary’s Karate Club network dataへの適用
各頂点はある大学の空手クラブのメンバ
リンクの有無は空手クラブ以外での交友の有無にたいえう
3種類のスペクトラルクラスタリングを適用した結果
クラスタリング精度の指標の一つ
修正Randインデックス (Adjusted Rand Index, ARI)
Lrw, Lの場合はARIが0.5725
Lsymの場合はARIが0.3291
性能が低い
2.6 実運用上の留意点と参考文献
2.6.1 どのアルゴリズムを選択するか
まずは酔歩正規化ラプラシアンを選ぶべき
2.6.2 Kの設定方法
ヒューリスティックな手法
1. 固有値の絶対値に閾値を定める
2. 隣あう固有値の差に閾値を定める
3. 隣あう固有値の比に閾値を用いる
モデル選択基準を用いてKの数を決める
赤池情報基準量(Akaike informato criterion, AIC)
ベイズ情報基準量(Bayes information criterion, BIC)
2.6.3 参考文献について
フォンラクスバーグによるスペクトラルクラスタリングのチュートリアル
84
ソフトウェアとしてはpyhtonのscimitar-learnモジュールに組み込まれている
2.6.4 スペクトルクラスタリングの限界と密結合クラスタリングの現状
無向かつ単一ドメイン関係データ(対称関係データ)にしか適用できない
対象でないケース
有向関係データ
複数ドメイン関係データ
memo
スペクトラルクラスタリング入門
スペクトラルクラスタリングとは
スペクトラルクラスタリングの特徴は、 データからグラフを生成し、 グラフの連結成分分解を応用してクラスタリングするところ
他の学習とは異なる特徴
行列の固有値問題によるグラフの連結成分分解
スペクトラルクラスタリングのアルゴリズムを説明する前に、 その基礎として、 行列の固有値問題によって グラフの連結成分分解を行う方法を説明
前提条件
正の重みつき無向グラフ G=(V = {v1, …, vn}, E)があるとする
隣接行列をW=(wij)とする
wijはwij=wji, wii=0, wij≥0
隣接行列W=(wij)に対して、 そのdegree matrixを上記で定義する
グラフGのunnormalized graph Laplacianを 上式で定義する
wijの値が大きいほど頂点i, jが近いようにする
例
頂点1と頂点2が重み100の辺でつながっていて、 頂点3と頂点4が重み20の辺でつながっているグラフG_1の場合
例
グラフG1の場合、degree matrix Dとunnormalized graph Laplacian Lは
隣接行列は
頂点1と頂点2が重み100の辺でつながっていて、 頂点3と頂点4が重み20の辺でつながっているグラフG_1の場合
この時次の定理が成り立つ
グラフG=(V = { v1, …., vn}, E)がK個の連結成分をもち、
V={v1,…,vc1}⋃{vc1+1},…,vc2}⋃….⋃{vcK-1+1,….,vn}と連結成分分解されるとする
このとき、Gのunnormalized graph LaplacianLの固有値はすべて非負の実数である
また、固有値0に対応する固有空間は、{ }を基底にもつ
例
Lの固有値0に属する固有ベクトルとすると
Lx=0xすなわち上式が成り立つ
このことから、固有値0に属する固有空間は(1,1,0,0)T,(0,0,1,1)Tを規定に持つことがわかる
ここからグラフの連結成分分解を得るには
まず、規定を並べて
この行列を構成する行ベクトルをそれぞれ見ると、
1行目と2行目がどちらも(1, 0)で同じ
3行目と4行目がどちらも(0, 1)で同じ
これはグラフの連結成分分解が {v1, v2}, {v3, v4}で与えられることと対応している
このようにグラフの連結成分を、行列の固有値問題を解くことによって得ることができる
グラフの連結成分分解アルゴリズム
グラフから隣接行列表現W、degree matrixDを計算し、unnormalized graph LaplacianLを計算する。
unnormalized graph LaplacianLの固有値0に属する固有空間の基底u1, …,uKを計算する。
基底ベクトルを並べて行列U=(u1u2 …. uK)を作る。
行列Uのi行目とj行目の行ベクトルが等しければ、頂点viとvjは同じ連結成分に属し、さもなくば違う連結成分に属している。
スペクトラルクラスタリングのアルゴリズムと計算例
スペクトラルクラスタリングのアルゴリズムは、 上に説明したグラフの連結成分分解を応用したもの
具体的なアルゴリズム
クラスタ数Kが与えられているとする
データからグラフを構築する
グラフから隣接行列表現W、degree matrixDを計算し、unnormalized graph LaplacianLを計算する
unnormalized graph LaplacianLの固有値が小さい順に固有ベクトルをK個u1,…, uKを計算する
固有ベクトルを並べて行列U=(u1 u2 … uK)を作る
行列Uの行ベクトルをK-meansなどを用いてK個にクラスタリングする
K-meansの結果、i行目が入ったクラスタを頂点viの入るクラスタとする
連結成分分解との相違点
データからグラフを構築するステップが必要である
連結成分分解のようにi行目とj行目が厳密に等しかったり異なったりしないために、K-meansなどを使う必要がある
例
上図のような連結グラフ
このグラフのunnormalized graph Laplacianは
その固有値は
対応する固有ベクトルは
固有値の小さいほうから固有ベクトルをK=2つとり、並べると
この4つの行ベクトルをK-meansでクラスタリングしたら
1行目と2行目の(1,1), (1, 0.99)のクラスタ
3行目と4行目の(1, -0.97), (1, -1.0)のクラスタに分かれる
v1, v2とv3,v4にクラスタリングされることに対応
K=3で計算すると、固有値の小さいほうから固有ベクトルをK=3個とって
4つ行ベクトルを3クラスタにK-Meansでクラスタリングする
(1,1,1), (1,0.99,0.59)のクラスタ
(1,-0.97, -64)のクラスタ
(1, -1.0, 62)のクラスタ
これはv1,v2とv3とv4にクラスタリングされる
図
関連する話題
スペクトラルクラスタリングではまずデータをグラフに変換します
グラフへの変換の仕方
Ε近傍法
K近傍法
全結合法

第3章　非対称関係データのクラスタリング技術：確率的ブロックモデルと無限関係モデル

3.1 非対称関係データの確率的「ブロック構造」クラスタリング
スペクトラルクラスタリングの問題点の分析
3.1.1 スペクトラルクラスタリングへの不満
問題点
適用可能な関係データが対象関係データ行列だけ
実際のデータはそうでないものが多い
抽出されるクラスタが密結合クラスタ(コミュニティ)のみである
クラスタが全て密結合なクラスタであるとは限らない
抽出するクラスタ数Kの設定が難しい
3.1.2 アプローチ:ブロック構造を仮定した確率モデル
目標: 非対称関係データに対しても適用可能で、 かつ密結合クラスタ以外も抽出可能な 関係データクラスタリングを行う
「潜在的なブロック構造」の存在を仮定する 「確率モデルアプローチ」
「潜在的なブロック構造」とは?
関係データ行列の 行インデックスiの集合と 列インデックスjの集合が それぞれ少数のクラスタに分割できる
関係データの観測値は、 これらのクラスタリングにより規定される
密結合および疎結合クラスタの双方を含む
密結合クラスタは、関係の強い(黒い)ブロックが対角要素に並ぶ
疎結合クラスタは逆に対角要素が弱く、非対角要素に関係が強いブロックが並ぶ
確率モデルアプローチ
潜在的なブロック構造が存在するものと仮定
そのブロック構造の特徴によって、 実際に観測される関係データが変化する
「ある特定のブロックを持つと、 どのようなデータが観測されるか?」 を計算するモデルを何らかのか仮説のもとに作成
右のモデルが確率的(stochastic)に振る舞う
確率的に振る舞うとは
モデルのパラメータを完全に固定しても 観測されるデータ(およびモデル内部の変数)に ランダム性が包含される
真の値とモデルには乖離があるので
確率モデルの持つランダム性は、 想定外の外部要因や乖離の影響を吸収してくれる
目的関数(objective function)
最小化すべきグラフカット関数
ステップ
仮定する確率ブロック構造モデルを決定
確率モデルのための目的関数の選択
最小化あるいは最大化するアルゴリズム
確率モデル一般の場合には、 確率モデルによる識別則の期待損失を 目的関数として最小化する
3.1.3 本章の対象:確率的ブロックモデルと無限関係モデル
確率的ブロックモデル (stochastic block model, SBM)
ブロック構造に基づく関係データクラスタリングの最も基本的な手法
行と列のブロック番号の組み合わせによって、 そのブロック内での関係データの値の分布が変わるモデル
無限関係モデル (infinite relational model, IRM)
SBMの拡張
ノンパラメトリックベイズ (Bayesian nonparametric)を利用
クラスタ数Kを指定せず、 理論上任意数(可算無限個)のクラスタを内包することが可能
無限関係データの簡単な実践例
3.2 確率的生成モデル
3.2.1 確率的生成モデルとは
系を構成する各変数の値が依存関係に従って 順番に生成される過程を確率で表現するモデル
モデル内の各変数の値の決定(生成)は、 指定された確率分布からのサンプリング(sampling)で実現される
単一の「真の値」の存在を考えるよりも、 「大体こういう範囲に散らばっているだろう」 という確率分布(probabilistic distribution)を重視
独立かつ同一に分布 (indépendant and identically distributed, i.i.d.)
例:サイコロを何回も振る
xi | π 〜 p(π) , I = 1,2,…,100
パラメトリックな確率分布
正規分布(normal distribution)
離散分布(discrete distribution)
3.2.2 確率モデルのベイズ推定
統計的機械学習で 確率生成モデルを使う際の、 興味の対象
真のパラメータや確率分布の形状は未知
多くの観測データから確率分布の形状やパラメータの値を推定する
推論(inference)
推論結果も確率分布で表現する
事後分布(posterior distribution)
観測データXを受け取った後の隠れ変数の集合Zおよびパラメータの集合{θ}
事後分布p(Z,θ| X) ∝ p(X | Z, {θ }) p(Z, {θ})
p(X | Z, {θ})
尤度(likelihood)
生成モデルで言うと観測量を生成する数式部
最終的に欲しいものが、 隠れ変数(のいずれか)のみの場合
パラメータ{θ}に基づく分布の不確かさを 全て織り込んで周辺化(marginalization)した事後分布
p(Z | X) ∝ ∫ p(X | Z, {θ}) p(Z, {θ}) dθ
ある変数のとりうる値とその発生確率を全て計算しつくす
「その変数の値がどのような値をとるか? その値によって全体の分布がどのように変化するのか?」 と言う影響を織り込み済みにする (関数の中から変数を消去)
解析的に求めることは困難である
マルコフ連鎖モンテカルロ(Markov chain Monte Carlo)法を使う
3.3 確率的ブロックモデル (stochastic blockmodel, SBM)
3.3.1 SBMの概要
SBMは 関係データ行列が潜在的なブロック構造を持つことを仮定して データ生成過程を構成する確率的生成モデル
各ドメインの「クラスタ間」で関係の強さを定義する
例
観察された関係データ行列
行にN1、列にN2
N1xN2個の関係の観測値を行K個と列L個のクラス間の関係の強さへ抽象化
行と列にそれぞれK,L個のクラスタが潜在していると仮定
少数のクラスタ同士の関係強さパラメータで、 関係データ行列全体を表現
3.3.2 SBMの定式化
前提条件
観測される関係データを{1,0}の2値関係とする
X={xii}, xii∈{1,0}
I=1,…,N1を第1ドメイン(関係データ行列の行)のインデックス
j=1,..,N2を第2ドメイン(関係データ行列の列)のインデックス
k=1,…,Kを第1ドメイン内のクラスタのインデックス
l=1,…,Lを第2ドメイン内のクラスタのインデックス
SBMの確率生成モデル
π1 | α1 〜 Dirichlet(α1)
第1ドメインに潜在するオブジェクトのクラスタの混合割合パラメータ
α1:K次元ベクトル
SBMのハイパーパラメータ
Dirichlet
ディリクレ分布 (Dirichlet distribution)
K次元のベクトルα=(α1,α2,…,αk)をパラメータとするK次元のベクトルの確率分布
例えると、K面のサイコロをランダムに生成する確率分布
具体的な確率分布
𝛤 (*)はガンマ関数(Gamma function)
π2 | α2 〜 Dirichlet(α2)
第2ドメインに潜在するオブジェクトのクラスタの混合割合パラメータ
α2:L次元ベクトル
SBMのハイパーパラメータ
Z1,i = k |π1 〜 Discrete(π1)
K面のサイコロ(各ドメインのクラスタの混合割合)に従って、 各ドメインのオブジェクトのクラスタ割り当てをサンプリングする
Discreteは離散分布を表す
K次元の「サイコロ」ベクトルをパラメータとして、 1からKのうちいずれかの値を一つ生成する
値を選ぶ確率はベクトルの各次元の値の大きさに比例する
具体的な確率分布
𝛱(*)は、 カッコ内の条件が成立するときには1、 成立しない時には0を返す関数
Z2,j = k |π2 〜 Discrete(π2)
L面のサイコロ(各ドメインのクラスタの混合割合)に従って、 各ドメインのオブジェクトのクラスタ割り当てをサンプリングする
同上
Θk,l | a0, b0 〜 Beta(a0, b0)
クラスタ間の関係の強さを定義するパラメータθの生成過程
2値の関係データを想定していることから、 後の計算の容易さを考えてベータ分布 (Beta distribution)を仮定
ベータ分布は0から1の間の実数θに対する確率分布
確率分布
xi, j | {θk,l}, z1,i, z2,j 〜 Bernoulli(θz1,i, z2,j)
クラスタの割り当てを表す隠れ変数z、 およびクラスタ間の関係の強さを表現するパラメータθを用いて、 観測値である関係データ行列X={0,1}N1xN2の 各エントリを生成する確率分布を規定
ベルヌーイ分布(Bernoulli distribution)を 観測値の確率分布として採用
ベルヌーイ分布は簡単に言うと、コイン投げの確率モデル
確率θで値x=1(表), 確率1-θで値x=0(裏)を生成する
確率分布
3.3.3 SBMの推論
周辺化ギブスサンプラーによる推論方法を具体的に導出
(A) 周辺化ギブスサンプラー (Collapsed Gibbs sampler, CGS)
MCMC法の一種
確率モデルは
隠れ変数Z={zi}, i=1,…,N パラメータの集合{θ} 観測データX から構成される
ベイズ推定の目的は 隠れ変数とパラメータの事後分布を計算すること
必要なのは隠れ変数の事後分布
パラメータは周辺化しておき、 Zだけに注目して推論を実施
CGSを用いて実施
CGSの手順
全ての隠れ変数の値を繰り返しサンプリング
一回のサンプリングでは、 一度に全ての隠れ変数の事後分布サンプリングを行うのではなく ある一つの変数だけのサンプリングを順番に行う
例
S回目のziの事後分布サンプリングでは、 それそれ以外の変数は 全て既知のものとした分布を計算し サンプルを生成する
全てのziについて、s=1,2,…,Sまで、S回のサンプリングを行う
CGSが実現可能
上記の積分が解析的に計算可能
前提条件
パラメータと隠れ変数の確率分布が “共役(conjugate)”の関係にある
(B) CGSによるSBMの推論式の導出
具体的なSBMのCGS計算式の導出
SBMで求めたいものはZ1,Z2
CGSを用いて SBMのクラスタ割り当て 隠れ変数Z={Z1,Z2}について推論を行う
CGSの事後分布サンプリング
i番目の隠れ変数のクラスタ割り当てを一時解除
他の隠れ変数のクラスタ割り当て、 クラスタの混合パラメータπ、 クラスタ間の関係強さパラメータ を考慮して事後分布の計算
改めてクラスタ割り当てを設定する (サンプリング)
(C) 十分統計量の定義と具体的な更新式
CGSのサンプリング操作は 割り当てz1,iに伴う 十分統計量の更新操作に帰着する
プログラム上は z1,iをサンプリングするために 十分統計量を足したり引いたりするだけ
SBMの第一ドメインの隠れ変数の更新に 必要な十分統計量の更新の式
各項について、k,lは k=1,2,…,K, l=1,2,….,Lの範囲を考える
式(3.14)は 第1ドメインであるクラスタkに属しているオブジェクトの数
式(3.15)および(3.16)は、 k番目の第1ドメインクラスタとl番目の第2ドメインクラスタで定義される 関係データ行列のブロック(k,l)において x=1(0)となる関係データの数を表す
実際のz1,iのサンプリングには以上の量を使う
右式の1つ抜きの十分統計量を使って、z1,iをサンプリング
その結果に従って本来の十分統計量を更新する
実際の計算
共役性
事前分布をディリクレ分布とするパラメータと 離散分布で定義される尤度項の組み合わせならば
パラメータの事後分布は必ずディリクレ分布になる
ディリクレ分布に関する恒等式
恒等式を入れて積分計算した結果
同様に2番目の式はベータ分布とベルヌーイ分布の共役性、およびベータ分布の恒等式を用いる
Z1,i=k, k∈{1,2,…,K}となる事後分布
右式で具体的なz1,iの値をサンプリング (右式で定義される「K面サイコロ」をふることで Z1,iを割り当てる
3.4 無限関係モデル (infinite relational model,IRM)
3.4.1 IRMの概要
IRMの最大の特徴は、 潜在するクラスタ数K,Lを自動的に決定できること
クラスタ数もある種の未知パラメータとして自動的に推定する
IRMは無限のクラスタの存在を仮定した生成モデル
3.4.2 IRMの定式化
IRMの生成モデル
Z1 | α1 〜 CRP(α1)
Z2 | α2 〜 CRP(α2)
Θk,l | a0,b0 〜 Beta(a0,b0)
Xi, j | {θk,l),z1,i, z2,j 〜 Bernoulli(θz1,i,z2,j)
SBMとの違い
SBMではπをパラメータとしてとる離散分布を使って、 クラスタの割り当てz1,i,z2,jを生成
IRMではCRPという分布から直接Z1,Z2を生成する
(A) 中華料理店過程(CRP)
CRP (Chinese restaurant Process, CRP)
ノンパラメトリックベイズ(Bayesian nonparametric)の 枠組みの確率過程(stochastic process)
確率過程と通常の確率分布との違い
通常の確率分布は値の確からしさの分布
確率過程は無限個の値(すなわち分布(あるいは関数))に 対する確からしさの分布
詳細説明[69] CRPは 無限個の要素集合に対する 分割(partitioning)に対して 確率を与える確率過程
IRMに適用したケースでは、 第1ドメインのN1個のオブジェクトおよび、 第2ドメインのN2個のオブジェクトを それぞれ適当なクラスタに割り当てる (分割した結果を生成する)
サンプリングのたびにクラスタ数は確率的に変動する
CRPの説明 (中華料理店の比喩)
中華料理店の仮定
各テーブルに着席できる客数に上限はない
中華料理店に配置できるテーブル数に上限はない
入店した客をどのテーブルに配置するか (与えられたデータ要素(客)を どのようにクラスタリング(テーブルへの配席)して分割するか)
各テーブルにはそれぞれ特定の料理が配膳されており その料理の内容でテーブルを区別する
各テーブルはそれぞれが異なる分割クラスを表し インデックスkで区別する
個々の客iはデータの要素 オブジェクトや特徴ベクトルなどを表す
I番目の客がk番目のテーブルに着席したことを クラスタ割り当ての変数zi=kで表す
各テーブルに配膳される料理はクラスタkのパラメータθkを表す
N人の客が入店して、総計K個のテーブルに着席している
CRPでのn+1番目に入店する客のテーブル選択確率
集中度ハイパーパラメータαに比例する確率で 新しいテーブル(k=K+1番目のクラスタ)を生成する
データ点数に合わせて適応的に期待値クラスタ数を変化させる
すでに多くの客が着席しているテーブルに 新たな客は集まりやすい
CRPによるデータ要素のクラスタリングは、 インデックスについて可換(exchangeable)である
どのような順番で客が入店しても 同じテーブル割り当てになる確率は変わらない
各要素のインデックスiの順番によらず 最終的に得られたクラスタ割り当ての構造だけで確率が決まる
(B) πの周辺化
SBMにおけるπに当たるものはどこへ行ったのか?
CRPは「無次元のディリクレ分布と無限次元の離散分布」を組み合わせたもの
無限次元のディリクレ分布との関連
π=(π1, π2,…,πk,…)の生成過程
棒折り過程(stick breaking process) に従えば、正しいディリクレ過程 (無限次元ディリクレ分布からのサンプリングを実現できる)
棒折り過程を利用した場合の、Zの生成過程
CRPは、上記の生成モデルにおいてπを周辺化したもの
3.4.3 IRMの推論
SBMとIRMの違いは、クラスタ混合割合パラメータπに関する部分だけ
SBMのCGSでもπは周辺化してサンプリングしている
推論アルゴリズムの計算式の導出と最終的な更新則の大部分はSBMの結果を借用できる
(A) CGSによるIRMの推論式の導出
第一ドメインの第i番目のオブジェクトのクラスタ割り当てzi,iのサンプリングに関して CGSの定義式を使って変数を配置する
定義式中の変数をIRMの生成モデルから当てはめる
隠れ変数に対しては SBMと同じくZ(i) = {Z1(i), Z2}
パラメータ集合{θ}については クラスタ間の関係の強さパラメータ {θk,l} k=1,2,…, l=1,2,…だけが対象
CRPによって クラス混合パラメータπ1,π2が 周辺化されているため
式
IRMではクラスタ数がCRPによって変動
現時点のZ1(i), Z2内の異なるクラスタ数をそれぞれK,Lとする
最終的な推論式
SBMより積分が一つ少ない
(B) 十分統計量の定義と具体的な更新式
IRMの第一ドメインの隠れ変数の 更新に必要な十分統計量を定義する
数式上はSBMの場合と同じ
Z1,iをサンプリングする際には、 Z1,iは「一度席を立ってもらい」 (つまり未知量として扱う)
「一つ抜き」の十分統計量を用いる点もSBMと同じ
右の一つ抜きの十分統計量を使ってZ1,iをサンプリング (その結果に従って新しい席にz1,iを「座らせる」)
事後分布
「K+1面サイコロ」を振ることで Z1,iの座る席を改めて決定する
サンプリング結果後の、十分統計量の更新
3.4.4 出力方法
具体的にIRMの推論結果を出力する手続き
事後分布の推定が目的
具体ては名クラスタリング結果には後処理が必要
CGSの出力値は、サンプリングの実現値
CGSをS周回実施すると、 全N1+N2個のオブジェクトに対して S個ずつサンプリング結果を得る
Z1,iについては(z1,i(1),z1,i(2),…,z1,i(S))
Z2,iについては(z2,i(1),z2,i(2),…,z2,i(S))
S個のサンプルのうち、最初の適当なB期間のサンプルを全て消去
適当なM個おきのサンプル実現値だけを抽出
サンプル間の独立を担保するための手続き
T=(S-B)/M個のサンプルを用いて 「あるオブジェクトiが、あるクラスタkに南海所属したか」 をカウントする
右のカウントを正規化すると 「オブジェクトのクラスタ割り当て隠れ変数zの事後分布」 が計算できる
確率分布の式
クラスタが消滅してもカウントするため、 クラスタにには一意の投資番号を付ける
事後確率から 「オブジェクト(1,i)はあるクラスタに属する」 というクラスタ割り当てCを求める
事後分布最大化(maximum a posterior, MAP)解を利用する
複数のクラスタに対する所属確率(重み)事後分布から 最大の所属確率を持つ(重みの大きい)クラスタを一つ選択する
式
ギブスサンプラーの周回中に何らかの評価値f(Z1,Z2)を計算して 最も評価値の高かった周回śのサンプリング結果を用いる方法もある
最も簡単な方法
S回目(最後のサンプリング周回)の結果をそのまま使う
3.5 IRMのまとめ
IRMの特徴
SBMから受け継いだブロック構造のモデル化
第一ドメイン(行)と第二ドメイン(列)を複数のクラスタにそれぞれ分解(Z1,Z2)
関係データ行列を クラスタ同士の結合パラメータ(θk,l)に従って 確率的に関係データを生成する
中華料理店過程(CRP)を用いることで 任意のクラスタ構造をモデル化可能とした
上記の説明では0,1の2値関係データだが 連続値、離散値の場合も Θk,lとxi,jの確率分布を適切に変更すれば 利用可能
文献17
推論方法
本書ではギブスサンプラーを適用
推論方法の導出、実装が比較的容易
スペクトラルクラスタリングでは適用不可能な非対称関係データに適用可能
密結合、疎結合クラスタ双方とも抽出可能
IRMの入出力
第一ドメインのノード数(オブジェクト数)N1∈ℕ、 第二ドメインのノード数(オブジェクト数)N2∈ℕ
N1xN2の行列関係データX(xi,j) 関係データの値については特に制約はない 2値でも、連続値でも、整数値でも、離散値でもOK
中華料理店過程のハイパーパラメータα1,α2 経験上ではNが大きい時はα1=α2=0.1 Nが小さい時はα1= α2=1.0と設定すると うまくいくことが多い
クラスタ間接続確率パラメータのためのハイパーパラメータ (本章ではa0,b0) ベータ分布の場合は例えばa0=b0=0.5が良い
周辺化ギブスサンプラーのための定数M,B,S 典型的な値はなくNや想定されるKの大きさ 利用可能な計算機リソースと時間を鑑みて決定
IRMの出力
主産物: クラスタへの割り当て隠れ変数の事後分布
p(Z={Z1,Z2}|X) Z1,i, z2,jは角度メインにおけるオブジェクトの クラスタ帰属(cluster assignment)事後確率を表す
副産物: クラスタ割り当ての代表値C={C1,C2}
c1,i,c2,jは 各ドメインにおけるオブジェクトの代表帰属クラスタを表す
副産物: 角度メインにおける推定クラスタ総数K1,K2
クラスタ総数はZあるいはCに絵けるクラスタインデックスのことなり数とする
アルゴリズム
3.6 実データへの適用例
Zachary’s Karate Club network dataへの適用
対象データ
スペクトルクラスタリングでは2つだが、4つのクラスタを見つけている
現実の分裂結果との一致度をARIで評価
行のクラスタリングでは0.6077
列のクラスタリングでは0.6404
スペクトルクラスタリングの値よりもよい
Enron社内のEメールデータセットへの適用
非対称データ
CGSでの結果(C)ではクラスタ数は安定して(K,L)=(6,3)に収束する
一部社員に職位情報があるためクラスタ結果の検証が可能
3.7 実運用上の留意点と参考文献
3.7.1 どのアルゴリズムを使うべきか
実運用を考えるとSBMの方が好ましい
実装が楽
クラスタすうやクラスタ番号の管理に気を使わなくても良い場合は適用しやすい
SBMでクラスタ数を大きめに設定すると、 データの複雑度に見合った程度の数のクラスタのみ お餅が残り、他のクラスタは小さくなる
3.7.2 IRMの限界と拡張
SBMやIRMの生成モデルは、 一つの行(列)オブジェクトは 必ず一つのクラスに所属するという仮定を置いている
不適当な例
友達のリンクで、 「どこで知り合ったか」のパターンによって リンクの生成確率を変えることは想定されていない
上記の例に対する対応
混合メンバーシップ確率ブロックモデル (Mixed membership stochastic block model)
トピックモデルと呼ばれる混合モデルの拡張モデルを利用した関係モデル
一つのオブジェクトがつながる相手に応じて所属クラスタを選択することができる
このモデルでは有限のクラスタ数を事前に固定する必要がある
未知数のクラスタを自動的に空いている拡張手法も存在する
文献54
交友ネットワークの時間の経過による変化
時系列関係データのモデリング手法
15,25,37,86
非常に巨大なネットワークを計算するための手法
並列計算による高速化
17
計算効率とスケーラビリティを重視したネットワークモデル
20,80,92
本書では、どれかのクラスタに所属するという暗黙の仮定で検討 (網羅的かつ排他的(exhaustive and non-overlapping)クラスタリング)
そうでないデータへの対応
部分IRM(subset IRM)法
27
IRMの拡張
行および列のオブジェクトがブロック構造に従うかどうかを推定する仕組みを持つ
ブロック構造に従わないものは「その他」クラスタに割り当てられる
非網羅的かつ非排他的(non-exhaustive and overlapping, NEO)
関係データ行列のうちブロック構造を持ち かつ「値が特に周囲と異なるブロック」のみを抽出
それ以外は「その他」クラスタに分類
13,26,44
3.7.3 参考文献について
Cによる実装はあるが、可変長配列(クラスタ数の増減)を考えると 通常の混合ガウシアンモデルを下敷きにして実装した方が良い

第4章　行列分解

4.1 準備
転置
行列Aの転置をATと表記する
縦ベクトルと横ベクトル
任意のベクトルa∈ℝ𝐼は縦に並ぶもの、すなわち𝐼x1行列として扱う
向きを変えたい場合 (1x𝐼行列として扱いたい場合) 転置記号を用いてaTと表記する
インデックス集合(添時集合)
ある自然数N∈ℕについて、 1からNまでの自然数の集合を [N]={1,2,…,N}と表記する
単位行列とゼロ行列
NxNの単位行列を𝐼Nと表記する
全ての要素の値が0で与えられるNxN行列をONと表記する
インデックス
ある𝐼x𝐽行列Aが与えられた時、
i番目(𝑖∈[𝐼])の行ベクトルを ai:T=(ai1,ai2,…ai𝑗)Tと表記する
j番目(𝑗∈[𝐽])の列ベクトルを a:jT=(a1j,a2j,…aIj)Tと表記する
まとめると上記のようになる
線形代数の表記
ベクトルの内積とノルム
ベクトル𝓧,𝓨∈ℝ𝐼の内積を上式のように定義
この内積から自然に定義されるℓ2ノルムを ∥x∥ = √⟨x, x⟩と表記する
行列の内積とノルム
行列X,Y∈ℝ𝐼x𝑱の内積を上式のように定義する
この内積から自然に定義される”フロベニウスノルムの概要とアルゴリズム及び実装例“で述べているフロべニウスノルムを ∥X∥Fro = √⟨X,X⟩と表記する
アダマール積
行列A,A’∈ℝ𝐼x𝑱のアダマール積A*A’は𝐼x𝑱行列を返す関数
その(i,j)成分(I∈[𝐼], j∈[𝑱])の値はaija’ijで与えられる
ベクトルと行列の積
K∈ℕに対し、 ベクトルa, a’∈ℝ𝐼, b∈ℝ𝑱および 行列A∈ℝ𝐼x𝑱, B∈ℝ𝑱x𝑲が与えられた時 これらの積は以下のようになる
aTa’はスカラを返し、その値は⟨a, a’⟩
abTは𝑰x𝐉行列を返し、その(𝑖,𝑗)主成分(𝒊∈[𝐼], 𝙟∈[𝐽])の値はaibj
Abは𝐼次元ベクトルを返し、その𝒊∈[𝐼]番目の値はa𝒊:Tb
ABは𝙄x𝑱行列を返し、その(𝒊,𝒋)成分(𝒊∈[𝑰], 𝒋∈[𝑱])の値はa𝒊:Tb:k
ランク
行列A∈ℝ𝐼x𝑱に対し、非ゼロな固有値の数をAのランクと呼ぶ
定義よりAのランクは0以上かつmin(𝑰,𝑱)以下
微分
Θ∈ℝを引数にとる一回微分可能な関数𝒇 : ℝ→ℝに対して
導関数を∂𝒇 / ∂θ
ある点θ’∈ ℝにおける微分係数を∂𝒇 / ∂θ|θ=θ’と表記する
勾配
Θ∈ℝ𝐼を引数にとる1解微分可能な 関数𝖌 : ℝ𝑰 → ℝに対し、
ある次元𝒊 ∈ [𝐼]に関する偏導関数を∂𝔤/∂θ𝒊と表記する
全ての次元に関する偏微分をまとめたものを勾配と呼び ∇g = (∂g/∂θ1, ∂g/∂θ2,…, ∂g/∂θ𝐼)Tと表記する
ある点θ’∈ℝ𝑰における勾配の値を ∇g(θ’)= (∂g/∂θ1|θ=θ’, ∂g/∂θ2|θ=θ’, …,∂g/∂θ𝑰|θ=θ’)と表記する
4.2 単純行列分解
はじめに
関係行列データの圧縮
例:映画推薦
顧客が行、映画が列として、要素の値が評価値となるような行列
顧客数をI ∈ ℕ, 映画数をJ∈ℕとする
I x J行列Xとして表現
Xijの値が高いほど評価が高い
架空の顧客であるa∈[I]さんを考える
行列分解
顧客iのR個の嗜好をまとめたものをui:∈ℝ それを行方向にまとめたものをU∈ℝR
映画jのR個の成分をまとめたものをvj:∈ℝ それを行方向にまとめたものをV∈ℝJxRとする
xij
X≃UVT
IJ個の成分を持つXはIR個の成分を持つUと、 JR個の成分を持つVの積により表現される
Xを「X≃UVT」の形に分解すること
XのランクR行列分解 (rank-R matrix decomposition)
UとVを因子行列 (factor matrix)と呼ぶ
4.2.1 目的関数
U,VをXからどのようにして求めるのか?
「Xをうまく表現できる」 ようなU,Vを求める
XとUVTの差が小さければUVTは 「Xをうまく表現している」
E = X -UVT
Eの各要素が0に近ければ近いほど近似している
Eの大小はどのようにして測れば良いのか?
機械学習における誤差の測り方
絶対誤差
全要素の絶対値の和
Xの一部に大きな値をとる異常値があっても影響を受けにくい
2乗誤差
全要素の2乗の和
異常値の影響を受けやすい
微分が常に定義されるので、微分を用いた恋最適化の技法を使える
最終的な計算式
単純行列分解 (Simple matrix decomposition)
係数1/2は微分の係数をキャンセルするために形式的に導入されたもの
4.2.2 最適化
上記で定義された 目的関数の最適化手法
特異値分解 (Singular value decomposition)
LAPACK(Linear Algebra PACKage)など、 十分に信頼性の高い、効率的な実装が複数提供される
Xに欠損値が含まれると、一般的には適用できない
汎用的ではない
勾配法による最適化
勾配さえ計算できればどのような目的関数でも適用できる
制約の追加も比較的容易
確率的勾配法と組み合わせて計算効率を高められる
Xが大規模な場合に有効
Rが2以上の場合、解には初期値依存性があり、毎回同じ解が求まるとは限らない
4.2.3 類似度としての解釈
行列分解で得られたU,Vは どのような解釈ができるのか
Xの要素ごとに分解
xij = ui:Tvj: + eij
簡単のために、uijとvijのノルムは 1に正規化されていると仮定
全てのI∈[I], j∈[J]について ∥ui:∥=∥vj:∥= 1
ui:Tvj:は正規化したベクトル同士の内積
映画Iと顧客jの類似度(similarity)
ui:とvj:が同じ方向を向いていたら1、 反対の方向を向いていたら-1
映画iに対する顧客jの評価 = 映画iと顧客jの類似度 + 誤差
映画Iと顧客jの類似度をsim(ui:,vj:)
誤差を小さくする = xijとsim(ui:,vj:)をできるだけ近くにする
顧客iが映画jを高く評価するようであれば ui:とvj:をできるだけ同じ方向になるようにする
UとVは顧客と映画の特徴を捉えたものとなる
潜在空間
4.3 さまざまな行列分解
4.3.1 2 正則化行列分解
単純行列分解では 最小化の対象として誤差のみに着目
誤差のみを最小化すると 悪い結果を生じる場合がある
Xに大きいノイズが加わっている場合
ノイズ成分をUやVに取り込んでしまう
ノイズの影響を取り除く一つの方法
正則化(regularization)
解の範囲をせばめ、過学習を防ぐ
ℓ2正則化行列分解 (ℓ2-regularized matrix decomposition)
正則化を加えた式
λ≥0は正則化の強さを決定する係数
λ＝0だと単純行列分解
λが大きいと、誤差よりもUとVを小さくする力が大きく働く
4.3.2 非負行列分解
上記までのアプローチはXの各成分が実数全体
解析対象によっては狭いクラスに限定できる
非負行列(全てのi,jについてxij≥0)
非負の制約を加えることで 恩恵を受けられる
例:
I個の文書からなるテキストデータ
文書内に登場する単語J種類
単語には1,…,Jまでのインデックスが割り当てられる
文書内に各単語が何回出現したかを数え、 その結果をまとめてZ∈ℕIxJとする
zijには文書iに単語jが出現した頻度が格納
全文書に出現した単語の総数(重複を含む)をT=∑∑zij
これでZの各要素を割ったものをX=Z/T
Xは何を表現しているのか?
行方向
Xのi番目の行ベクトルxi:は 文書iに出現する単語の(正規化されていない)確率分布
文書と単語の同時確率を表す
例
Iが和食に関する本
「だし」や「醤油」や「卵」といった和食に関連する単語が頻出
該当するxi:の値は高くなる
I’が機械学習に関する本
「行列」や「ネットワーク」や「学習」といった技術系の単語が頻出たんごの上位を占める
概念的な意味
Rはパターンの総数
p(パターンr)はパターンrが全パターンに占める割合
単語の確率分布の式の展開
uir=p(パターンr)p(文書I|パターンr)
vjr=p(単語j|パターンr)
とすると行列分解の式になる
非負行列分解 (non-negative matrix decomposition
自然言語処理では「パターンr」のことを「トピック」
行列X(あるいはZ)からトピックょ抽出することを「トピック抽出」と呼ぶ
非負行列分解では2乗誤差最小化だけではなく、カルバック・ライブラー擬距離なども使われる
4.4 アルゴリズム
4.4.1 1次交互勾配降下法
交互勾配降下法 (Alternating gradient descent)
2つ以上の変数からなる関数の最適化法の一つ
K∈ℕ個の変数θ1,θ2,…,θkを引数にとる 関数h(θ1,θ2,…,θk)∈ℝに対し
交互勾配降下法は 以下のように変数を更新する
θknew = θk – η∇θkh(θ1new,…,θk-1new,θk,θk+1,…,θK)
η>0は学習率(learning rate) と呼ばれる最適化パラメータ
𝛈 は一回の更新において その勾配方向にどれだけ進むか、 その幅を決めるパラメータ
𝛈が小さいと変数が少しづつしか進まず 固定店にたどり着くまでに多くの反復数がかかる
大きすぎると解が発散あるいは、振動する
一般的には小さい値(例えば𝛈=0.01)から初めて 徐々に大きくしていくのが良い
ℓ2正則化行列分解の 交互勾配降下法による最適化
まず勾配を導出
フロべニウスノルムの定義より Xの各成分に関する和として書き直せる
途中の計算式1
途中の計算式2
ℓ2正則化行列分解の交互最適化アルゴリズム
変数が行列になっていてそれに関する微分を考える時
行列のまま微分をとろうとすると必要以上に難しくなる
行列石などを要素ごとの和に変換し目的関数をスカラ変数で書き直す
4.4.2 疑似2次交互勾配降下法
1次交互勾配降下法では、更新式が単純で実装も比較的容易 目的関数の1階微分(勾配)の情報しか用いないので収束が遅い場合がある
1階微分に加えて2階微分の情報を用いて収束を早くする
ヘッセ行列(Hessian matrix)の利用
ヘッセ行列を用いた逐次的な最適法 ニュートン法
ヘッセ行列をブロック対角で近似すると 逆行列を効率的に計算できる
ブロック対角行列の性質を利用した逆行列計算
𝜱と𝜳の逆行列だけを計算すれば良い
𝜱と𝜳自身もブロック行列を持つため、 さらに効率化できる
途中計算省略
擬似2次交互勾配降下法アルゴリズム (Quasi second-order alternating gradient descent method)
4.4.3 制約つきの最適化
制約を入れたい場合は交互最適化と 射影勾配降下法(projected gradient descent method) を組み合わせる
射影勾配降下法によるUの更新
ℓ2正則化つき非負行列分解の擬似2次交互射影勾配降下法アルゴリズム
4.4.4 確率勾配降下法
Xが密でかつI,Jが大きくて、 Xがメモリに乗り切らない場合
確率勾配降下法 (Stochastic gradient descent method) を用いて行列分解を計算
全サンプルを用いて勾配を厳密に計算する代わりに 少数のサンプルで勾配を近似的に計算し最適化する方法
直感的解釈
目的関数fが「サンプルnのみに依存する関数fn」の和で書ける時 (すなわち、f(x1,x2,…) = ∑ fn(xn)のように定義される時)
fに対する勾配の代わりにfnに対する勾配を使う
サンプルの取り方により 異なる確率勾配アルゴリズムが定義される
(A) 要素ごとの確率勾配法
Xの各要素xijをサンプルとする、 最も一般的な方法
要素ごとの和としてまとめると、 ui:とvj:のみに依存する関数の和として書き直せる
これに対するui:とvj:の勾配
この勾配から更新式が求まり、 これを全てのi,jについて繰り返すことでU,Vを最適化する
要素ごとの確率勾配アルゴリズム
Rを大きくして 解きたい場合は確率勾配降下方を誓うのが良い
(B) 行列ごとの確率勾配法
X全体にアクセスするのは高コストだが Xの各成分には低コストでランダムアクセスできる設定の場合に 威力を発揮する
Xの各成分のへのアクセスに時間がかかる場合に どうすれば良いのか?
例 : 単語の共起行列の行列分解
共起行列とは
複数の文書が与えられたときに 同じ文書内で同時に出現した単語の数を 数えることによって作られる行列
これを分解することで 単語の関係性を抽出することができる
M個の文書があり、全体の語彙数をI∈ℕとする
共起行列 X ∈ ℕIxIは各列と各行が単語のインデックスに対応し Xijには単語iと単語jが同時に出現した文書の数が格納されている
全文書m=1,…,Mを読み込み、 xij=∑cij(m)を計算する必要がある
IやMが大きい場合には その計算量が大きくなる
Xが複数の行列の和で与えられるとき 一つ一つの行列をサンプルとする確率勾配法を導出できる
サンプルの平均をとる作用素としてℙm[*]を導入
単純行列分解の目的関数を fstd(X;U,V) = 1/2∥X – UVT∥2Froと表記すると
C(m)を文書m∈[M]の共起行列
Xを分解する目的関数fstd(Z;U,V)は、 C(m)を分解する目的関数の和fstd(MC(m);U, V)で書ける
行列ごとの確率勾配降下法アルゴリズム
4.5 欠損値がある場合の行列分解
4.5.1 欠損値を除外
欠損値を完全に未知だとし、学習の際に一切使わない方法
勾配の式
4.5.2 欠損値を補完
欠損値を何らかの値で補完し、 目的関数に含めて計算する
Xの欠損値を補完した行列をẊとすると
補完値は前反復の値を使う
更新式は除外の場合と全く等価になる
4.6 関連する話題
4.6.1 確率モデルとしての解釈
行列分解は確率モデルとして解釈することも可能
行列分解を 確率モデルおよびそのパラメータ推定問題として 解く場合の利点
ノイズモデルの拡張
Xの各要素が実数の値をとる場合、 同じく実数の値をとるガウスノイズの過程は雨程度自然
映画推薦のような離散値を取る場合には、 ガウスノイズは適切でない
確率モデルとして考えること でXの分布を変更すること(例:正規分布からポアソン分布に変更)で この種の問題に対応可能
Xの分布を正規分布から 指数型分布族に一般化したモデル
一般化行列分解
10
事前知識の統合
木台によってはX以外にも情報が得られる場合がある
例 映画推薦の場合
顧客に対する情報(年齢や性別)
映画に対する情報(公開日やジャンル)
これらの情報は階層ベイズモデリングの考え方を用いて U,Vに対する事前分布として統合可能
U,Vのベイズ推定
事前分布を導入しベイズ推定することで U,Vの情報を分布として捉えることができる
U,Vを点として推定するのではわからなかった 推定量の確からしさなど、 より詳細な情報を得ることができる
Rの決定
尤度関数をパラメータに関して周辺化した 周辺尤度を計算することでRを決定できる
行列分解の確率的解釈の文献
57
4.6.2 非線形モデルへの拡張
本章ではXとU,Vの間に線形性を仮定
微分が容易に取れて、最適化が容易になる
非線形への対応
カーネル主成分分析
カーネル法によって行列分解を非線形化したもの
隠れ変数ガウス過程モデル
カーネル主成分分析をベイズ的に拡張したもの
Kvkをガウス過程でモデリングする方法
計算量の増大が問題
モデルの複雑度はRに依存する
4.6.3 Rの決めかた
Rの実用的な決め方
検証誤差を見る
Xの観察部分のいくつかの要素をランダムで選び検証データとする
検証データとして選ばれたデータを欠損値として学習を行う
学習後に検証データと予測値の誤差(検証誤差)を計測する
上記のステップをいくつかのRで行う
検証誤差が最小になるRを選ぶ
4.6.4 実データ
MovieLensデータ
映画推薦のデータ
中規模
Yahoo!音楽推薦データ
KDD cup 2011で使われたデータ
最大級の大きさ
KONECT
ノードすうが数百〜数千万と様々な規模の実ネットワークデータが250個
SNAP
比較的大規模なネットワークデータ

第5章　高次関係データとテンソル

はじめに
オブジェクトの種類が3種類以上の場合には、行列では軸が足りない
テンソルは行列を拡張した概念で、行と列に加えてさらなる「軸」を追加するもの
5.1 用語の定義
テンソルとして表現されるデータの例
これまでの2項から3項の関係へ
例
「顧客、映画、時間」
「経度、緯度、風速」
「主語、述語、目的語」
例:映画
I人の顧客、J個の映画、K個の時点(I,J,K∈ℕ)
評価はIxJxKテンソル(tensor)𝓧で表現できる
𝓧は3つのインデックスを決めたときに 評価の値を返す3次元配列
各要素は xijk = 時点kでの顧客iによる映画jの評価
オブジェクトの種類の数(この例では3)を テンソルでは次数(number of order) あるいは階数とよぶ
次数L次のテンソルをL次テンソルとよぶ
テンソルの各軸をモード(mode)と呼ぶ
各軸の長さを次元と呼ぶ
映画の例だ𝓧のモード「顧客」の次元はI
全ての次元をまとめたもの(IxJxK)をテンソルのサイズと呼ぶ
利便性のため、 あるモードの次元が1となるようなテンソルは、 そのモードについて省略して 一つ次数が低いテンソルとして扱う
IxJx1テンソルはIxJ行列として扱う
5.2 テンソルにおける線形演算
はじめに
ベクトル、行列、テンソルにおける操作の比較
内積に相当する「全モードを縮約」
行列やベクトルの席に相当する「モードを一つ減らす」
行列や行列席に相当する「モードの次元を変更」
外積に相当する「モードを一つ増やす」
5.2.1 加算，スカラ倍，内積
テンソルにおいて、足し算(引き算)とスカラ倍は行列と同様に定義される
IxJxKテンソル𝓧と𝓨があるとき
テンソルの足し算𝓧+𝓨は 要素ごとに𝓧と𝓨を足し合わせる操作
全ての𝒊∈[𝑰],𝑗∈[𝑱],𝒌∈[𝑲]に対して [𝓧+𝓨]𝒊𝒋𝒌=𝒙𝒊𝒋𝒌+𝐲𝒊𝒋𝒌
ある実数αに対して𝓧の要素を全てα倍する操作を スカラ倍α𝓧([α𝓧]𝒊𝙟𝒌=α𝒊𝙟𝒌
テンソルの内積
2つのテンソルの要素同士の積の総和の定義
テンソルによるフロべニウスノルムを定義
5.2.2 モード積
前準備
ベクトルや行列の積について
あるI次元ベクトルxとIxJ行列AのATx=∑xiaiを考えると
xとの積は IxJ行列を入力として J次元ベクトルを返すような作用素として 捉えることができる
上記の一般化
L∈ℕ次のテンソルを入力として (L-1)次のテンソルを返すような xとの演算を考える
ベクトルと行列の積では 「行と列のどちらに対し和を取るか」 で2通りの異なる演算が考えられる
ベクトルとテンソルの積では 「後のモードに対して和を取るか」 でL通りの異なる演算が得られる
モード𝒍∈[𝐿]に対し和を取る積を モード積(mode-l multiplicaton)と呼ぶ
I1xI2x..xILテンソルと I1次元ベクトルxのモード積A xl B = Bを 上式のように定義する
Bは(L-1)次のテンソル
行列とテンソルの積
ベクトルとテンソルの積では 特定のモードが「つぶれる」ような演算
行列とテンソルの積では 特定のモードが「伸び縮みする」ような演算
L次テンソルAとIlxJ行列Yとのモードl積
L=3の場合の式
イメージ
例
IxJ行列A、IxK行列B、JxN行列Cに対し
A x1 B x2 C =BTAC
5.2.3 スライス
行列の場合
IxJ行列Xに対し
I∈[I]番目の行を止めて列方向に インデックスを動かすことで 行ベクトルxi:が得られる
j∈[J]番目の行を止めて列方向に インデックスを動かすことで 行ベクトルx:jが得られる
テンソルの場合
i番目のモード𝒍スライス
𝒳i(l)
I1xI2x…,ILテンソル𝓧のモード𝒍をi∈{I𝒍]番目で固定し その他のインデックスを動かすことで定義される I1xI2x…xI𝒍-1xI𝒍+1x…ILテンソル
テンソルのスライスイメージ
5.2.4 外積
ベクトルの場合
ベクトルxとyの外積xyTは、一次のテンソル(ベクトル)から 2次のテンソル(行列)を作り出す操作と言える
テンソルでの一般化
L次のテンソルからL+1次のテンソルを作り出す
任意のテンソル𝓧と任意のベクトルyに対して 𝓧x(𝓧の次数+1)yを𝓧⚪︎𝓨と書く
行列とベクトルの外積のイメージ
5.3 行列演算への変換
はじめに
テンソル演算の多くは、 ベクトル化作用素とクロネッカー積を使うことで 等価な行列演算として書き直すことができる
書き直された行列演算は BLAS(Basic Linear Algebra Subprograms)などの 標準の線形演算ライブラリが利用できる
5.3.1 ベクトル化作用素とクロネッカー積
準備としてテンソルの ベクトル化作用素vecを導入する
テンソルの要素を番号が小さいモードから ベクトルに並べ替える作用素
例
IxJxKテンソル𝓧に対して vsc(𝓧)はIJK次元のベクトルを返す
[vec(𝓧)]i+Ij+Ijk = xijk
クロネッカー積
A⊗B
5.3.2 モード積，外積，内積
モード積
I1xI2x…xILテンソルと Il x Il’行列Al(l∈L)のモード積
外積で定義されるテンソル(すなわちランク1テンソル)を ベクトル化したものはクロネッカー積で書くことができる
テンソル同士の内積は、 定義よりベクトル化したもの同士の 内積として書くことができる
5.3.3 計算量の注意
モード積と行列演算化の違い
素直なテンソル演算はテンソルを「膨らませる」作業と「圧縮する」作業を同時並行的に行う
行列演算化は「膨らませる」作業をまとめて行い(クロネッカー積)、 その後に「圧縮する」素行を行う(ベクトルとの積)
中間処理に必要な変数のきぼが大きくなる
実際には線形演算ライブラリの恩恵を受けられる
5.3.4 テンソル演算の式展開のコツ

第6章　テンソル分解

6.1 テンソルの次元圧縮
テンソルで次元圧縮を行うにはどうすれば良いか?
例:3次元での検討 データとしてIxJxKの大きさを持つ テンソル𝓧を考える
6.2 CP分解
はじめに
与えられたテンソルをランク1テンソルの輪に分解する方法
6.2.1 動機
6.2.2 目的関数
6.2.3 アルゴリズム
6.3 タッカー分解
はじめに
CP分解に比べてより複雑な表現が可能な方法
6.3.1 動機
6.3.2 目的関数
6.3.3 アルゴリズム
6.4 CP分解とタッカー分解の違い
はじめに
CP分解はタッカー分解の特殊ケース
6.4.1 類似度としての解釈
使い分け
ある程度独立であればCP分解
相互にいづんするのであればたかー分解
6.5 補足
6.5.1 実応用例
国同士の関係性を抽出するためにタッカー分解を拡張した方法を提案
国x国x時点x行動
72
言語データへの CP分解の効率的な計算による応用
主語x動詞x目的語
32
脳情報データの解析
脳波 (Electroencephalogram EEG)
チャネルx周波数x時間
どの脳波の部位がどのようなタスクに関連して活発化するか
6.5.2 さまざまな最適化と学習アルゴリズム
確率降下勾配法
凸最適化
6.5.3 2 次の交互作用のみからなる分解
6.5.4 複数のテンソルや行列が与えられた場合
6.5.5 その他のテンソル分解