カーネル法のバイオインフォマティクスへの応用

配列および化合物データ解析
のためのカーネル法
阿久津達也
京都大学 化学研究所
バイオインフォマティクスセンター
サポートベクターマシン




カーネル法の一つ、ニューラルネットワークと類似
1990年代に、Cortes と Vapnik が発明
トレーニングデータとして与えられた正例と負例から、
それらを分離する超平面を計算
機械学習、統計学、人工知能、パターン認識、バイオ
インフォマティクスなど様々な分野に応用






配列分類
タンパク質フォールド予測、二次構造構造
遺伝子発現データ解析
タンパク質相互作用予測
化合物の性質推定
c.f. Kernel Methods in Computational Biology, MIT Press,
2004
サポートベクターマシン


正例と負例を与
えて、それらを最
適(マージンを最
大)に分離する超
平面を学習
カーネルを適切に
定義することによ
り超平面以外で
の分離が可能
テストデータ
margin
SVMによるテストデータの分類


学習データより超平面
を学習(SVM)
テストデータは、対応す
る点の超平面に対する
位置(上下)で判定

テストデータ
テストデータとサポート
ベクター間のカーネル関
数値の重み付き和でテ
ストデータを類別
margin
f (x ) 
  K (x,x )  
i:xi  
i
i
j: x j  
j
K (x , x j )
カーネル




サポートベクターマシン:基本的には超平面で分離
Φ(x) (特徴ベクトル):「非線形曲面⇒超平面」に写像
カーネル K(x,y)=φ(x)・φ(y)
x と y の類似度が高い ⇔ K(x,y)が大
φ(x)
カーネルの例




線形カーネル:
K(x,y) = x・y
多項式カーネル: K(x,y) = (x・y + c)d
RBFカーネル: K(x,y) = exp (-||x - y||2 /2σ2 )
シグモイドカーネル(厳密にはカーネルではない):
K(x,y) = tanh (κx・y - δ)
カーネルとなるための条件


カーネルの定義: K(x,y)=φ(x)・φ(y)
Mercer条件を満たす ⇒ カーネル

連続値の場合
K
(
x
,
z
)
f
(
x
)
f
(
z
)
d
x
d
z

0


離散値の場合 ( x1,x2,…,xn が入力データ)
K  (K (x i , x j ))
n
i , j 1
が半正定値行列
カーネルの作り方
データから特徴ベクトル(feature vector)を
作るのが一般的、かつ、
多くの場合に実用的
 特徴ベクトル: 実数値の列
 例えば、各化合物 x に対し、


Φ(x) = (分子量, 容積, 表面積, logP,…)
とすれば、化合物 x,y に対するカーネルは
Φ(x) と Φ(y) の単なる内積
アライメントカーネルによる構造予測
1.
2.
3.
4.
5.
6.
SCOPとスーパーファミリー予測
既存カーネル
配列解析手法(アライメント、HMM)
新カーネル
計算機実験結果
結論と課題
タンパク質立体構造予測



アミノ酸配列から、
タンパク質の立体
構造(3次元構造)
をコンピュータによ
り推定
実験よりは、精度が
悪い
だいたいの形がわ
かれば良いのであ
れば、4~5割の予
測率
タンパク質配列(アミノ酸配列)
T C A V F G L G G V R L S D
V
コンピュータ
タンパク質
立体構造
フォールド予測 (Fold Recognition)
タンパク質配列(アミノ酸配列)


精密な3次元構造
ではなく、だいたい
の形(fold)を予測
立体構造は数千
種類程度の形に分
類される、との予
測(Chotia, 1992)
に基づく
T C A V F G L G G V R L S D
V
1000個のテンプレート構造
SCOPデータベース

タンパク質立体構造を形状を中心に、人手で、
階層的に、分類したデータベース
SCOP
Root
‥‥‥‥‥
Class.1 Class.2
‥‥‥‥‥
Fold.1 Fold.2
Super
Super
Family.1 Family.2
Family.1
‥‥‥‥‥
Family.2
mkkrltitlsesvlenlekmaremglsksam
ispqarafleevfrrkqslnskekeevakkcg
isvalenykkgq
itplqvrvwfinkrmrs
Family.3
Super Family 予測

