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

配列および化合物データ解析
のためのカーネル法
阿久津達也
京都大学 化学研究所
バイオインフォマティクスセンター
サポートベクターマシン




カーネル法の一つ、ニューラルネットワークと類似
1990年代に、Cortes と Vapnik が発明
トレーニングデータとして与えられた正例と負例から、
それらを分離する超平面を計算
機械学習、統計学、人工知能、パターン認識、バイオ
インフォマティクスなど様々な分野に応用






配列分類
タンパク質フォールド予測、二次構造構造
遺伝子発現データ解析
タンパク質相互作用予測
化合物の性質推定
c.f. Kernel Methods in Computational Biology, MIT Press,
2004
サポートベクターマシン


正例と負例を与
えて、それらを最
適（マージンを最
大）に分離する超
平面を学習
カーネルを適切に
定義することによ
り超平面以外で
の分離が可能
テス トデー タ
m argin
SVMによるテストデータの分類


学習データより超平面
を学習(SVM)
テストデータは、対応す
る点の超平面に対する
位置（上下）で判定

テス トデー タ
テストデータとサポート
ベクター間のカーネル関
数値の重み付き和でテ
ストデータを類別
m argin
f (x )

i: x i   

i
K ( x ,x 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 )dx dz  0

離散値の場合 ( x1,x2,…,xn が入力データ)
K  ( K ( x i , x j ))
n
i , j 1
が半正定値行列
カーネルの作り方
データから特徴ベクトル(feature vector)を
作るのが一般的、かつ、
多くの場合に実用的
 特徴ベクトル： 実数値の列
 例えば、各化合物 x に対し、


Φ(x) = (分子量, 容積, 表面積, logP,…)
とすれば、化合物 x,y に対するカーネルは
Φ(x) と Φ(ｙ) の単なる内積
アライメントカーネルによる構造予測
1.
2.
3.
4.
5.
6.
SCOPとスーパーファミリー予測
既存カーネル
配列解析手法（アライメント、HMM）
新カーネル
計算機実験結果
結論と課題
タンパク質立体構造予測



アミノ酸配列から、
タンパク質の立体
構造（３次元構造）
をコンピュータによ
り推定
実験よりは、精度が
悪い
だいたいの形がわ
かれば良いのであ
れば、4～５割の予
測率
タ ン パ ク質 配 列 （ア ミ ノ酸 配 列 ）
T
C
A
V
F
G
L
G
G
V
R
L
S
V
D
コンピュー タ
タンパ ク質
立体構造
フォールド予測 (Fold Recognition)
タ ン パ ク質 配 列 （ア ミ ノ酸 配 列 ）


精密な３次元構造
ではなく、だいたい
の形（fold)を予測
立体構造は数千
種類程度の形に分
類される、との予
測(Chotia, 1992)
に基づく
T
C
A
V
F
G
L
G
G
V
R
L
1 000 個 の テ ン プ レー ト 構 造
S
V
D
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)
S p ectrum カ ー ネ ル
A A
ACCT A
( 0
φ (x)
A C
A G
CC
CG
CT
T A
1
0
1
0
1
1
)
配列アライメント




バイオインフォマティクス
の最重要技術の一つ
２個もしくは３個以上の配
列の類似性判定に利用
文字間の最適な対応関
係を求める（最適化問題）
配列長を同じにするよう
に、ギャップ記号（挿入、
欠失に対応）を挿入
A
L
G
F
G
S
L
Y
A
L
G
G
V
S
V
G
A
L
G
F
G
S
L
A
L
G
S
V
G
V
G
Y
G
G
スコア行列（置換行列）

残基間（アミノ酸文字間）の類似性を表す行列

