Skip to content

最尤法

MEGAによる最尤法系統樹推定

MEGAは視覚的にわかりやすいGUI形式ですので,

特に解説は必要ないと思いますが,備忘録的に.



アライメントした配列を用意します.COI配列です.

Data→Phylogenetic Analysisを選択.



アミノ酸のコーディング領域かどうかを聞かれます.

今回はCOIなのでYesです.



うまくいけばこのようなウィンドウが表示されます.

MEGAの基本ウィンドウ(ここではMEGA 6.06(6140226))でPhylogeny→Construct/Test Maximum Likelihood Tree…を選択.



ここで解析のパラメーターを設定します.



Phylogeny Test :どのように系統を検定するかです.各枝の信頼性が得たい場合はここでBootstrap Methodを選択しましょう.

No. of Bootstrap Replications:Bootstrap検定を何回行うかです.100回でも構いませんが,1000回くらいはしたほうが無難です.2-3000回やってる論文もたまに見かけます.

Substitution Type:塩基(Nucleotide)かアミノ酸(Amino Acid)が選べます.選べないと困ります.

Model/Method:塩基(アミノ酸)置換モデルが選べます.

Rates among Sites:座位ごとの置換の頻度の違いを,ガンマ分布に基づいて分類するか否かを設定できます.

No of Discrete Gamma Categories:ガンマカテゴリ数を設定できます.

Gaps/Missing Data Treatment:ギャップの扱いの設定です.Complete deletionとすると,一つでもギャップがあるサイト(列)を解析からのぞけます.

ML Heuristic Method:最尤法系統樹の探索方法を選べます.系統樹探索方法についてはまたどこかで書くかもしれませんが,局所解がたくさんあるようなデータセットだと,ここでの方法選びは結構重要だと思います.いくつか試してみても良いかもしれません.

Initial Tree for ML:系統樹探索を行う際の,最初の出発点の系統樹の作成法を選びます.上記した通り,これも結構重要だと思います.デフォルトのNJでも問題はないと思いますが,もし既にそれらしい系統樹があるのであれば,それを設定する事もできます.

Blanch Swap Filter:系統樹探索の際の枝の入れ替えの大胆さを決めます.Strongにすると枝長の入れ替えがより消極的になり,解析時間は短くなりますが,考慮する系統樹は少なくなります.よりWeakにすると,枝長の入れ替えが大胆になり,解析時間は長くなりますが,考慮する系統樹が多くなるようです.



パラメーター設定が終わったら,Computeをクリックして解析開始です.

Progressが100%になるまで,気長に待ちましょう.



解析が終わると,このようにTree Explorerに系統樹が表示されます.

各枝の上にブートストラップ確率(2桁の整数)と枝の下に枝長(有理数)が示されています.



このExplorerで色々系統樹をいじれます.例えば特定の枝を選択して,左上のコマンドの中からPlace root on Branchを選ぶと



外群を指定できます.



他にもCompress/Expand Subtreeで,枝を一つにまとめたり,



Fit Tree to Screenで,ウィンドウ内に系統を収めたり,



Flip Subtreeで枝の上下を入れ替えたりできます.

良い時代になりました.



この他にもたくさんのコマンドがありますが,感覚的に理解できると思いますので,いろいろ試してみてください.



今日はここまで.


iqtreeによる最尤法系統樹推定

※もっと良いやり方をご存知の方はそっと教えてください.

iq treeは,パーティション分けが可能な最尤法系統樹推定ソフトです.

コマンドプロンプト対応で非常に使い勝手が良く,しかも解析速度がべらぼうに速いのが特徴です.

RAxMLがコマンドプロンプトに対応した今,やり方自体はRAxMLとほぼ同じです.
 

塩基配列データ作成(MEGA, SeaVeiwなど)
  • 事前準備として解析用フォルダをつくりましょう(iqとします).

アライメントソフトでの作業

  • アライメントを行ったセッションを.phy形式で出力する(iq.phyとしましょう).
  • パーティション分けとそれぞれのモデル指定のファイルを作る(iq.nexとしましょう).以下は例と解説(赤字)です.

#nexus

begin sets;

このコマンドで読み込みを開始します.

 

charset part1 = 1-430 431-1940\3 432-1940\3;

charset part2 = 433-1940\3;

charset part3 = 1941-2889;

ここで,モデルごとのパーティション分けを指定します.

ハイフンでつないだ配列が一つの領域で,それぞれをスペースで分けます.

アミノ酸指定領域の場合は,”日本円マーク”か”\”で各コドンをします.