入力配列が SCOP のどのスーパーファミリー
に属するかを予測
Super Family.1
タンパク質配列
madqlteeqiaefkeafslfdkdgdgtittkel
gtvmrslgqnpteaelqdminevdadg
Super Family.2
ngtidfpefltmmark
:
:
Super Family.3
既存手法の主なターゲット
Class
Secondary Structure Prediction
Fold
Threading
Super
Family
HMM, PSI-BLAST, SVM
Family
SW, BLAST, FASTA
タンパク質配列解析のための既存カーネル

HMMから特徴ベクトルを抽出



配列から直接特徴ベクトルを抽出



Fisher カーネル (Jaakkola et al., 2000)
Marginalized カーネル (Tsuda et al., 2002)
Spectrum カーネル (Leslie et al., 2002)
Mismatch カーネル (Leslie et al., 2003)
他の配列とのスコアを特徴ベクトルとして利用

SVM pairwise (Liao & Noble, 2002)
Spectrumカーネル
AA AC AG
A CCT A
φ(x)
( 0
1
0
CC CG CT
1
0
1
TA
1
)
配列アライメント




バイオインフォマティクス
の最重要技術の一つ
2個もしくは3個以上の配
列の類似性判定に利用
文字間の最適な対応関
係を求める(最適化問題)
配列長を同じにするよう
に、ギャップ記号(挿入、
欠失に対応)を挿入
A L G F G S L Y G
A L G G V S V G
A L G F G
A L G
S L Y G
G V S V
G
スコア行列(置換行列)

残基間(アミノ酸文字間)の類似性を表す行列

PAM250, BLOSUM45 など
A
A
R
N
D
C
Q
E
G
H
I
L
K
M
F
P
S
T
W
Y
V
5
-2
-1
-2
-1
-1
T
W
Y
V
-2 -1 -2 -1 -1 -1 0 -2 -1 -2 -1 -1 -3 -1 1 0
7 -1 -2 -4 1 0 -3 0 -4 -3 3 -2 -3 -3 -1 -1
-1 7 2 -2 0 0 0 1 -3 -4 0 -2 -4 -2 1 0
-2 2 8 -4 0 2 -1 -1 -4 -4 -1 -4 -5 -1 0 -1
-4 -2 -4 13 -3 -3 -3 -3 -2 -2 -3 -2 -2 -4 -1 -1
1 0 0 -3 7 2 -2 1 -3 -2 2 0 -4 -1 0 -1
R
N
D
C
Q
E
G
H
I
L
K
M
F
P
-3
-3
-4
-5
-5
-1
-2
-1
-2
-3
-3
-1
0
3
-3
-4
-1
-3
BLOSUM50 スコア行列
(置換行列)の一部分
S
ペアワイズ・アライメント


配列が2個の場合でも可能なアライメントの個数は
指数オーダー
しかし、スコア最大となるアライメント(最適アライメ
ント)は動的計画法により、O(mn)時間で計算可能
(m,n:入力配列の長さ)
入力配列
AGCT, ACGCT
アライメント
AGCT ACGCT
スコア -3
AG - CT
ACGCT
1
最適アラ
イメント
A - GCT
ACGCT
3
- AGC - - T
AC - - GCT
-5
(同じ文字の時: 1、違う文字の時: -1、ギャップ1文字: -1)
動的計画法による大域アライメント(1)
(Needleman-Wunschアルゴリズム)



入力文字列から格子状グラフを構成
アライメントと左上から右下へのパスが一対一対応
最長経路=最適アライメント
G
G
F
V
D
5
-5
K
-2
-5
Y
-5
7
D
1
-6
-1
-2
-3
-2
1
0
-4
4
-7
-7
-7
-7
アライメント
スコア
-7
-7
-7
-7
GKY
D
G F V D
5 -7 +7
-7 +4 = 2
GK Y D
GF V D
-7 -7 -1 +0
-7 -7 = -29
GKY D
-7 -7 -5 -7
-7 -7 -7 = -47
G
F V D
動的計画法による大域アライメント(2)
F(0,0)
G
F(1,0)
=-d
K
F(2,0)
=-2d
=0
G
F(0,1)
=-d
F(i-1, j-1)
F(i, j-1)
DP (動的計画法)による
最長経路(スコア)の計算
F (0, j )   jd , F (i,0)  id
 F (i  1, j  1)  s ( xi , y j )

F (i, j )  max
F (i  1, j )  d

F (i, j  1)  d

s(K,F)
F
-d
F(0,2)
=-2d
F(i-1, j)
-d
F(i, j)
⇒ O(mn)時間
行列からの経路の復元は、
F(m,n)からmaxで=となっている
F(i,j)を逆にたどることに行う
(トレースバック)
ローカルアライメント(1)
(Smith-Watermanアルゴリズム)