PAM250, BLOSUM45 など
A
A
R
N
D
C
Q
E
G
H
I
L
K
M
F
P
S
T
W
Y
V
R
N
D
C
Q
E
5 -2 -1 -2 - 1 -1 -1
-2
7 - 1 -2 -4
G
P
S
T
Y
V
0 -2 - 1 -2 -1 -1 -3 - 1
H
I
L
K
M
F
1
0 -3 -2
W
0
3
1
0 -3
0 -4 -3
3 - 2 -3 -3 -1 -1 - 3 -1
-1 -1
7
2 -2
0
0
1 -3 -4
0 - 2 -4 -2
-2 -2
2
8 -4
0
2 - 1 -1 -4 -4 -1 - 4 -5 -1
0
1
0 - 4 -2 -3
0 -1 - 5 -3 -4
-1 -4 - 2 -4 13 -3 -3 - 3 -3 -2 -2 -3 - 2 -2 -4 -1 -1 - 5 -3 -1
-1
1
0
0 -3
7
2 -2
1 -3 -2
2
0 -4 -1
B LO S U M 5 0 ス コ ア 行 列
（置 換 行 列 ）の 一 部 分
0 -1 - 1 -1 -3
ペアワイズ・アライメント


配列が２個の場合でも可能なアライメントの個数は
指数オーダー
しかし、スコア最大となるアライメント（最適アライメ
ント）は動的計画法により、O(mn)時間で計算可能
（m,n:入力配列の長さ）
入力配列
A G C T, AC G C T
イメント
ア ラ イ メント
AGCT AC G C T
スコア
-3
最 適 アラ
AG - CT
AC G C T
1
A - GCT
AC G C T
- AGC - - T
AC - - G C T
3
-5
（ 同 じ 文 字 の 時 : 1、 違 う文 字 の 時 : -1 、ギ ャ ッ プ 1 文 字 : - 1）
動的計画法による大域アライメント(1)
(Needleman-Wunschアルゴリズム)



入力文字列から格子状グラフを構成
アライメントと左上から右下へのパスが一対一対応
最長経路＝最適アライメント
G
G
F
5
K
-2
Y
D
-5
1
ア ラ イメント
-7
G K Y
-5
-5
7
-6
-7
G
D
F V D
G K Y D
V
D
-1
-2
-3
-2
1
0
-4
4
-7
-7
-7
-7
スコア
-7
G
-7 +4 = 2
-7 -7 -1 +0
-7 -7 = -29
G F V D
-7 -7 -5 -7
G K Y D
-7
5 -7 +7
F V D
-7 -7 -7 = -47
動的計画法による大域アライメント(2)
F (0 ,0 )
G
F (1 ,0 )
=-d
K
F (2 ,0 )
= -2 d
DP (動的計画法)による
最長経路(スコア)の計算
F ( 0 , j )   jd ,
=0
G
F ( i -1 , j -1 )
F (0 ,1 )
=-d
F (i, j- 1 )
F ( i , 0 )   id
 F ( i  1, j  1)  s ( x i , y j )

F ( i , j )  max 
F ( i  1, j )  d

F ( i , j  1)  d

s (K ,F )
F
-d
F (0 ,2 )
= -2 d
-d
F ( i -1 , j )
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 i , y 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)

K 1  K 2 ( x, y ) 
K1(x , y )K 2 (x , y )
x x  x, y y  y

ギャップなしブロックが n 個ある場合のスコア
K ( n ) ( x , y )  K 0  K a  K g





n 1

 Ka  K0
ただし、
K
0
 1
LAカーネル

K LA ( x , y ) 


i0

K (i) ( x, y )
F V - - E K L G A V - - T
F L L D D R L - - V L L T
β
Ka
β
Kｇ
β
Ka
β
Kｇ
β
β
β
Ka Kｇ Ka
LAカーネルとSWスコアの関係



π：(ローカル)アラ
イメント
S(x,y,π): アライメン
トπのスコア
Π：可能なアライメ
ントの集合
SW ( x , y )  max
  ( x , y )


1

S ( x, y, )
ln  max exp(  S ( x , y ,  )) 
   ( x , y )

K LA ( x , y ) 
 exp(
 S ( x , y ,  ))
  ( x , y )
