クラスタリングから系統樹

(1)

配列解析アルゴリズム特論渋谷配列解析アルゴリズム特論渋谷

クラスタリングと系統樹

渋谷

東京大学医科学研究所ヒトゲノム解析センター（兼）情報理工学系研究科コンピュータ科学専攻 http://www.hgc.jp/~tshibuya

(2)

配列解析アルゴリズム特論渋谷

今日の話題

¦

クラスタリング・アルゴリズム

¿

K-meansのバリエーション

uK-medoids uCLARA

¿

SVD

¿

probabilistic clustering

¿

DBSCAN

¿

階層的クラスタリング

¿

系統樹関連アルゴリズム

(3)

学習アルゴリズムの分類

¦

教師無し学習

¿

正解セットも何もなく、何らかの情報を得る

uクラスタリングなど（今回）

¦

教師有り学習

¿

学習用の正解セットをもとに、学習を行う

uナイーブベイズ、ＨＭＭ、ＳＶＭなど u実は圧縮も「教師有り学習」のひとつ

(4)

クラスタリング

¦ データを類似しているものごとにグループ化すること ¿ データ u 数値的な特徴ベクトル u 非数値的な特徴ラベル ’ 生物の表現型など ¿ 類似度の定義はいろいろ u ミンコフスキー距離 ’ {Σ_i |x_i-y_i|d}1/d ’ ユークリッド距離(d=2)、マンハッタン距離(d=1) u ベクトル同士の角度 ’ あるいはそのcosine u 二次形式距離 ’ ( (x−y)T M (x−y) )1/2 » 各波長の強度からなる高次元ベクトルで表された色の間の距離など » M=I ならユークリッド距離と同じ u ハミング距離、編集距離、アラインメントスコア ’ 文字列間 u Ｊａｃｃａｒｄ(タニモト)距離 1−#(X∩Y)/#(X∪Y) ’ （ラベルなどの）集合間 ’ _{1 − (𝑥𝑥 � 𝑦𝑦)/(|𝑥𝑥|}2_{+ |𝑦𝑦|}2_{− 𝑥𝑥 � 𝑦𝑦)} 一般ベクトル拡張版 u 相関係数 ’ 時系列の変化量など

(5)

クラスタリングの種類

¦

クラスタリングの種類

¿

非階層型クラスタリング

uいきなりいくつかのクラスタに分割

¿

階層型クラスタリング

uボトムアップ ’ 局所的なグルーピングから順番に uトップダウン ’ 全体の分割から順番に

(6)

最遠点クラスタリング（１）

¦