配列の一部のみ共通部分があることが多い
⇒共通部分のみのアラインメント
配列検索において広く利用されている
例えば、HEAWGEH と GAWED の場合、
AWGE
A W -E
というアライメントを計算
ローカルアライメント(2)
動的計画法
の式
0
 F (i  1, j  1)  s ( x y )

i, j
F (i, j )  max
 F (i  1, j )  d
 F (i, j  1)  d
LAカーネル



SWアルゴリズムをカーネルとして利用したい
⇒ MAX 操作のためカーネルとならない
一方、ペアHMMはカーネルとなることが既知
本研究



SWアルゴリズムを模倣するペアHMMを構成
SWアルゴリズム: 最適なパスのみ
LAカーネル:


全てのローカルアライメントの(重みつき)和
両者ともに時間計算量はO(mn)だが、LAカーネルの
方が数十倍、遅い
LAカーネルの定義(1)

文字(残基)ペアのスコア: Kaβ (x,y)
if | x | 1 or | y | 1
0
K a ( x, y)  
 exp(s( x, y)) otherwise


ギャップのスコア: Kgβ (x,y)
K a ( x, y )  exp ( g (| x |)  g (| y |) )


ただし、 g (0)  0, g (n)  d  e(n  1)
LAカーネルの定義(2)

カーネルの畳み込み(convolution)
K1  K2 ( x, y) 


x x  x, y y  y
ギャップなしブロックが n 個ある場合のスコア




 n 1
K ( n ) ( x, y )  K 0  K a  K g

K1 ( x , y ) K2 ( x , y )
 K a  K 0
ただし、 K 0  1
LAカーネル


K LA
( x, y)   K (i ) ( x, y)
i 0
F V- - E K L GAV- - T
F L L DDRL - - VL L T
Kaβ Kgβ Kaβ Kgβ Kaβ Kgβ Kaβ
LAカーネルとSWスコアの関係



π:(ローカル)アラ
イメント
S(x,y,π): アライメン
トπのスコア
Π:可能なアライメ
ントの集合
定理
SW ( x, y)  max S ( x, y,  )
  ( x , y )
 1 ln max exp( S ( x, y,  ))
  ( x, y )


K LA
( x, y) 

exp( S ( x, y,  ))


( x , y )
lim  ln(K LA ( x, y))  SW ( x, y)
1
 
LAカーネルとSWスコア


SWスコア: 1個の最適なアライメントのみを考慮
LAカーネル: すべての可能なアライメントを考慮
配列 x
HAWGEG
配列 y
AGEHV
SWスコア
AWGE
A - GE
LAカーネル
π1
AWGE
A - GE
S(x,y,π)=24
π2
AWGE
AG - E
S(x,y,π)=19
π3 HAWGE
S(x,y,π)=14
π4 HAWGE - G
S(x,y,π)=5
A -G -E
A -G EHV -
SVMによるスーパーファミリー予測 (1)


各スーパーファミリーごとにSVMを作成
最も高いスコアを出したSVMに対応するファミリーを出力
Super Family.1
SVM.1
Super Family.2
SVM.2
Super Family.3
SVM.3
タンパク質配列
madqlteeqiaefkeafslfdkdgdgtittkel
gtvmrslgqnpteaelqdminevdadg
ngtidfpefltmmark
SVMによるスーパーファミリー予測 (2)

黄色の○


テストデータ
赤い×


スーパーファミリーに所
属するタンパク酸配列
データ(正例)
スーパーファミリーに所
属しないタンパク質配列
データ(負例)
緑の□

どのスーパーファミリーに
所属するかを予測したい
タンパク質配列データ(テ
ストデータ)
margin
カーネル計算の並列化

LAカーネルの計算

K LA
( x, y)
1回あたりO(mn)時間(|x|=m, |y|=n)
 配列データの個数を N とし、配列の平均長を n
⇒ 全部で O(N2n2) 時間 ⇒ 並列化が必要

x1 x2 x3 x 4
x1
x2
x3
x4
並列化
CPU1
CPU2
CPU3
並列計算機の利用

LAカーネルの計算

1回あたりO(mn)時間だが大量の計算が必要



並列計算機



SGI ORIGIN 3800 (R14000(500MHz) × 256CPU)
PCクラスタ HPC (2.8GHz Xeom × 8CPU)
並列化