定理
lim
 
1


ln( K LA ( x , y ))  SW ( x , y )
LAカーネルとSWスコア


SWスコア： 1個の最適なアライメントのみを考慮
LAカーネル： すべての可能なアライメントを考慮
配列 x
H A W G EG
配列 y
AGEH V
LAカ ー ネ ル
π1
AWGE
A - GE
S (x,y, π )= 2 4
π2
AWGE
AG - E
S (x,y, π )= 1 9
π3
HAWGE
A -G -E
S (x,y, π )= 1 4
π4
HAWGE - G
A - G EH V -
S (x,y, π )= 5
SWス コア
AWGE
A - GE
SVMによるスーパーファミリー予測 (1)


各スーパーファミリーごとにSVMを作成
最も高いスコアを出したSVMに対応するファミリーを出力
Super Family.1
SVM.1
Super Family.2
SVM.2
Super Family.3
SVM.3
タンパク質配列
madqlteeqiaefkeafslfdkdgdgtittkel
gtvmrslgqnpteaelqdminevdadg
ngtidfpefltmmark
SVMによるスーパーファミリー予測 (2)

黄色の○


テス トデー タ
赤い×


スーパーファミリーに所
属するタンパク酸配列
データ（正例）
スーパーファミリーに所
属しないタンパク質配列
データ（負例）
緑の□

どのスーパーファミリーに
所属するかを予測したい
タンパク質配列データ（テ
ストデータ）
m argin
カーネル計算の並列化

LAカーネルの計算
x1 x2 x3 x 4
x1
x2
x3
x4

K LA ( x , y )
 １回あたりO(mn)時間（|x|=m, |y|=n）
 配列データの個数を N とし、配列の平均長を n
⇒ 全部で O(N2n2) 時間 ⇒ 並列化が必要
並列化
CPU1
CPU2
CPU3
並列計算機の利用

LAカーネルの計算

