|
|
|
|
RAxML (ラックスエムエル) は最尤法を用いて系統樹探索を行うプログラムです．fastDNAml の解析方法を改良したそうです．塩基配列とアミノ酸配列の解析が可能です．最尤法による系統樹探索では，2016 年 2 月の時点でもっとも速いプログラムと考えられています．ここでは，主に Mac のターミナルを用いた操作を解説します．
まずは，インターネット版 (RAxML Web-Servers: このページの下にあるコラム参照) を試してみるといいです．
疑問点などがあったら，RAxML Google Groups を参照するように，マニュアルで推奨されています．
|
|
簡単なチュートリアルを作りました．こちらです (2016 年 10 月)．
|
|
2014 年になって，ダウンロード方法が変更されました．RAxML のウェブサイトにある指示通り github というページに行ったら，右端にある「Download ZIP」というボタンを教えてください． ダウンロード&解凍されて「standard-RAxML-master」というディレクトリできます．13MB なので，やや時間がかかります．
現在のダウンロード方法はわかりやすいとは言えないので
，もしかしたら近いうちに変更されるかも知れません
(2014 年 9 月)．
|
|
解凍して作成された RAxML のフォルダーに入って以下を入力して下さい．
raxmlHPC 用:
make -f Makefile.gcc
raxmlHPC-PTHREADS 用:
make -f Makefile.SSE3.PTHREADS.gcc
コンパイル方法によって違うバージョンのプログラムが得られます．raxmlHPC-PTHREADS-SSE3 は通常の raxmlHPC-PTHREADS よりも速いそうです．詳しくはマニュアルをご覧下さい．こちらにある step-by-step tutorial も参考にしてください [2014 年 9 月]．
|
|
ターミナルから解凍してできたフォルダ「RAxML_WebExample5spp」に入ってください．intel Mac でコンパイルした raxmlHPC と raxmlHPC-PTHREADS-SSE3 いうアプリケーションが入っています．このままでも動くかもしれませんが，自分でコンパイルして作成されたアプリケーションに変更してください．script ファイルにある以下の行をコピーしてターミナルに入力してください．このとき，コンパイルして得られたアプリケーションの名前が raxmlHPC でなければ，適宜この部分を変更してください．数秒で解析が終わるはずです．
./raxmlHPC -f a -x 12345 -p 12345 -# 20 -m GTRGAMMA -s Coe5_500.PHYLIP -q part -o Scca -n outfile
なお，
sh script.sh
でも，script ファイルに書いてあるコマンドが実行されます (2014 年 9 月)．
RAxML_WebExample5spp.tar
|
|
ブートストラップと ML tree search を同時に行う
ブートストラップ解析を行って，同時に ML tree も探索してくれます．今後はこのやり方が RAxML ではスタンダードになると思います (2014 年)．
-f a -x 12345 -p 12345 -# 300 -m GTRGAMMAI -T 2 -s Seq.phy -q partition.txt -o OutGroups -n outfile
これで 300 回のブートストラップ解析を行ってくれます．-x と -p は Random Seed Number です．
raxmlHPC-PTHREADS -f a -x 12345 -p 12345 -# 300 -m GTRGAMMAI -s Coe12tr.PHYLIP -q part -o Scca -n outfile -T 8 > out_RAxML &
大規模データの解析
GTRCAT_FLOAT のように，FLOAT をつけると，memory の使用を減らすことができます．その分，解析速度は遅くなります．ver 7.2.3 のマニュアルですが，こちらを参考にしました．ver 7.4.4 でも使えるオプションです．
"-m GTRCAT_FLOAT" : Same as above but uses single-precision floating point arithemtics instead of double-precisionUsage only recommened for testing, the code will run slower, but can save almost 50% of memory.If you have problems with phylogenomic datasets and large memory requirements you may give it a shot.Keep in mind that numerical stability seems to be okay but needs further testing. [2010 年 3 月]
-e オプションによって，より厳密な ML tree search を行えるようです．私もよくわかっていませんが，tree search でより良い尤度を求めて hill climbing する際に，系統樹間の差異がどれぐらいになったら収束したかを決める値が e のようです．デフォルトでは 0.1 ですが，0.000001 を勧められたこともあります．
なお，-e オプションを使うと，RAxML_bipartitions.* などのファイルが作成されないことがあります．どのような場合なのかは，わかりません．
その他
・ アミノ酸と DNA の混合データ解析も可能なようです．partition 設定のファイルを変更します．
・ -f h オプションで SH テストができるようになったみたいです．
・ -M オプションで，partition ごとの branch length を出力してくれるようです．
・ V223 にもありましたが，V7 のマニュアルでも -c と -i の正しい設定を探すように書いてあります．デフォルトのままでも，MrBayes とほぼ同じ系統樹が得られるので，どれぐらい重要なのか不明です．
・ 解析に用いる種が少なすぎると動きません．4 種から解析可能なようです．
|
|
マルチコアの解析を簡単に行うことができるので，ぜひ PTHREAD version を使った方が良いです．MPI などによって並列解析を行うのは手続きが面倒ですが，PTHREADS はコンパイルの仕方がちょっと違うだけで，複数の CPU を同時に走らせることができます．Mac の場合は「この Mac について」などで搭載されている CPU の数がわかります．「-T 2」が CPU を 2 つ使うことを意味しています．上記にある例題を使って，-T 2 や -T 8 などとしてスピードを確認してから実際の解析を行った方が良いです． あやまって大きな数を -T に入れてしまうと，解析速度が非常に遅くなります．
raxmlHPC-PTHREADS -f a -x 12345 -p 12345 -# 300 -m GTRGAMMAI -s Coe12tr.PHYLIP -q part -o Scca -n outfile -T 2 > out_RAxML &
Mac の場合は，アプリケーション > ユーティリティ にあるアクティビティモニタ.app で，CPU が動いているか確認できます．私の Mac はプロセッサが「 3.06 Ghz Intel Core 2 Duo」となっていますが，raxmlHPC-PTHREADS-SSE3 を -T 2 で動かすと，以下のようになりました [2013 年 3 月]．
|
|
|
my $raxmlNumThreads = 2;
Overall execution time for full ML analysis: 6.456854 secs or 0.001794 hours or 0.000075 days
my $raxmlNumThreads = 4;
Overall execution time for full ML analysis: 5.048122 secs or 0.001402 hours or 0.000058 days
my $raxmlNumThreads = 6;
Overall execution time for full ML analysis: 4.650692 secs or 0.001292 hours or 0.000054 days
my $raxmlNumThreads = 8;
Overall execution time for full ML analysis: 7.548798 secs or 0.002097 hours or 0.000087 days
|
|
|
Infile で用いる OTU 名
以前 10 文字以内でないと正常に作動しないように書きましたが，そのようなことはありませんでした．さらに，OTU 名に「_」が入っていても問題なく動くようです．RAxML に限らず様々なソフトウェアを用いる場合は，OTU 名に気を付ける必要があります．アルファベットと数字のみ，にすると，問題になることが少ないと思います．複数のファイルで OTU 名の変更を一気に行う場合は，こちらを参考にして下さい．
Match First も読み込む
infile で MacClade の「Match First.」を用いても，問題なく動きます．「Match First..」とは，2 OTU 目以降に同じ配列を「.」で表示する方法です．
|
|
|
パーティションの設定は「-q」コマンドです．ここでは V7 のパーティション設定を紹介しています．
./raxmlHPC -f d -m GTRMIX -s infile.phy -q partition_12 -# 100 -n MultipleOriginal
「partition_12」という名前のファイルにパーティションの設定を書いています．内容は以下です．
DNA, gene1=1-200, 800-1000
DNA, gene2=201-799
アミノ酸の解析では，パーティションごとに異なるモデル (rate matrix) を指定することができます．「F」オプションを付けると，データからアミノ酸組成を推定 (最尤推定ではなさそう) して解析に取り入れます．
JTT, gene1 = 1-500
WAGF, gene2 = 501-800
DNA, gene1codon1=1-10608\3
DNA, gene1codon2=2-10608\3
DNA, gene1codon3=3-10608\3
DNA, gene1codon1=1-10608\2
DNA, gene1codon2=2-10608\2
としましたが，問題なく動いているようです．
DNA とアミノ酸の両方が入ったデータも解析可能です．
DNA, gene1 = 1-500
WAGF, gene2 = 501-800
|
|
|
Monophyly の constraint は「-g」オプションです．
./raxmlHPC -f d -m GTRMIX -s test9.phylip -q part_12tr -g 2mono -n 2monoRes
「2mono」という名前のファイルに制約を書いています．
ここではコマンドが 2 行に渡っていますが，実際はリターンを入れずに 1 行で書いてください．他のコマンドも同様です．
((PolErpCal ,PolPoSeSe ,AciPolSpa ,AciAciTra ,CoeLaChTa ,CoeLatMen ,(PolPolOrn,AciScaAlb)), ScyScyCan);
根幹の分岐は 2 分岐でも大丈夫ですが，たまに 3 分岐にするよう error message が出て解析が進まないことがあります．
以前は r オプションの記述と混同して，外群を外した制約を与えていましたが，g オプションでは外群も制約に書き込んだ方が良いです．g オプションで外群を外した制約を行うと，外群直下の分岐に制約を与えても結果に反映されないことがありました．
系統関係の制約は「-r」オプションで与えられます．これはマニュアルで backbone tree と呼ばれていて，この tree に含まれていない taxa (data にはある) が tree に付け加えられていきます．このため constraint となる tree は完全二分岐で，data set よりも少ない OTU 数であった方が良いみたいです．
|
|
|
スピードと精度の観点からも，ブートストラップ解析と最尤樹推定を同時に行う -f a コマンドを行うのが効率が良いと思います．詳しくは，上記 「ブートストラップと ML tree search を同時に行う」参照してください．コマンドは以下です [2010 年 12 月]．
-f a -x 12345 -p 12345 -# 300 -m GTRGAMMAI -T 2 -s Seq.phy -q partition.txt -o OutGroups -n outfile
|
DNA 配列の解析
./raxmlHPC -f d -m GTRMIX -s infile.phy -q partition -# 100 -n MultipleOriginal
-#100 のように，できるだけ多くの解析を行います．作者はマニュアルで，(おそらく 1000 以上の OTU からなるデータセットを用いて)，-# 200 と設定した，と書いています．
計算が始まると，自動的に最節約樹が推定され (RAxML_parsimonyTree.NondMultipleOriginal_0)，これを starting tree として最尤法に基づいて系統樹探索を行っているようです．
-d オプションを用いると，starting tree がまったくの random tree になるようです(RAxML_randomTree.2MultipleOriginal_0)．
計算が終了すると，「RAxML_info.〜」というフォルダ内に，どの系統樹の尤度が最も高いか書いてあります．
アミノ酸配列の解析
モデルは mtREV+F+Gamma (F はデータセットからアミノ酸組成を推定するオプション) を使用しています．
./raxmlHPC -f d -m PROTGAMMAMTREVF -s ActCYB8.PHYLIP -q partition -# 100 -n ActCYB.out
|
|
|
./raxmlHPC -f d -m GTRCAT -s Nuc12Act57.PHYLIP -q partition_12 -# 100 -b 12345 -n bsout
-f d
系統樹探索のアルゴリズム
-m GTRCAT
モデルの選択．ブートストラップ解析では CAT-based のモデル (尤度を最大化するのではなく，ベストの系統樹を探索することに的を絞ったモデル) がマニュアルでは推奨されています．作者の Stamatakis さんに問い合わせたところ，BP 解析をCAT ではなく MIX モデルで行っても大きな違いはないということでしたが，やはり計算時間は尤度最大化と樹長の計算をする分だけ (10〜20%) 時間がかかる，ということでした．
以上のように書きましたが，GTRMIX で 8 時間かかった解析が，GTRCAT では 4 時間になるようです (宮正樹先生に教えて頂きました)．このためマニュアルに従って，ブートストラップ解析は GTRCAT で行った方がよいと思います．
-s Act57CYBR.PHYLIP
シーケンスファイル
-# 100 b 12345
100 回反復のブートストラップ解析を，12345 という seed number から始める．
-n MultipleBootstrap
Outfile
アミノ酸配列の場合
./raxmlHPC -f d -m PROTCATMTMAMF -s Act57CYBRM.PHYLIP -# 100 -b 12345 -n bsout
ここでは mtMAM+F (アミノ酸組成をデータから推定) モデルを用いています．
ML tree に BS 値をのせる
最尤樹はあらかじめ推定しておく必要があります．ML tree は樹長がないとうまく動きませんでした．BS 値をのせた後でも樹長は保存されています．
./raxmlHPC -f b -m GTRCAT -s Coe5_500.PHYLIP -z BS.trees -t ML.tree -o Scca -n BS_TREE
RAxML_bipartitions.BS_TREE に BS 付き系統樹が得られています．
この機能はわりと便利で，MrBayes の consensus tree などに RAxML で得られた BS 値をのせることもできます [2014 年 12 月]．
RAxML_Consensus_fol.tar
* 以前は以下のように考えていましたが，ML tree および BS trees の根幹が二分岐でなくても，解析が動きました [2010 年 9 月]．
BS trees と ML tree は根幹も含めて完全 2 分岐にしておいた方が良さそうです．なお，複数の tree について根幹の 3 分岐を 2 分岐にする場合はこちらを参考にしてください．
Consensus tree の作り方
-J オプションでできます．詳しくはマニュアルを見てください (2014 年 12 月)．
RAxML_Consensus2_fol.tar.gz
|
|
|
raxmlHPC725 -f e -m GTRGAMMA -t brown.tree -s brown1.nuc -n outfile
-f e: 樹長を推定します．
例題．
raxmlHPC-PTHREADS-SSE3 -f e -e 0.0000001 -m PROTGAMMAWAGF_FLOAT -t mammal36s.Unrooted.tre -s mammal36sAA.txt -n out -T 26
アミノ酸配列からなる巨大なデータでは，メモリが大きくなりすぎて解析できないことがあります．Alexis に聞いたところ，上記のように「_FLOAT」をモデルにつけることによって single prevision にすれば，メモリの使用を半分にすることができます．single prevision の意味は私にもよくわかりませんが，こちらに多少説明があります．ここでは 26 コアを用いて 94GB のメモリを使ってようやく解析が動きました [2010 年 3 月末]．
|
|
RAxML の Web-Servers は，よく整備されていて使いやすいです．ミトコンドリアゲノム配列程度 (あるいはそれ以下) の長さの配列を用いるのであれば，特別な解析を行わない限り，Web_Servers による解析で十分だと思います．
Swiss Institute がサービスを提供しています．塩基とアミノ酸の両方で解析が可能です．
BS 解析を行った後に ML tree を求め，樹長を ML tree に載せた tree を自動で出してくれます．ML の starting tree を代えた tree 探索の回数を設定することはできませんが，BS の解析は指定できます．解析スタートと終了時点で，登録したアドレスにメールが届きます．
さすがに 8 CPU を用いた PTHREADS 解析などにはかないませんが，かなり解析が速くて，詳細な設定も論文出版用の解析に十分対応しています．
例題はこちらです．
|
|
| Perl script: cDNA 配列を 1st と 2nd に分ける |
|
| raxml_1st2nd.tar.gz (2016 年 3 月) |
|
|
まだ作成途中です．
RAxML_perSiteLn.tar.gz
|
|
v728 で動いたファイルが，v813 では動きませんでした．
[jun-inoue:raxmlSgeTest]$ raxmlHPC-PTHREADS-SSE3 -f a -x 12345 -p 12345 -# 100 -m GTRGAMMA -s 010_sequenceFileTestDir/ENSP00000216034.txt -q 010_partitionFileTestDir/ENSP00000216034.txt -o Drosophila_FBpp0073033_mge -n ENSP00000216034.txt -T 4
This is the RAxML Master Pthread
This is RAxML Worker Pthread Number: 1
This is RAxML Worker Pthread Number: 2
This is RAxML Worker Pthread Number: 3
This is RAxML version 8.1.3 released by Alexandros Stamatakis on August 13 2014.
RAxML was called as follows:
raxmlHPC-PTHREADS-SSE3 -f a -x 12345 -p 12345 -# 100 -m GTRGAMMA -s 010_sequenceFileTestDir/ENSP00000216034.txt -q 010_partitionFileTestDir/ENSP00000216034.txt -o Drosophila_FBpp0073033_mge -n ENSP00000216034.txt -T 4
Time for BS model parameter optimization 0.013184
raxmlHPC-PTHREADS-SSE3: evaluateGenericSpecial.c:3433: evaluateIterative: Assertion `partitionLikelihood < 0.0' failed.
Aborted (core dumped)
|
|
rate category 補正なしの解析
作者に問い合わせたところ，ソースコードから関数 optimizeRateCategories() をすべて取り除いて compile する必要があるそうです．この関数は複数のファイルに存在するので，BBEdit の検索 (コマンド F) で others を選択して (複数のファイルを検索してくれる) 探すと良いでしょう．とくに optimizeModel.c というファイルでは
void optimizeRateCategories(tree *tr, int categorized, int _maxCategories, analdef *adef)
{
int i, k;
double initialRate, initialLikelihood, v, leftRate, rightRate, leftLH, rightLH, temp, wtemp;
double lower_spacing, upper_spacing;
．．．．．
}
のように void の行とこれに続く {．．．．} のすべてを削りました．作者に確認はしていませんが，ブートストラップ解析をガンマ補正した場合と比較したら違いが出たので，問題なく動いていると思います．ちなみにソースコードを書き換えずに「-c 1」(カテゴリーを 1 つとする) としてはどうかと作者にたずねたら，どうなるかわからないと言うことでした．
|
どれぐらい速いか
以下の論文では MrBayes と計算速度を比較して，その速さを強調しています． |
Stamatakis A, Ludwig T, Meier H. RAxML-II: a program for sequential, parallel and distributed inference of large phylogenetic. CONCURRENCY AND COMPUTATION-PRACTICE & EXPERIENCE 17 (14): 1705-1723 Sp. Iss. SI DEC 10 2005
|
とくにアミノ酸の解析速度が MrBayes など他のプログラムに比べると非常に速いです．ミトコンドリアゲノムの解析では RAxML で精度が劣ることはないです．50 OTU 以上を用いて，広く一般的に認められている脊椎動物における高次の系統関係を再現できるか調べたところ，MrBayes と同じような答えを返しています [2007 年]．
|
|