学習データ中のすべての配列ペアに対して計算
1CPUだと数十日を要する
LSF (Load Sharing Facility) と script の組み合わせ
単純なデータ分割(分割されたデータごとに別CPUで計算)
半日程度でカーネル計算が終了

並列化手法は単純だが、非常に有効
提案手法の評価法
ROCによる性能評価
カーブが
上にある
ほど良い
性能
mRFPによる性能評価
カーブが
上にある
ほど良い
性能
結論

タンパク質ホモロジー検出のための新たなカーネル



Smith-WatermanアルゴリズムとペアHMMの組み合わせ
ベンチマークテストにおいては最高クラスの性能
単純な並列化が非常に有効
課題

タンパク配列の個数(学習データ数)が少ないスー
パーファミリーの予測
グラフカーネルによる化合物の
性質予測
グラフ・カーネル

グラフ

情報科学において幅広く利用されているデータ表現法



頂点と辺で構造を表す(点と線で構造を表す)
V: 頂点の集合
E: 辺の集合
バイオインフォマティクスにおいても幅広い利用


G(V,E)
化学構造、遺伝子ネットワーク、代謝ネットワーク
グラフカーネル

二つのグラフ G1(V1,E1) 、G2(V2,E2) 間の類似性の指標
G(V,E)
Marginalized カーネル