ここでは,431番目から始まるアミノ酸指定の第一コドンと第二コドンは同じパーティションですが,

第三コドンは違うパーティションとしています.

 

charpartition mine = GTR+I+G:part1, TN93+I+G:part2, K2P+G: part3;

ここで,各パーティションごとのモデルを指定します.

よくモデルテストで見る形式をそのまま入力すればよいので楽です.

主な塩基置換モデルのリストはコチラ

 

end;

おしまい.


iqtreeを走らせる(コマンドプロンプト,iqtree)
  • 事前準備として,iqフォルダにiqtreeからDLしたiqtree.exeとiqtree-click.exeを入れておきます.

コマンドプロンプトでの作業

  • コマンドプロンプトを立ち上げ,iqフォルダまでのパスを通します.
  • iqtreeを走らせるためのコマンドを入力します.以下,例とコマンドの解説(赤字)です.

iqtree -s iq.phy -spp iq.nex -m TEST -bb 1000

“-s”で配列ファイルの読み込みを行います.

“-spp”でパーティション分けファイルの読み込み,

“-m” で検定方法を指定,

“-bb”でブートストラップ検定とその回数の指定を行います.


 

  • 解析が終わったら,iqフォルダ内に解析のログファイルや,系統樹ファイルiq.phy.treefileが生成されるはずです.
  • 参考までに,iqtreeでは他にも様々なコマンドがあります.コチラをご参照ください.

 

今日はここまで.


RAxMLによる最尤法系統樹推定

※もっと良いやり方をご存知の方はそっと教えてください.

塩基配列データ作成(MEGA, SeaView)

MEGAでの作業

  • MEGAでアライメントを行ったセッションを出力する.配列の結合を行う場合は,その後SeaViewで読み込むので,”.nexus” or “.nex”でも”.fasta” or “.fas”とかでもOK.行わない場合は”.phylip or “.phy”で出力.
  • ↑のファイルのOTU名は,それぞれのファイルで同じにしておく(”Ophi281_COI”とか”Ophi281_16S”じゃダメ).また,配列内のギャップ(”-“で指定する事が多い)と混同する事が多いため,OTU名には”-“は使わないほうがいい.

SeaViewでの作業(配列を結合する場合)

  • FileタブのOpen ***(***はファイルのフォーマット)を選び,結合させたい配列その1を選ぶ.
  • 同じようにその2,その3を選ぶ.それぞれ別ウィンドウで表示されればOKだが,不安な場合はFileタブのNew windowを選び,新しいウィンドウで配列を読み込めばOK.
  • FileタブのConcatenateを選ぶと,結合先の配列選ぶウィンドウが表示されるので,結合したい配列を選ぶ.by nameを選び,”OK”をクリック.
  • 上手くいけば,各OTUの後に配列が付加される.
  • これで準備完了.最終的な結合配列をFileタブのSave asからphylip形式で保存する(RA.phyとしましょう).
最尤法系統樹作成(RAxML)
  • 解析用の適当なフォルダ(RAとしましょう)を作り,その中にRA.phyを入れておく.
  • RA内にパーティション分け用のファイルを作っておく.例えば1-300 bpまでが16S, 301-1200までがCOIの場合は,以下のようなテキストデータファイルを作ります(partition.txtとしましょう).

DNA, 16S=1-300
DNA, COI=301-1200


  • RA内に”raxmlHPC.exe”を入れておく.これがRAxMLの起動ファイルです.これで準備完了.
  • コマンドプロンプトを起動(ウィンドウズメニューから探すか,ファイル検索で”cmd.exe”を検索)し,”cd 「RAのアドレス」”を入力してRA内に移動します.「RAのアドレス」はフォルダのアドレスバーの部分を左クリックすると”C:\Users\****\Desktop\研究\分子系統解析\RA”みたいなのが選択できるので,コピーしてコマンドプロンプト内で右クリックするとペーストできます.
  • “>”の後ろにRAxMLの起動コマンドを入力します.このコマンドは色々あるのですが,例えば

raxmlHPC -f a -m GTRGAMMA -p 12345 -x 12345 -# 1000 -s RA.phy -q partition.txt -n out

  • と入力すれば解析が始まり,GTRGAMMAモデルでブートストラップ1000回行った最尤法系統樹が得られます.
  • このコマンドなどについては井上潤さんのページにすごく詳しいです.
  • 出力されたファイルのうち,”RAxML_bipartitions.out”が系統樹のファイルです.TreeViewやFigtreeなどで閲覧して見てみましょう.

更新:2017年3月5日