macでインフォマティクス

macでインフォマティクス

HTS (NGS) 関連のインフォマティクス情報についてまとめています。

コード領域のアミノ酸配列を考えてマルチプルアライメントを行うMUCSE

 

塩基配列からコード領域のアミノ酸配列を予測してマルチプルアライメントを行う場合、従来はギャップやミスを補正せず全ての配列をアミノ酸に変換してアライメントを行なっていた。しかしこのような一義的に変換する方法だと、シーケンスエラーや擬遺伝子のstopコドンに引っ張られてアライメント精度がおかしくなる場合があった。

オーサーの発表したMUCSEは、stopコドンなどにペナルティを設定して、擬遺伝子を含むリストでも従来のNeedleman-Wunsch algorithmと同等の精度とスピードでマルチプルアライメントを行う方法論である。

 

 

公式サイト

http://bioweb.supagro.inra.fr/macse/index.php?menu=intro

チュートリアルファイル

http://bioweb.supagro.inra.fr/macse/index.php?menu=download

ロシア語のチュートリアルが多く分かりにくい。

 

インストール

brewで導入できる。

brew install MUCSE

 

実行方法

brewで導入したなら"macse -prog <commnad>"の形で使う。

 

alignSequences

元のfastaファイル(先頭5配列)。

f:id:kazumaxneo:20170909003555j:plain

 

このような短いコード領域の配列なら、デフォルトのパラメーターでランする。

macse -prog alignSequences -seq tmem184a.fasta

数分で解析は終わる。

塩基出力(tmem184a_macse_NT.fasta)。

f:id:kazumaxneo:20170909003655j:plain

 アミノ酸出力(tmem184a_macse_AA.fasta)。

f:id:kazumaxneo:20170909003711j:plain

アミノ酸ベースのマルチプルアライメントなので、塩基で確認するとアライメント精度が低い部位もある。

出力ファイル名は-out_NT-out_AAで変更できる。

 

 

 genetic codeがスタンダードではない生物の場合、コードを明示して解析する。例えばクロロプラストの遺伝子だと

macse -prog alignSequences -seq poacea_matk.fasta -gc_def 11

 NCBI genetic code(http://www.ncbi.nlm.nih.gov/Taxonomy/Utils/wprintgc.cgi

 

 

擬遺伝子も混じっている可能性がある(stopコドンが入る可能性がある)データセットで、擬遺伝子がわかっているなら、ファイルを別々に入力し、ペナルティスコアを個別に与える。この方が精度が高くなる。

macse -prog alignSequences -seq ENAM_genes.fasta -seq_lr ENAM_pseudos.fasta -fs_lr 10 -stop_lr 10
  • -seq "reliable" sequences
  • -seq_lr "less reliable" sequences

 上ではpseudo geneのペナルティスコアを-fs_lr 10 -stop_lr 10にし、信頼できる遺伝子のペナルティスコア(-fs-stop)をデフォルトでランしている。

 

 

NGSリードやcontigを使う場合

NGSのデータを直接使うとシーケンスエラーの可能性があるので、公式サイトではペナルティの与え方をはじめは-fs_lr 10 -stop_lr 15にして、後でチューニングするよう指示している。

macse -prog alignSequences -seq TMEM214.fasta -seq_lr TMEM214.fasta -fs_lr 10 -stop_lr 15

マルチプルアライメントはかなり計算が重くなる。100以上の配列があったりすると、計算が終わるのに非常に時間がかかる。

 

 

 splitAlignment

アライメント結果から、リストで与えた配列を別ファイルに出力する。

macse -prog splitAlignment -align tmem184a_NT.aln -subset primate_seqId.list

入力したリストは以下のようなものである。

user$ cat primate_seqId.list 

>Pongo

>Homo

>Gorilla

>Macaca

>Pan

 

 

他にもいくつかのコマンドがある。詳細は公式ページで確認してください。

 

2016年にversion2にアップデートされ、2018年夏にv2の論文が出ました。

https://academic.oup.com/mbe/advance-article/doi/10.1093/molbev/msy159/5079334

引用

MACSE: Multiple Alignment of Coding SEquences Accounting for Frameshifts and Stop Codons

Vincent Ranwez , Sébastien Harispe, Frédéric Delsuc, Emmanuel J. P. Douzery

PLoS One. 2011;6(9):e22594

 

Needleman-Wunsch algorithm

グローバルアライメント | グローバルアライメントを求める Needleman–Wunsch アルゴリズム