各クラスタの直径(最も遠い2点間の距離）の最大を最小化する

¿ クラスタ数 k はgiven

¦

そのバリエーション

¿ 各クラスタを含む最小の円（球）の半径のうち最大のものを最小化する u ユークリッド空間でなく、球を考えにくい場合は代表要素からの最大距離

¦

しかし残念ながら、2次元以上の数値データではNP完全

(7)

最遠点クラスタリング（２）

¦

2倍近似アルゴリズム

1. 適当な要素を選んで要素集合 P に入れる 2. P の要素からの距離（一番近いもの）が、もっとも大きい要素を選ん で P に入れる 3. |P| が k になるまで繰り返し 4. P の各要素を代表点として、P のどの要素に近いかによって全要素 を分類

¦

ちなみにこのアルゴリズムにおける最大直径は最適解の2倍

以内

¿ 半径も最適最大半径最小化の2倍以内 ¿ 2倍というのはあくまで最悪比で、実際にはもっとよい値をとることが多い

¦

実は、1.969近似以下の最遠点クラスタリングはNP完全

¿ 同様に1.822近似以下の最大半径最小化クラスタリングもNP完全 ¿ ユークリッド空間でない一般距離の場合は2近似未満がNP完全

(8)

最遠点クラスタリング（３）

¦

2倍であることの証明

¿ 下のPから作成したクラスタの各直径は2d 以下 u半径はd 以下 ud より距離の離れた点はもう残っていないため ¿ k+1個だけを考えた最適クラスタリングの直径はd 以上 uどの2点間をとっても d 以上の長さであるため、（１つだけある）大きさ２ のクラスタの直径もd 以上 ¿ したがって、全体の最適クラスタリングの直径も当然 d 以上 P : グリーディーに 最初に選んだk個 k+1個目 d どのクラスタをとっても半径d以下 k+1個目を含むクラスタは直径d以上

(9)

分散和最小化クラスタリング（１）

¦

最遠点はあまり密度分布を考慮していない

¦

かわりにクラスタ内の分散の和を最小化したい

¿

(X-C(X))

2

の平均を最小化することを考える

uC(X) は X を含むクラスタの重心 u分散の定義: V(X) = E((Ｘ-E(X))2)

¿

ただし

（あるグループ内の各点間距離の二乗の平均） ∝ （そのグループの重心からの距離の二乗の平均）

が成立するため、点間の距離が計算できれば、必ずしも

重心が計算可能である必要はない

1 𝑛𝑛 �_𝑖𝑖 𝑛𝑛 (𝑣𝑣_𝑖𝑖−𝐺𝐺)2 = _𝑛𝑛1₂ (𝑛𝑛 � 𝑖𝑖 𝑛𝑛 𝑣𝑣_𝑖𝑖 2 − (� 𝑖𝑖 𝑛𝑛 𝑣𝑣_𝑖𝑖)2) = _2𝑛𝑛1₂(∑_{𝑖𝑖,𝑗𝑗}(𝑣𝑣_𝑖𝑖 − 𝑣𝑣_𝑗𝑗)2)

(10)

分散和最小化クラスタリング（２）

¦

分散和が最小であるための必

要条件

¿

各クラスタの重心を考えたとき、

ある要素から各重心の中で一番

近い要素は、その要素を含むク

ラスタの重心である

¦

しかし、、、

¿

クラスタを求めないと重心は定ま

らない。

¿

重心がさだまらないとクラスタが

求められない。

¿

しかも十分条件、というわけでも

ない

(11)

ボロノイ図

¦

分散和最小化クラスタリングは、各クラスタの重心

によって作られるボロノイ図的な分割になっている

¦

ボロノイ図

¿

与えられた点集合のどの点に近いかによって空間を分割

した図

(12)

データが一次元であれば、

¦

DPで簡単に解ける

1. データを数値でソート 2. Cluster(i, k)  i 番目までの数値に対するクラスタ数kの最適クラスタリング結果  Cluster(j, k−1) (k−2 < j < i) のすべての値があれば計算可能

¦

O(kn

3

)

¿ ただし、頑張ればO(n·k)あるいはO(n log n)にできる

uMonge性(既出）を利用 [Aggarwal et al, '87, '88, '94]

¦

しかし、2次元以上ではNP完全

¿ 定数近似アルゴリズムも存在しない a _b c d ここでのMonge 性： a+b ≤ c+d

(13)

定数次元で、

k

=定数であれば多項式時間

¦

2次元上の2-means

¿

直線による２次元上の点集合の分け方すべての中から、分

散が最小のものを探す方法

¿

直線による集合の分け方の数

u高々O(n2)通りしかない（＝O(２点を通る直線の数)） u全部調べても多項式時間

¦

d=const, k=const ならやはり（巨大だが）多項式時間

この「間」のいかなる直線による「集合の分け方」も、左右どちらかの２点をとおる直線（ただし２点の「ぎりぎり」手前）で同じ分け方が可能

(14)

k

-Means法

1. 適当にクラスタを作る

2. 各クラスタの重心を求める

3. どの重心に近いかによって、それぞれの要素を各クラスタに分

けなおす。

¿ 必ず分散和は改善される！ → （局所最適解に）収束する u このアルゴリズムはEMアルゴリズムの一種と見ることもできる

4. クラスタに変化がなくなれば終わり。そうでなければ２へ

各重心のボロノイ図で分割 [Forgy '65]

(15)

クラスタ数の決定

¦

明確な指標はない

¦

１つの方法としては、グループ内での各点間距離の

二乗の総和の変化で判断

¿

そもそもこれを最適化したかった

1 2 3 4 5 6 グループ内距離の二乗総和クラスタ数総和ががくんと落ちた直後は、「良い」クラスタ数であると考えられる

(16)

ちなみに、凝集性と分離性

¦

凝集性（cohesion）

¿

クラスタ内の各点間距離の二乗平均

¦

分離性（separation）

¿

クラスタ間の各点間距離の二乗平均

uクラスタが３つ以上の場合、各要素に対しその要素に一番近いクラスタのみ考える

¦

シルエット係数(siｌhouette coefficient)

¿

１－（凝集性／分離性）

¿

ひとつひとつの要素についても、ひとつひとつのクラスタに

ついても、全体についても計算可能

(17)

バリエーション１：

k

-Medoid法 (PAM: partition around medoids)

¦

重心のかわりにMedoidを用いる

¿

クラスタを代表する要素

uクラスタ中の他の要素との距離の和が最小の要素 uクラスタ中の他の要素との共通ラベル数が最大の要素（k-mode法） uetc.

¦

特徴

¿

数値以外のデータにも適用が可能

¿

outlier の影響を小さくすることができる

¿

ランダムネスを導入することも可能

(18)

バリエーション２：

CLARA (Clustering Large Applications)

¦

k-medoidの高速化

1. 2k 個程度の要素をランダムに選ぶ

2. k-medoidを計算する

3. 得られた代表データを各クラスタの代表点だと考え

て、全体をクラスタリングを行う

4. 以上を異なるランダム要素の集合から何回か繰り

返して一番よさそうなものを選ぶ



あまりよくない解を避けるため

(19)

バリエーション３：

BIRCH

(Balanced Interactive Reducing and Clustering using Hierarchies)

¦ k-Meansを階層的クラスタリングに拡張 ¦ CF (cluster feature)木 ¿ 子供の数は高々 B ¿ 葉がクラスタを表す u 要素の分散平均が T 以下 u クラスタの要素数は L 以下 ¿ 各ノードでは、その子孫の葉に属する要素の重心と分散和を管理する ¦ 作成方法 ¿ 最初はルートのみ ¿ 要素をひとつずつ加える u 重心の近い子供を辿っていく u 葉のクラスタに入れても分散平均がT以下、要素数L以下のままであれば、そのまま入れる u そうでない場合は、２つに分割する ’ 通常の2分割のクラスタリング(k-meansなど)を行う ’ 親の次数がBを超えたら親も分割 » 重心に対するクラスタリング(k-meansなど) ’ 全体をバランスさせる必要があれば、バランスさせる

(20)

バリエーション４：

Fuzzy C-Means

¦

データが複数のクラスタに属することを許した

k-Means

¿

それぞれの要素は、すべてのクラスターに距離に反

比例した比率で帰属するものとする

u

ある要素から３つの重心までの距離が１，2，2であれば、

それぞれのクラスタに0.5, 0.25, 0.25の比率で属しているも

のとする

¿

重心は重みをつけて計算

¦

こういうクラスタリングのことをソフトクラスタリン

グという

¿

通常のクラスタリングはハードクラスタリング

(21)

バリエーション５：

SOM (Self-Organizing Maps)

¦

高次元空間→低次元空間のマッピング + クラスタリング

¿ 同じ次元の空間へのマッピングでもOK ¿ もともとはニューラルネットワークのお話 u ニューロンがどのようにして空間を認識するのか、のモデル（コホーネン・ネットワーク） ’ ニューロンに刺激があると、「近く」のニューロンも刺激を受けてしまう

¦

アルゴリズム

1. 写像先空間に適当な大きさの格子上の格子点（ユニット）を置く 2. 写像先空間の各格子点にそれぞれ対応する代表点（プロトタイプ）を原空間に置く 3. 各データ x_i に対し、最も近いプロトタイプ p_iを探し、それに対応するユニットu_iの「近傍」にあるユニットの集合をU_iとする。 U_iの各ユニットに対応するプロトタイプそれぞれをx_iの方向へ動かす  p_j= p_j+ α∙w(|u_j − u_i|) ∙(x_i − p_j)  w: 対応するユニットがu_iに近いほど大きく  近傍が自分自身のみであれば、k-meansと近い挙動をする 4. αをだんだん小さくしながら3を繰り返し、あまりプロトタイプが動かなくなるか、既定回数以上計算をしたならば終了各データ x_i に最も近いプロトタイプ p_iに対応するユニットが、そのデータの新たな空間での座標となる（近いデータは近い場所にマップされる）

(22)

特異値分解（

SVD）による次元削減（１）

¦ 特徴数が膨大な時など、特徴ベクトルの次元数がきわめて大きくなる ¿ すると計算時間がきわめて大きくなる ¦ 次元を減らしたい ¿ 特徴ベクトルのrankが r ならば、r 次元に減らせる ¿ 特異値分解をする u A= USVT u B=UT_{A = SV}T _{の列ベクトルの} _{1～r 番目の要素が求める r 次元ベクトル集合} 特徴ベクトル１＝ _× n m m n _× ・・・ U S VT σ₁ σ₂ σ₃ ・・・ O O 正規直交特異値：σ₁≥ σ₂≥ …≥ 0を満たす m m m m Aのランクがｒであれば、i > r の場合 σ_i= 0 A r r 特徴ベクトル２・・・正規直交基底ベクトル１基底ベクトル２・・・変換ベクトル１変換ベクトル２意味のある基底ベクトルはr個のみ r 非負対角

(23)

特異値分解（

SVD）による次元削減（２）

¦

その意味

¿

特徴ベクトルが

r 次元に乗っている

特徴ベクトルは（r 個の基底ベ クトルで張られた）この r 次元 空間内に存在する

(24)

特異値分解（

SVD）による次元削減（３）

¦

ヒューリスティックでよければ、さらに次元は減らせる

¿

特異値が小さい次元は無視してもよい

u特異値が大きければ大きいほど、影響度が大きい u特異値が大きいものから k 個選べば、 k 次元空間に落とすことが可能 特徴ベクトル１＝ _× n m m n _× ・・・ U S VT σ₁ σ₂ σ₃ ・・・ O O m m m m A r r 特徴ベクトル２・・・基底ベクトル１基底ベクトル２・・・変換ベクトル１変換ベクトル２ r k 次元分のみを用いる k k k

(25)

特異値分解（

SVD）による次元削減（４）

¥

こうして

k次元超平面に落としたものは、k次元回帰超平面に

なっている

¿

各データの超平面との距離の二乗誤差の和が最小

¦

高次元のベクトルを低次元で可視化することができる

¦

妥当な次数

¿

∑

_i=1..k

σ

_i2

/ ∑

_i=1..r

σ

_i2

（すなわち行列のフロベニウス・ノルムの

比率）が適当な値を超えるように設定

(26)

SVDによるクラスタリング

¦

直接的には、

¿

次元を落とした空間上で、通常のクラスタリングを行う

u計算時間の削減 uノイズの削減にも効果

¦

あるいは

¿

SVDによって得られた（特に特異値の大きな）基底ベクトル

は、データの代表的ベクトルだと考えることができる

uそれらは直交したベクトルで、互いにほとんど似ていない！

¿

それらの基底ベクトルとの類似度によって、クラスタリングを

行うことができる

uハードクラスタリングとすることも、ソフトクラスタリングとすることもできる uあるいはそれを初期解として、別のアルゴリズムを走らせることもできる

(27)

MDS (Multi-Dimensional Scaling)

¦

（特徴ベクトルではなく）要素間の距離が与えられてい

る場合に低次元にマッピングする方法

¿

cf. カーネル法（SVM)

¦

以下の値を最小化する

¿

Σ(d

_ij

− D

_ij

)

2

/ Σ(D

_ij

)

2 uD_ij: もとの距離 d_ij: マップされた距離

¿

一般的な場合、この最適化は難しい

uヒューリスティックで解く ’ 最初は適当に配置して改善していく、など

¿

さまざまなバリエーション

¦

SVDの時と同様に、この後通常の（低次元用）クラスタ

リングをすればクラスタリングができる

(28)

確率的クラスタリング（１）

¦

正規分布である等、数値の分布がわかっている場合

にはもっと統計的にきちんとした扱いをすることが好

ましい

¦

混合分布モデル

¿

複数のモデルに基づいてデータが存在しているものと仮定

¿

それぞれのモデルのパラメータを推定したい

1次元上の分布 2次元上の分布

(29)

確率的クラスタリング（２）

¦ EMアルゴリズム（確率パラメータの一般的学習法）をそのまま適用する ¿ ただし、ＥＭはあくまでヒューリスティックによる局所解を得る解法なので過信は禁物 ¦ パラメータを推測できたら、後は属する確率が最大であるクラスタに各要素を分配するだけ ¿ 確率をそのまま出力すれば、ソフトクラスタリング

¦

クラスタ数も学習の対象とできる

¿ 情報量規準 A B C この要素がAに属する確率が50%、Bに属する確率が４５％、Cに属する確率が５％、ならば、それぞれ0.5個、0.45個、 0.05個のA,B,Cのデータがあると考えて、もう一度A,B,Cのパラメータを推測する。

(30)

EM Algorithm (1)

¦

データが一部しかない場合の最尤推定アルゴリズム

¦

やりたい最尤推定

¿

x: 観測データ

y: 観測できないデータ

¿

θ: 知りたいパラメータ

¿

P(x|θ) = Σ

_y

P(x,y|θ) を最大化するθを求めたい

¦

EM アルゴリズム

¿

適当な初期解を与えて、それを改善して行くアルゴリズム

u必ず改善することと、何らかの解に収束することが証明されている uただし、厳密な最適解を得るアルゴリズムではなく、あくまで局所最適解を求めるヒューリスティックにすぎない ’ 厳密解は（多くの場合）極めて難しい

(31)

EM Algorithm (2)

¦

例題

¿

Input: 男女あわせて

n 人の身長のデータ

¿

Output: 男女比率および男女それぞれの平均身長

¦

本来はこんな問題をとける筈はないが、適当な仮定

のもと、推測することはできる（かもしれない）

¦

仮定

¿

男性の身長の平均＞女性の身長の平均

¿

それぞれ正規分布に従う

¿

尤度がより大きいパラメータがより正しいパラメータである

男性女性身長分布

(32)

EM Algorithm (3)

¦

アルゴリズム

¿ まず適当に男女比、男女それぞれの平均・分散を仮定する。 uある程度正確な初期値であることは必要 ¿ そのモデルに基づき、各データがどれくらいの確率で男性であるか、あるいは女性であるかを予測する u男性 r , 女性 1−r ¿ 男性である確率が r であるデータを男性 r 人、女性 1−r 人分のデー タとして扱って、男女比・平均・分散を計算する ¿ 収束するまでこれを繰り返す男性女性身長分布男女比と男性・女性それぞれの平均・分散がわかっていれば、この身長の持ち主が男性である確率は85%、というような推定が可能

(33)

EM Algorithm (4)

¦

直感でなく式で書くと、、、

¿

θ

t+1

：

Q(θ|θ

t

)=E

_y

(log P(x,y|θ)|

x, θ

t

) を最大化する θ

ux: 観測データ、y: 潜在データ（観測できない）、 θ ：パラメータ ’ 先の例では x が身長データ、y が各データの男女別データに相当

¿

この時、

0 ))

,

|

(

),

,

|

(

)

|

(

)

|

(

)

|

(

log

)

|

(

log

≥

+

−

=

−

t t t t t

x

y

P

x

y

P

KL

Q

x

P

x

P

θ

∑

− = ≥ ⋅ = = y t t t y t t y t t y x y P x y P Q x y P y x P x y P x y P y x P x y P y x P x P )} , | ( log ) , | ( { ) | ( )}] , | ( / ) | , ( log{ ) , | ( [ )} , | ( / ) | , ( ) , | ( { log )} | , ( log{ )} | ( log{( θ θ θ θ θ θ θ θ θ θ θ θ ちゃんと書くと Jensen's inequation x x これを大きくしたいこれを大きくすればよい（するとP(x|θ)の値が必ずよくなる） x でも、この項（正の値）もあるので、選んだθが最適なわけではない

(34)

Density-based Clustering

(35)

CLIQUE (Clustering In Quest)

¦

まず（

d次元）格子で特徴空間を分割

¦

密度の高い部分をつなげる

(36)

DBSCAN

(Density-based spatial clustering of applications with noise)

¦

要素xの近傍

¿

xからの距離がε以内（εは適当）

¦

密度の高い近傍

¿

近傍内に

p個以上の要素がある（pは適当）

¦

近傍関係のグラフ

¿

yがxの近傍に、yがxの近傍に入っていれば、入っていれ

ば枝を張る

¦

そのグラフの連結成分分解をしたものが出力

p=3

(37)

階層的クラスタリング

¦

トップダウン型は、非階層的クラスタリングを繰り返し

行うことで可能

¦

ボトムアップ型

¿

距離の近いものをまとめて、クラスタを作ることを繰り返す

¿

クラスタ間、要素-クラスタ間の距離の定義が問題

¿

途中で止めることで、非階層クラスタリング結果を得ること

も可能

(38)

クラスタ間の距離

¦

最短距離 (Single linkage)

¿ D (C₁, C₂) = min { d (x,y) | x∈C₁, y ∈C₂ } ¿ ノイズに非常に弱い

¦

最長距離 (Complete linkage)

¿ D (C₁, C₂) = max { d (x,y) | x∈C₁, y ∈C₂ } ¿ ノイズに非常に弱い

¦

平均距離 (UPGMA, Average linkage）

¿ D (C₁, C₂) = <d (x,y) | x∈C₁, y ∈C₂ > （平均）

¦

代表値間の距離

¿ 中央値、重心、etc.

¦

Ward法

¿ D(C₁, C₂) = E (C₁∪C₂) − E (C₁) − E (C₂) u_{E (C) = ∑}_x∈C{d (x, c)}2 _（_{c は Cの重心）} ’ クラスタ内の二乗誤差 u 分離性を考慮にいれたクラスタ内分散

(39)

Single / Completeの意味

¦

Unit-diskグラフ

¿

距離がある閾値以下ならば要素と要素を枝で結ぶ

uちなみに構造生物学ではこれのことをContact map graphとよぶ

¦

ある閾値で切った時の（非階層）クラスタリング結果

¿

最短距離(single link)クラスタリング

u連結成分分解

¿

最遠距離(complete link)クラスタリング

u極大完全部分グラフへの分解

¦

グラフ理論を用いたクラスタリングアルゴリズムは多い

¿

CLICK(Cluster Identification via Connectivity Kernels)

uMinimum weight cut でグラフを分割

’ 分割した時に切られる枝の重みを最小化

(40)

系統樹とは

¦

系統樹とは

¿

生物的な対象を「進化」の視点から階層

的クラスタリングを行ったもの

¦

系統樹の種類

¿

rooted

¿

unrooted

¿

network (recombinationを考慮)

¦

系統樹のアルゴリズム

¿

Distance-based methods

u距離が与えられている場合

¿

Maximum-Parsimony methods

uラベルから ’ 生物の「性質」など

(41)

レポート課題

¦

以下の中から一題選び、レポートを提出せよ。

1. SETHの元で計算量の下限が示されている問題について、その下限の証明を示すとともに、その問題に対して知られている最良のアルゴリズムを説明せよ。 2. 講義では扱っていないWavelet treeの応用を紹介し、それについて論じよ。 3. Tandem Repeatを求める線形時間アルゴリズムについてわかりやすく説明せよ。 4. 以下より１つ選び、実装・実験・考察せよ。 a. ψ配列による圧縮接尾辞配列を実装し、通常の接尾辞配列と性能比較せよ（作成アルゴリズムは問わない）。 b. FFTを用いたパタンマッチングアルゴリズムと、そうでないものの性能を比較せよ（対象は文字列でも画像でも時系列データでもなんでも） c. Tandem Repeatを求めるアルゴリズムを2種類（以上）実装し、性能を比較せよ（対象は文字列でもDNAでも01ビット列でもなんでも）。

締め切り：

8月5日（木）

提出先：

ITC-LMSへ