１回あたり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 , G 2 )  
 p ( h ) p ( h ' ) K ' (l ( h ), l ( h ' ))
*
*
hV1 h 'V 2




h: グラフ G1 におけるパス
h’: グラフ G2 におけるパス
l(h): パス h のラベル（原子名）の列
K’(x,y): ラベル列間のカーネル関数
（例： K’(x,y)=1 if x=y, otherwise 0
）
Marginalized グラフ・カーネル(2)
u3
G1
u1
u2
H
C
O
G2
u4
Cl
v1
H
h  ( u 1 , u 2 , u 3 ), h '  ( v1 , v 2 , v 5 ) 
l ( h )  ( H, C, O ), l ( h ' )  ( H, C, O )
K ' ( l ( h ), l ( h ' ))  1
H
v2
v4
v3
v5
C
O
v6
H
H
h ' '  (u1 , u 2 , u 3 , u 2 , u 4 ) 
l ( h ' ' )  ( H, C, O, C, Cl )
h ' ' '  ( v1 , v 2 , v 5 , v 2 , v1 ) 
l ( h ' ' ' )  ( H, C, O, C, H )
K ' ( l ( h ' ' ), l ( h ' ' ' ))  0
Marginalized グラフ・カーネル(3)
x
H
O
C
φ (x)
Cl
H
( 0.03
H
H
H
C
C
C
N
C
O
H
H
0.03
0.0
0.02
0.0
0.01
0.002
)
Marginalized グラフ・カーネル(4)
p ( v1 , v 2 , v 3 ) 
0 . 25  0 . 9  0 . 3  0 . 1
G1
p (v2 , v4 , v2 , v3 ) 
0 . 25  0 . 3  0 . 9  0 . 3  0 . 1
p 0 ( v i )  0 . 25
p q ( v i )  0 .1
START
0 .2 5
p a (v3 | v 2 )  1 / 3
(1  p 0 ) p a ( v 3 | v 2 )  0 . 3
0 .2 5
O
0 .3
v1
H
0 .3
v2
C
0 .3
0 .9
0 .1
v3
0 .9
0 .9
p a ( v 2 | v1 )  1 . 0
(1  p 0 ) p a ( v 2 | v1 )  0 . 9
0 .2 5
0 .2 5
0 .1
v4
Cl
0 .1
EN D
0 .1
Marginalized グラフ・カーネル(5)

グラフカーネルの定義
K (G1 , G 2 )  
 p ( h ) p ( h ' ) K ' (l ( h ), l ( h ' ))
*
*
hV1 h 'V 2



h: グラフ G1 におけるパス
h’: グラフ G2 におけるパス
すべてのパスを考慮するので、そのまま実行した
のでは指数時間かかるが、工夫すると逆行列の
計算に帰着できる
(|V1||V2|×|V1||V2|サイズの逆行列の計算)
Marginalized グラフ・カーネル(6)
p s (v )  p0 (v ) p q (v )
p t (u | v ) 
1  p q (v )
p q (v)
p a (u | v ) p q (u )
n
p ( v1 ,  , v n )  p s ( v1 )  p t ( v i | v i  1 )
i2
Marginalized グラフ・カーネル(7)
 s ( u 1 , v1 )  p s ( u 1 )  p s ( v1 )
(1 )
(2)
 t (( u i , v i ) | ( u i 1 , v i 1 ))  p t ( u i | u i 1 )  p t ( v i | v i 1 )
(1 )
(2)
n
 (( u 1 , v1 )( u 2 , v 2 )  ( u n , v n ))   s ( u 1 , v1 )   t (( u i , v i ) | ( u i 1 , v i 1 ) )
i2
G 1 u1
H
G 1× G 2 ( u 1 , v1 ) ( u 1 , v2 ) ( u 1 , v3 )
H ,K
H ,O
H ,H
Cl
u2
G2
v1
v2
v3
K
O
H
C l,K
C l,O
C l,H
( u 2 , v1 ) ( u 2 , v2 ) ( u 2 , v3 )
Marginalized グラフ・カーネル(8)
K (G1 , G 2 )  
 p ( h ) p ( h ' ) K ' (l ( h ), l ( h ' )) 
*
*
hV1 h 'V 2
 s  ( s ( v )) vV
(V  V  V )
1
2
 t  ( t ( u ' | u )) u 'V , u V
h (V1 V 2 )

 (h)

*
 (h )   s  ( t )
n 1
*
hV , | h |  n

 
n 1


K ( G1 , G 2 )  
  ( h )    s (  t ) 1


*
n 1 hV , | h | n

 n 1
2
1
1

x

x
   1 /( 1  x )
  s  ( I   t ) 1


1


Marginalized グラフ・カーネル⇒逆行列の計算
Marginalized グラフカーネルの問題点

パス（の集合）だけを用いて化学構造を表現


反応中心などの情報を十分に取り入れることが困
難？
行列のサイズが大きく（化合物xの原子数×化
合物yの原子数）2）なるため、逆行列の計算に
時間がかかる

すべてのトレーニングデータのペア（化合物のペア）
について、それぞれ、逆行列を計算することが必要
⇒ 構造情報(Morgan Index)との組み合わせ
により行列のサイズを減らす
Morganインデックス

化学構造の一意名を計算機により計算するために
1960年代に考案



CAS(Chemical Abstract Service)で利用
等価な原子に同じ番号（整数値）が与えられるような、
各原子への番号づけを計算
簡単な繰り返し計算による番号づけ

等価で無い原子にも同じ番号がつく可能性（でも、低い）
⇒ Marginalized グラフカーネルにおいて、原子名ととも
に、モーガンインデックスを利用
原子名およびモーガンインデックスの両者が一致するパス
のみを考慮
⇒ 部分構造に関する特徴も、ある程度、取り入れられる

Morganインデックスの計算法


すべての原子に番号１を割り当てる
すべての原子 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個の正例、６３個の負例を利用
各例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.