配列アラインメント Sequence alignment

Transcription

1 Traceback 以下の図に示すように 最大重みパスを D[ を用いて頂点 (m,n) から逆にたどりながら アラインメントを後ろから前へと計算していけばよい トレースバック i-,j-)=9 j-)= d=- トレースバック i-,j-)=8 j-)= d=- Traceback Procedure Traceback ( ) k ; i m; j n; while i or j do if i or D[ D[ j ] d do ' '; s [ ; j j ; else if j or D[ D[ i j, d do s [ i]; ' '; i i ; else do s [ i]; s [ ; i i ; j j ; k k ; t (... k ]) return ( t, t ); t (... k ]) ; ; i-,=8 d=- = d=- i-,= = ここで s は s を反転させた配列を表す 例えば s=gtの場合 s TG とする 9+= -= 8/4/4 8/4/4 Traceback Procedure Traceback ( ) k ; i m; j n; while i or j do if i or D[ D[ d do ' '; s [ ; j j ; else if j or D[ D[ i j, d do s [ i]; ' '; i i ; else do s [ i]; ' s [ '; i i ; j j ; k k ; t (... k ]) t (... k ]) return ( t, t ); 計算時間 : ループを 回実行するごとに i か j の少なくともどちらか一方は だけ減少することがわかる 8/4/4 3 ; ; Traceback Procedure Traceback ( ) k ; i m; j n; while i or j do if i or D[ D[ d do ' '; s [ ; j j ; else if j or D[ D[ i j, d do s [ i]; ' '; i i ; else do s [ i]; ' s [ '; i i ; j j ; k k ; t (... k ]) t (... k ]) return ( t, t ); 8/4/4 4 よって ループをたかだか n+m 回実行することによって i=j= となることがわかる ループ 回あたりの計算時間は明らかに定数時間であり ループ全体でかかる時間は O(m+n) となる ; ;

2 Traceback Procedure Traceback ( ) k ; i m; j n; while i or j do if i or D[ D[ d do ' '; s [ ; j j ; else if j or D[ D[ i j, d do s [ i]; ' '; i i ; else do s [ i]; ' s [ '; i i ; j j ; k k ; t (... k ]) t (... k ]) return ( t, t ); ; ; 動的計画法の計算時間 定理 線形ギャップスコアの下での 個の配列に対する最適大域アラインメントは O(mn) 時間及び O(mn) 領域で計算できる 配列の反転は明らかに線形時間で実行できるので 全体の計算時間は O(m+n) となる このように 動的計画法を用いて計算したデーブルの値から それを逆方向にたどって最適解を復元する手続きは一般に トレースバック (traceback) とよばれる 8/4/4 5 問題を定義する : 二本配列 SEQX and SEQY, * 単純スコアリングシステム (MTH, MISMTH, and gap INDEL 変数を定義する #define SEQX "TTT" #define SEQY "TGTGT" #define MTH #define MISMTH - #define GP - #include <stdlib.h> #include <stdio.h> #include <string.h> int main(void) { char *x,*y; /* 二本配列, x_..m and y_..n */ int M,N; /* 配列の長さ x,y */ int j; /* 座標 in x and y */ int **S; /* 動的計画法の行列 */ char *ax,*ay; /* 結果 : アラインメントした二本配列 */ int K; /* pairwise alignment の長さ */ int k; /* pairwise alignmentの変数,..k- */ int sc; /* tmp best pathを探すの変数 */ char tmp; /* tmp swap charsするための変数 */

3 再帰法 : 初期化 M = strlen(seqx); N = strlen(seqy); x = malloc(sizeof(char) * (M+)); /* +: leading blank, and trailing */ y = malloc(sizeof(char) * (N+)); strcpy(x+, SEQX); /* x_..x_m now defined; x_ undefined */ strcpy(y+, SEQY); /* y_..y_n now defined; y_ undefined */ /*, /* 動的計画法の行列..M rows by..n columns.*/ S = malloc(sizeof(int *) * (M+)); for (i = ; i <= M; i++) S[i] = malloc(sizeof(int) * (N+)); S[][] = ; for (i = ; i <= M; i++) S[i][] = i * GP; for (j = ; j <= N; j++) S[][ = j * GP; 再帰法 : DP lgorithm for (i = ; i <= M; i++) for (j = ; j <= N; j++) { /* case #: j are aligned */ if (x[i] == y[) S[i][ = S[i-][j-] + MTH; else S[i][ = S[i-][j-] + MISMTH; sc = S[i-][ + GP; /* case #: i aligned to - */ if (sc > S[i][) S[i][ = sc; sc = S[i][j-] + GP; /* case #3: j aligned to - */ if (sc > S[i][) S[i][ = sc; /* The result (optimal alignment score) is in S[M][N]. */ Traceback ( see also above) /* llocate for our aligned strings (which will be..k-)*/ ax = malloc(sizeof(char) * (M+N+)); ay = malloc(sizeof(char) * (M+N+)); /* Start at M,N; work backwards */ i = M; j = N; k = ; while (i > && j > ) { /* ase : best was x_i aligned to x_j? */ if (i > && j > ) { /* watch out for boundaries */ sc = S[i-][j-]; if (x[i] == y[) sc += MTH; else sc += MISMTH; if (sc == S[i][) { /* residue i aligns to j: */ ax[ = x[i]; ay[ = y[; k++; /* record the alignment */ i--; j--; /* move to next cell in trace */ continue; 3

4 /* ase : best was x_i aligned to -? */ if (i > ) { /* watch out for boundaries */ if (S[i-][ + GP == S[i][) { ax[ = x[i]; ay[ = '-'; k++; /* record the alignment */ i--; /* move to next cell in trace */ continue; /* ase 3: best was - aligned to y_j? */ if (j > ) { /* watch out for boundaries */ if (S[i][j-] + GP == S[i][) { ax[ = '-'; ay[ = y[; k++; /* record the alignment */ j--; /* move to next cell in trace */ continue; /* Done with the traceback. * Now ax[..k-] and ay[..k-] are aligned strings, but backwards; * and k is our alignment length. * Reverse the strings and we're done. */ K = k; /* K = remember the alignment length */ for (k = ; k < K/; k++) { /* a sly, efficient in-place reversal by swapping pairs of chars: */ tmp = ax[k-k-]; ax[k-k-] = ax[; ax[ = tmp; tmp = ay[k-k-]; ay[k-k-] = ay[; ay[ = tmp; Output printf("sequence X: %s n", x+); printf("sequence Y: %s n", y+); printf("scoring system: %d for match; %d for mismatch; %d for gap n", MTH, MISMTH, GP); printf(" n"); printf("dynamic programming matrix: n"); for (j = -; j <= N; j++) printf(" %c ", j >? y[ : ' '); printf(" n"); for (i = ; i <= M; i++) { printf(" %c ", i >? x[i] : ' '); for (j = ; j <= N; j++) printf("%4d ", S[i][); printf(" n"); printf(" n"); printf("lignment: n"); printf("x: %s n", ax); printf("y: %s n", ay); printf("optimum alignment score: %d n n", S[M][N]); exit(); 4

5 G G T G 宿題 (3): G G T T G T T G G T T G T T G G T G G G T T G T T G G T -8 - G - -4 GGTTGTT GG- T- G = +7-8= GGTTGTT GGT- - G = +7-8 = - G T G Example: T T G T

6 G T G トレースバック T T G T 配列アラインメント Sequence alignment -Smith-Waterman (SW) algorithm- Local sequence alignment ローカルアライメント 配列アライメントから分かること 二つの DN の有機体は相似点と相違点を特定するのに重要である 以下のことを理解しやすくする : どのように 二つの有機体が進化するのか? ( 例 : 人間とチンパンジーの相違点 ) なぜ? どのように? インフルエンザのウィルスが変化するのか? なぜ病気になりやすい人がいるのか? ( 例 :DN 塩基の変化 ) 遺伝子をターゲットにした特定の薬を特定する DN 大域アラインメント DN ローカルアライメント 6

7 グローバルアライメントとローカルアライメント ローカル対グローバル ( 大域 ) 前回学んだアルゴリズム グローバルアライメント では あまり似ていない配列の比較は? ローカルアライメント ( よく似ているところだけそろえてみる ) 大域アラインメントは二本の配列に対して最適アラインメントを見つける DLGVFLDRYFQ DLGRTQN-DRYYQ 類似領域だけを見つける ローカルアラインメントは部分配列中に類似領域 ( ドメイン ) を見つける DLG DRYFQ DLG DRYYQ 例 : Global lignment --T ---GT -TTGT-GGGGG -GTGG-G TTGG-GTGT-T-TTG-----GTTTG T-GT-- Local lignment better alignment to find evolutionary-conserved segments ローカルアラインメント : 進化的保存部分列を発見するためのよりよいアラインメントである aattgccgccgtcgttttcaggtttgtgatc tccgtttgtggggacacgagcatgcagagac Smith-Waterman アルゴリズム 初期化,)=,,=, )= d はギャップのペナルティ ( スコア ) d= - 再帰式で 最大スコアを行列中から探し そこから があるところまでトレースバック ローカルアライメントのアルゴリズム New! i, j ) w( s[i], t[) max i, d j ) d w( x, x) x yのときw( x, y) New! New! 7

8 Example: グローバルアライメント 対角 対角 T T ?? - -4?? マッチノの重み + ミスマッチの重み - ギャップの重み T Example: ローカルアライメント 対角 対角 T T???? マッチの重み + ミスマッチの重み - ギャップの重み - T MX + = - - = = - 4 MX - = = -6 + = - MX + = - = - - = - MX = - - = - = - 例 : ローカルアライメント 例 : ローカルアライメント T G T G T T 8

9 例 : ローカルアライメント 例 : ローカルアライメント T G T + + T G T Local alignment 例 : ローカルアライメント 例 : ローカルアライメント G T G T T G T 初期化 T T G T T G T G T

10 例 : ローカルアライメント 初期化 行列中の最大スコア T T G T T G T G T 例 : ローカルアライメント 初期化 があるところに到達するまでトレースバック GT GT Local alignment T T G T T G T G T DN { G T { 4 文字 :DN 塩基 今まで 大域アラインメントローカルアライメント 今後 たんぱく質のアライメント 今まで 二本のDN 配列を比較して来ました 今から たんぱく質配列の場合を研究しましょう そんなに違わない たんぱく質 { 文字 : アミノ酸 大域アラインメント ローカルアライメント たんぱく質を比較するとき 単純スコアリングシステム ( 例えば : マッチの時は + ミスマッチの時ー ) は十分ではない

11 タンパク質の大域アラインメント : 動的計画法 New! DN 大域アラインメント DN ローカルアライメント i max, j ) M ( s[i], t[) i, d j ) d M Matrix G T G T M: アミノ酸の類似性行列 ( 例えば :BLOSUM6, or PM5 ) DN similarity matrix. たんぱく質大域アラインメンたんぱく質アライメント G T G T M Matrix PM5 BLOSUM6

12 Example: PM 5 Example : Blosum6 derived from block where the sequences share at least 6% identity 配列はすくなくとも 6% の類似度がある 類似なアミノ酸がより大きいスコアになる タンパク質のローカルアラインメント : Smith-Waterman (local alignment) i, j ) M ( s[i], t[) max i, d j ) d M: アミノ酸の類似性行列 (BLOSUM6, or PM5 for example) New! タンパク質の大域アラインメント : 動的計画法 S P E R E S H K E These values are copied from the PM5 matrix (see earlier slide) PM5, Gap =6 S -6 H - -8 K -4 E -3 S P E R E SPERE S-HKE

13 タンパク質の大域アラインメント : 動的計画法 S P E R E S H K E These values are copied from the PM5 matrix (see earlier slide) PM5, Gap =6 S P E R E S H K E SPERE S-HKE - タンパク質の大域アラインメント : 動的計画法 S P E R E S H K E These values are copied from the PM5 matrix (see earlier slide) PM5, Gap =6 S P E R E S H K E SPERE S-HKE - たんぱく質はアミノ酸で構成される 進化中のアミノ酸の置換は生化学的性質による PM 進化はたんぱく質を比較するのに重要である 3

14 PM: Percent of ccepted Mutation matrix Example: PM 5 PM 中の化学的類似性 : PM は似かよった化学的性質を持ったアミノ酸には より大きなスコアを与える cids & mides D,E,N,Q (sp, Glu, sn, Gln) Basic H,K,R (His, Lys, rg) romatic ( 芳香族 ) F,Y,W (Phe, Tyr, Trp) Hydrophilic ( 親水性 ),,G,P,S,T, (la, ys, Gly, Pro, Ser,Thr) Hydrophobic ( 疎水性 ) L,M,V,I (Ile, Leu, Met, Val) Similar amino acids have greater score 類似アミノ酸がより大きいスコアになる PM Matrices PM 5 Family of matrices PM 8, PM, PM 5 行列のグループ The numbers with PM matrices represent evolutionary distance. 数字は進化距離を示す Greater numbers are greater distances より大きな数はより進化距離 () スコアが 以上のペアアミノ酸は 両方とも似た機能を持つ () スコアが 以下のペアアミノ酸は 両方とも違う機能を持つ 4

15 Example: sp & Glu OO - OO - + H 3 N HH H + H 3 N HH H BLOSUM O O - spartate (sp, D) Score = 3 O HH O - Glutamate (Glu, E) BLOSUM: Blocks Substitution Matrix Based on BLOKS database BLOKS databaseによる ~ blocks from 5 families of related proteins 5の関係たんぱく質のグループの中に blocks がある Families of proteins with identical function たんぱく質のグループは同じ機能を持っている Blocks are short conserved patterns of 3-6 long without gaps ブロックは3から6ののギャップ無しの長さのパターンです BD BBD DBD..BBBB BBBDB.B D.DBDB BDB.DBBD BB BLOSUM Matrices () BLOSUMn は少なくとも類似度が n % の配列による () BLOSUM6 は BLOSUM45 よりさらに配列は似ている 5

16 Example : Blosum6 Entry (i) is greater than any entry (, ji. (i) エントリは最大です derived from block where the sequences share at least 6% identity すくなくとも 6% の類似度がある BLOSUM6 Symmetric matrix 対称行列 of x entries Dynamic programming PM5, Gap =6 S P E R E S H K E These values are copied from the PM5 matrix (see earlier slide) タンパク質のグローバルアライメント S -6 H - -8 K -4 E -3 S P E R E Global alignment Dynamic programming PM5, Gap =6 S P E R E S H K E These values are copied from the PM5 matrix (see earlier slide) S P E R E S H K E SPERE S-HKE - 6

17 This matrix is constructed using the PM5 matrix S P E R E S H K E Global alignment between proteins Dynamic Programming PM5, Gap =6 S P E R E S H K E SPERE S-HKE - 問題 : K P R T P E K end Protein Scoring Systems Sequence Sequence PTHPLSKTQILPEDLSEDLTI PTHPLGERIGLRLEEDFGM Blosum6 Scoring matrix S T P G N D.. 9 S - 4 T - 5 P G N D T:G = - T:T = 5 Score = 48 7

18 Dynamic programming for global alignment: Proteins. Dynamic programming for protein sequences: Smith-Waterman (local alignment) i max, j ) M ( s[i], t[) i, d j ) d d=- M: Similarity matrix for amino-acids (BLOSUM6, or PM5 for example) i, j ) M ( s[i], t[) max i, d j ) d d=- M: Similarity matrix for amino-acids (BLOSUM6, or PM5 for example) Global alignment Dynamic programming PM5, Gap =6 S P E R E S H K E These values are copied from the PM5 matrix (see earlier slide) S -6 H - -8 K -4 E -3 S P E R E Global alignment Dynamic programming PM5, Gap =6 S P E R E S H K E These values are copied from the PM5 matrix (see earlier slide) S P E R E S H K E SPERE S-HKE - 8

19 This matrix is constructed using the PM5 matrix S P E R E S H K E Global alignment between proteins Dynamic Programming PM5, Gap =6 S P E R E S H K E SPERE S-HKE - 問題 : K P R T P E K 9