Tsudaらが2002年に提案
定義
K ( x, y)    p(h | x) p(h'| y)K '(( x, h), ( y, h'))
h h'


h,h’: 隠れ変数群、K’:カーネル
配列解析やRNA二次構造解析に応用
Marginalized グラフ・カーネル(1)

Kashimaらが2003年に提案
K (G1,G2 )    p(h) p(h')K '(l (h), l (h'))
hV1* h'V2*




h: グラフ G1 におけるパス
h’: グラフ G2 におけるパス
l(h): パス h のラベル(原子名)の列
K’(x,y): ラベル列間のカーネル関数
(例: K’(x,y)=1 if x=y, otherwise 0
)
Marginalized グラフ・カーネル(2)
G1
u3
u1
H
u2
C
O
u4
Cl
G2
H
v1 v2
H
h  (u1 , u2 , u3 ), h'  (v1 , v2 , v5 ) 
l (h)  (H, C, O), l (h' )  (H, C, O)
K ' (l (h), l (h' ))  1
v4
C
v3
v5
O
v6
H
H
h' '  (u1 , u 2 , u3 , u 2 , u 4 ) 
l ( h' ' )  ( H, C, O, C, Cl)
h' ' '  (v1 , v2 , v5 , v2 , v1 ) 
l ( h' ' ' )  ( H, C, O, C, H )
K ' (l ( h' ' ), l ( h' ' ' ))  0
Marginalized グラフ・カーネル(3)
x
H
O
C
φ(x)
Cl
H
C
N
( 0.03 0.03 0.0
H
C
H
O
C
H
H
C
H
0.02
0.0
0.01
0.002
)
Marginalized グラフ・カーネル(4)
p(v1 , v2 , v3 ) 
0.25 0.9  0.3  0.1
p(v2 , v4 , v2 , v3 ) 
0.25 0.3  0.9  0.3  0.1
p0 (vi )  0.25 pq (vi )  0.1
pa (v2 | v1 )  1.0
(1  p0 ) pa (v2 | v1 )  0.9
pa (v3 | v2 )  1 / 3
(1  p0 ) pa (v3 | v2 )  0.3
G1
START
0.25
0.25
0.25
0.25
O
0.3
v1
H
0.3
v2
0.9
C
0.3
0.9
0.9
0.1
v3
0.1
v4
Cl
0.1
END
0.1
Marginalized グラフ・カーネル(5)

グラフカーネルの定義
K (G1,G2 )    p(h) p(h')K '(l (h), l (h'))
hV1* h'V2*



h: グラフ G1 におけるパス
h’: グラフ G2 におけるパス
すべてのパスを考慮するので、そのまま実行した
のでは指数時間かかるが、工夫すると逆行列の
計算に帰着できる
(|V1||V2|×|V1||V2|サイズの逆行列の計算)
Marginalized グラフ・カーネル(6)
p s ( v )  p0 ( v ) p q ( v )
pt (u | v) 
1  pq (v )
pq ( v )
pa (u | v) pq (u )
n
p(v1 ,, vn )  ps (v1 ) pt (vi | vi 1 )
i 2
Marginalized グラフ・カーネル(7)
 s (u1 , v1 )  ps(1) (u1 )  ps( 2) (v1 )
 t ((ui , vi ) | (ui 1 , vi 1 ))  pt(1) (ui | ui 1 )  pt( 2) (vi | vi 1 )
n
 ((u1 , v1 )(u2 , v2 )  (un , vn ))   s (u1 , v1 )  t ((ui , vi ) | (ui 1 , vi 1 ))
i 2
G1 u 1
H
G2
Cl
K
O
v1
u2 G1×G2 ( u1 , v1 ) ( u1 , v2 ) ( u1 , v3 )
v2
v3
H
H,K
H,O
H,H
Cl,K
Cl,O
Cl,H
( u2 , v1 ) ( u2 , v2 ) ( u2 , v3 )
Marginalized グラフ・カーネル(8)
K (G1,G2 )    p(h) p(h')K '(l (h), l (h')) 
hV1* h'V2*
 s  ( s (v))vV (V  V1  V2 )
 t  ( t (u' | u))u 'V ,uV


 (h)
h(V1V2 )*
 (h)   s  ( t ) n 1  1
hV * , |h|  n



 
n 1


K (G1,G2 )  
 (h)    s ( t ) 1



*
n 1 hV , |h| n
 n1
2

1
1

x

x
   1 /(1  x)
  s  (I  t ) 1


Marginalized グラフ・カーネル⇒逆行列の計算
Marginalized グラフカーネルの問題点

パス(の集合)だけを用いて化学構造を表現


反応中心などの情報を十分に取り入れることが困
難?
行列のサイズが大きく(化合物xの原子数×化
合物yの原子数)2)なるため、逆行列の計算に
時間がかかる

すべてのトレーニングデータのペア(化合物のペア)
について、それぞれ、逆行列を計算することが必要
⇒ 構造情報(Morgan Index)との組み合わせ
により行列のサイズを減らす
Morganインデックス

化学構造の一意名を計算機により計算するために
1960年代に考案



CAS(Chemical Abstract Service)で利用
等価な原子に同じ番号(整数値)が与えられるような、
各原子への番号づけを計算
簡単な繰り返し計算による番号づけ

等価で無い原子にも同じ番号がつく可能性(でも、低い)
⇒ Marginalized グラフカーネルにおいて、原子名ととも
に、モーガンインデックスを利用
原子名およびモーガンインデックスの両者が一致するパス
のみを考慮
⇒ 部分構造に関する特徴も、ある程度、取り入れられる

Morganインデックスの計算法


すべての原子に番号1を割り当てる
すべての原子 x について以下を実行

x に結合している原子の番号を総和を、x の番号とする
1
2
1
1
1
1
1
1
1
1
N
O
1
2
2
1
4
3
2
2
3
2
2
1
5
4
3
7
5
5
7
5
6
7
1
1
N
1
3
N
3
O
O
3
O
O
5
O
計算機実験

MUTAG データを利用





標準的ベンチマークテストの一つ
化合物のサルモネラ菌の変異性への影響データ
125個の正例、63個の負例を利用
各例1個のみをテストデータとし、他を学習データと
したテストを繰り返した
ソフトウェア
SVMソフトとして、GIST
(http://microarray.cpmc.columbia.edu/gist)
を利用
 他は C++ で記述

計算機実験の結果: 予測精度
Marginalized
カーネル
+
モーガン法
他手法
計算機実験の結果: 計算時間
結論

モーガンインデックスの利用により以下を達成


Marginalizedカーネルと、同様の精度
数十倍以上、高速
今後の課題



他のインデックス手法の利用、開発
他手法との比較
大規模な計算機実験(現在、実行中)

並列処理が必要
参考文献

SVMおよびカーネル一般


バイオインフォマティクスにおけるカーネル


Kernel Methods in Computational Biology, MIT Press, 2004.
Marginalized Kernel + Morgan Index


N. Cristianini & J. Shawe-Taylor: An Introduction to Support Vector
Machines and Other Kernel-based Learning Methods, Cambridge Univ.
Press, 2000.
P. Mahe, N. Ueda, T. Akutsu, J-L. Perret, J-P. Vert: Extensions of
marginalized graph kernels, Proc. 21st Int. Conf. Machine Learning,
552-559, 2004.
LAカーネル

H. Saigo, J-P Vert, N. Ueda, T. Akutsu: Protein homology detection
using string alignment kernels, Bioinformatics, 20:1682-1689, 2004.