PowerPoint プレゼンテーション

生命情報解析 第3回
シグナル配列の統計解析(2)
慶應義塾大学先端生命科学研究所
塩基の偏りの程度を
エントロピーで測る
エントロピー H  
 P log
i a,c,g, t
i
2
Pi
但しPiは塩基iの頻度(0≦ Pi≦1), 0log0 = 0と定義
•0≦H≦2
• 塩基の偏りが強い → Hが0に近づく
• 塩基の偏りが弱い → Hが2に近づく
A,C,G,Tが1/4ずつなら
A,Cが1/2, C,Gが0なら
全てAなら
H=2
H=1
H=0
大腸菌開始コドン周辺の
塩基のエントロピー
2.5
エントロピー
2
1.5
1
0.5
0
-100
-50
0
開始コドンからの相対位置
50
100
シグナル配列が見られるところで、値が
下がるというのは感覚に合わない?
• Schneider T et al. 1986はRseqを導入した
Rseq = エントロピーの最大値 – 対象位置のエントロピー
• 塩基の偏りが強い → Rseqが2に近づく
• 塩基の偏りが弱い → Rseqが0に近づく
• 不確定性の減少の度合い
最大エントロピーから
対象位置のエントロピーを引く
と…
2.5
2
R seq
1.5
1
0.5
0
-100
0
開始コドンからの相対位置
100
配列ロゴ
http://www.lecb.ncifcrf.gov/~toms/gallery/ribo.logo.gif
タンパク質のリン酸化部位の解
析
http://en.wikipedia.org/wiki/Protein_kinaseより抜粋
•
•
•
•
細胞は、その機能を維持するため、細胞内のタンパク質をリン酸化、
脱リン酸化する反応を繰り返している。
このリン酸化によってタンパク質は酵素活性、細胞内での局在や他の
タンパク質との会合状態を変化させる。
細胞内の30%ものタンパク質がキナーゼによる変化を受け、細胞内に
おける様々なシグナル伝達や代謝の調節因子として機能している。
キナーゼ遺伝子はヒトゲノム中に約500種類があり、また真核生物の全
遺伝子の約2%を占める。
P
P
P
ADSLQMWSA
P
MAALLSL
ADSLQMWSWLLW
タンパク質リン酸化部位の
シグナル配列
Schwartz D and Gygi SP(2005)
マイコプラズマ菌の
開始コドン周辺のエントロピー
2.5
エントロピー
2
1.5
1
0.5
0
-100
0
開始コドンからの相対位置
100
ゲノムの塩基組成を反映させるに
は?
• マイコプラズマ菌などはゲノム全体のAT含量
が高い
• シグナルがないところでもエントロピーが低
くなってしまう
• ゲノムの塩基組成を考慮したい
増加情報量(1)
• ゲノム全体の塩基iの組成をBiとする。
• シグナル配列の塩基iの組成をPiとする。
• 与えられた位置がシグナル配列であることが判明、
その情報量I(Pi//Bi) は?
増加情報量(2)
• 情報量の加法性より、
塩基iが判明したときに得られる情報量
= シグナル配列であることが判明したときに得られる情
報量
+その上でさらに塩基iが判明したときの情報量
• 式で表すと、-log Bi = I(Pi//Bi) + -log
Pi
従って、 I(Pi//Bi) =-log Bi --log Pi = log Pi / Bi
増加情報量(3)
• シグナル配列上で期待値をとって、
増加情報量I ( P // B ) 
 P  I ( P // B )
i  a,c,g, t
i
i
i
Pi
  Pi log2
Bi
i  a,c,g, t
但しBiは塩基iの対照となる頻度(0≦ Bi≦1)と考える
大腸菌開始コドン周辺の
塩基の増加情報量
2.5
増加情報量
2
1.5
1
0.5
0
-100
-50
0
開始コドンからの相対位置
50
100
3
マイコプラズマ菌開始コドン周
辺の
塩基の増加情報量
2.5
増加情報量
2
1.5
1
0.5
0
-100
0
開始コドンからの相対位置
100
演習問題
Ba=0.3,Bc=0.3,Bg=0.3,Bt=0.1
として、与えられた位置における塩基iの個
数Niが
Na = 40, Nc =40, Ng = 10, Nt = 10
のときの増加情報量を求めよ。
log23 ≒1.585
統計的有意性の検定
• これまでの議論はシグナル配列中の塩
基がどれくらい偏っているかだった
• しかし、このままではその偏りが偶然
に生まれている可能性は否定できない
• 偶然に起きる確率の計算 → 統計的有意性
偏りのないサイコロ?
(歪みのないサイコロ)
• サイコロを10回振っ
たとき、6の目が9回
出た
• 6の目が10回中、9回
以上も出る確率は
(1/6)9(5/6)1×10 + (1/6)10
≒8.4344675608 × 10-7
-7
P
<
8.44
×
10
確率計算によって有意性を調べる
確率分布と有意性(1)
サイコロを10回ふったときの6の目が出る回数とその確率との関係
0.35
0.30
P = 0.00243815649926
0.25
P < 0.003
確率
0.20
0.15
0.10
0.05
0.00
0
1
2
3
4
5
6
7
8
9
10
6の目が出る回数
棒グラフの右側部分の面積の合計が確率、すなわち有意性を表す
確率分布と有意性(2)
サイコロを100回ふったときの6の目が出る回数とその確率との関係
1.20E-01
1.00E-01
確率
8.00E-02
試行回数が増えると、正規分布と
呼ばれる左右対称の分布に近づく
6.00E-02
4.00E-02
2.00E-02
0.00E+00
0
10
20
30
40
50
60
6の目が出る回数
70
80
90
100
正規分布
Gaussian Distribution
p
μ
1
 2 ( x )2
1
p( x |  ) 
e 2
2 
σ
μ:平均、σ:標準偏差
x
• 多くの自然現象、社会現象は正規分布になる?
–
–
–
–
サイコロを1000回振ったときに、6の目が出る回数
工業製品の品質(重さ,長さなど)、多くの要因が重なる誤差
大学1年生男子の身長
一日の降水量?(対数をとる)
• サンプル数が多ければ、その平均は正規分布に従う(中心極限定
理)
Z Scoreによる分布の変換
Z Score =
N
Nobs
p
N obs  Np
Np(1  p)
:試行回数
:対象となる現象の観測数
:対象となる現象が起こる確率
Z Scoreの特徴
0
Z Scoreの分布
1.96
1.20E-01
• 平均が0、分散(データの散らばり)が1
になる
• 元の分布が正規分布なら、そのZ
Scoreは標準正規分布となる
• どんな正規分布でも、Z Scoreに直せ
ば同じ土俵(?)で確率計算ができる
• Z Scoreが1.96を超える確率は0.025
1.00E-01
確率
8.00E-02
6.00E-02
4.00E-02
2.00E-02
0.00E+00
-4.47 -3.13 -1.79 -0.45 0.89
2.24
3.58
4.92
6.26
7.60
8.94 10.29 11.63 12.97 14.31 15.65 16.99 18.34 19.68 21.02 22.36
Z S core
Z Scoreの計算
サイコロを100回振って、90回”6”の目が出るときのZ Scoreは
N obs  Np
90  1001 / 6
Z Score 

 21.17
Np(1  p)
901 / 6  5 / 6
塩基の方も…
• ここでは簡単のため、1塩基の偏りだ
けを考える
• ゲノム全体の塩基組成を考えて、塩基i
が対象となる場所において観測される
確率はpiとする
• 今、N本の配列のうち、 Ni個について、
対象となる位置に塩基iが観測された
• この条件では通常、Niは正規分布に従
う
頻出塩基の統計的有意性
Z Score =
N
Ni
Bi
N i  NBi
NBi (1  Bi )
:解析する配列数
:観測された塩基iの数
:ゲノム中における塩基iの割合
Z Scoreは標準正規分布に従う
Z Score > 1.96なら、P < 0.05
演習問題
ゲノム全体でGの組成は0.2であるとする。
今、200本の配列の特定の位置において、
Gが20本の配列で観測された。このとき
のGの個数のZ Scoreを求めよ。
複数の塩基の有意性を
同時に検定するには?
• ゲノム全体の塩基組成をBa=0.3, Bc=0.3,
Bg=0.3, Bt=0.1として、与えられた位置
における塩基iの個数Niが Na = 40, Nc =40,
Ng = 10, Nt = 10のとき、偏りは有意か?
• Z Scoreを4つも計算すると…
– 4つも値が出て、取り扱いが煩雑になる
– 偶然に高い値を示すものが出やすくなる
• Χ2値を使う
Χ2値
• n個の互いに独立なZ Score:Z1, Z2, Z3, …,
Znがあるとき、
Χn2値 = Z12+Z22+Z32+…+Zn2
• Χn2値は自由度nのΧn2分布に従う
Χ2分布
1.4
1.2
1
自由度
1
2
3
4
5
確率
0.8
0.6
0.4
0.2
0
0
0.5
1
1.5
2
2.5
3
3.5
4
4.5
5
2
Χ 値
5.5
6
6.5
7
7.5
8
8.5
9
9.5
10
Χ2分布に従う値を求める
• しかし、Na, Nc, Ng, Ntは互いに独立ではないた
め、これらのZ Scoreを足しても自由度4のΧ2
分布には従わない。そこで…
2
(
P

B
)
2
i
3 値  N  i
Bi
i a,c,g, t
但しNは解析対象の配列数
• 塩基がそれぞれBiの頻度で出現するとき、上記
Χ2値は自由度3のΧ2分布に従う
• Χ2値>12.84ならP < 0.005
Χ2計算の例
• ゲノム全体の塩基組成をBa=0.3, Bc=0.3, Bg=0.3, Bt=0.1
として、与えられた位置における塩基iの個数Niが Na
= 40, Nc =40, Ng = 10, Nt = 10のとき、偏りは有意か?
(0.4  0.3) 2 (0.4  0.3) 2 (0.1  0.3) 2 (0.1  0.1) 2



0.3
0.3
0.3
0.1
 0.033 0.033 0.133 0  0.2
• この例では偏りが有意とは言えない。
14000
大腸菌開始コドン周辺の
2
塩基のχ 乗値
12000
Χ 自乗値
10000
8000
6000
4000
2000
0
-100
-50
0
開始コドンからの相対位置
50
100
2つの数式の関係(1)
2
(
P

B
)
2
i
 n1 値  N  i
Bi
i 1
n
n 値  Z12  Z22  Z32   Zn2
2
• 4塩基を2種類に分類して考える
– プリン(A,G)、ピリミジン(C,T)
• ゲノム中のプリン、ピリミジンの頻度をそれぞれBpur,
Bpyrとする
• 対象となる位置で観測されたプリン、ピリミジンの
頻度をそれぞれPpur, Ppyrとする。但し、Ppur + Ppyr = 1
2つの数式の関係(2)
 ( Ppur  Bpur ) 2 ( Ppyr  Bpyr ) 2  ( NPpur  NBpur ) 2
2

N

 Z pur

 NB (1  B )
B
B
pur
pyr
pur
pur


• 2種類の塩基の数をもとに計算したΧ2値は自
由度1のΧ2分布に従う
• 4種類の塩基の数をもとに計算したΧ2値は自
由度3のΧ2分布に従う
• 自由度は自由に動ける変数の数を意味する
演習問題
Ba=0.3,Bc=0.2,Bg=0.2,Bt=0.3として、
与えられた位置における塩基iの個
数Niが
(1) Na = 50, Nc = 30, Ng = 10, Nt =
10
(2) Na = 500, Nc = 300, Ng = 100, Nt
= 100
のときの増加情報量、χ2値を求めよ。
log23 ≒1.585、log25 ≒2.322
演習問題 解答
• (1), (2)ともに
1
5 3
3 1
1 1
1
増加情報量 log2  log2  log2  log2
2
3 10
2 10
2 10
3
 0.368 0.175 0.1  0.159  0.285
• (1)
 ( 12  103 ) 2 ( 103  15 ) 2 ( 101  15 ) 2 ( 101  103 ) 2 

  (50  30  10  10)




3
3
1
1
5
5
10
 10

 100 (0.13  0.05  0.05  0.133)  36.7
• (2)
3 2
3
3 2 
1
1 2
1
1 2
1

(

)
(

)
(

)
(

)
 2  (500 300 100 100) 2 310  10 1 5  10 1 5  10 3 10 
5
5
10
 10

 1000 (0.13  0.05  0.05  0.133)  367
2
特定のシグナル配列の存在頻度
• 様々な塩基配列の偏りを調べるのではな
く、特定のシグナル配列の存在頻度を調
べたい (ex. SD配列 “AGG”)
• 最も単純なのは、頻度=あるシグナル配
列が観測される配列数÷解析対象の配列
数
0.14
大腸菌開始コドン周辺の”AGG”の頻
度
0.12
"AGG"の頻度
0.1
0.08
0.06
0.04
0.02
0
-100
-50
0
開始コドンからの相対位置
50
100
頻出塩基配列パターンの
統計的有意性
Z Score =
N obs  Np
Np(1  p)
Nobs :パターン観測数 N:解析する配列数
p:ゲノム中におけるパターンの割合
パターンの出現頻度がpのとき、Z Scoreは
標準正規分布に従う
Z Score > 1.96なら、P < 0.05
大腸菌開始コドン周辺の”AGG”のZ-Score
80
70
60
Z Score
50
40
30
20
10
0
-100
-50
0
-10
開始コドンからの相対位置
50
100
翻訳開始シグナル抽出結果
Escherichia coli
16S rRNA 3- terminal: gcggttggatcacctcctta3
Expected SD Sequence: 5taaggaggtgatccaaccgc
Pat.
agga
ggag
aagg
gagg
gaga
Z-Sc.
94.97
82.94
58.15
53.08
42.23
Pos.
-11
-10
-11
-11
-9
シグナル配列出現の評価
• 塩基の偏り
– 偏りの程度 … エントロピー、増加情報量
– 偏りの有意性 … Χ2値
• 配列パターン
– 出現の程度 … 頻度
– 出現の有意性 … Z Score