BEAST
1. ソフトウエアの入手先
2. Infileの準備
BEASTには、infile作成専用のプログラム「BEAUTi」が同梱されている。
BEAUTiのタブに従って入力していけば、基本的に躓くところはない。
2.1. sequence input
BEAUTiではnexuxあるいはFASTA形式が使用可能。
手持ちのデータフォーマットが上記と異なる場合、以下のページ等でconvert。
http://www.hiv.lanl.gov/content/sequence/FORMAT_CONVERSION/form.html
2.2. STR input
Extended Bayesian Skyline Plot (EBSP)を使うことにより、STRデータに基づく集団サイズ推定が可能となった。
formatに関しては、testdataを用意したので参照していただきたい。
また、このperlスクリプトを使えば、GENPOPデータも簡単にBEAST形式に変換してくれる。
3. BEAUti
BEAUtiを立ち上げると、まず"Partitions"というタブが選択されている。
ここで、準備したinfileをimportする。
menu>File>Import data
以降、それぞれのTabについて、役割を記述する。
[Partitions tab]
Partitionsでは、それぞれの配列データセットに関して（たとえばCOI, 18S, Rhodopsin等）、
Site model (塩基置換モデルの設定)、Clock model (進化速度の設定)、Partition tree (系統樹作成)を個別に行うか、
あるいは単一の設定をすべてに適応するかを選択できる。
パーティションごと、あるいは単一の設定を用いたい場合、
partitionリスト上部にあるUnlink Subt. Model./Link Subt. Model等をクリックする。
[Taxa tab]
Taxa tabでは、データセットの中で特定のグルーピング (Taxon Set)を規定することができる。
たとえばデータセットが世界中から採集されたお魚のmtDNA配列だったとすると、
Set1は太平洋産、Set2はさらにその内群である日本近海産、といった具合。
それぞれのグループの単系統性を規定値として与えることができるが、
それゆえ事前に系統関係が判明している必要がある（BEASTの系統探索は弱いらしい）。
BSP/EBSPはインプットデータ全体に対してしか描けない。
つまり、root以降のすべてのLineageに対してDemographic modelが適応されるので、
特定の内群についてDemographic reconstructionを行いたいときは、データセットをその都度作り直してBEAUtiに読み込ませる必要がある。
[Tips]
省略
[Trait]
省略
[Sites]
Sites tabでは各パーティションの塩基置換モデルを設定する。
Partition tabでLink/Unlinkしたグループごとに設定することになるため、
パーティションごとに別のSite modelを適用したいときは、Unlink Subst. Modelsを選択しておく。
事前にMEGAやJ-Modeltest、Kakusan4等で最適なモデルを知っておく必要がある。
HKY/GTR/TN93/SDR06が使用可能。もちろん、G + I (3 codon position)も可。
[Clocks]
分子進化速度を各パーティションごとに設定できる。
Siteの場合同様、パーティションごとに別の分子進化速度/Clock Modelを適用したいときは、Unlink Clock Modelsを選択しておく。
Modelカラムでは、Strict/Lognormal/Exponential/RandomLocalの各ClockModelが選べるが、
BSP/EBSPを描く際にはStrict Clockを選択し、塩基置換速度を固定して推定を行うのが常。
各パーティションに関して文献値がある場合はよいが、そうでない場合は分からないパーティションの進化速度を推定する。
Estimateチェックボックスにチェックをつけると、尤もらしい進化速度を推定してくれる。
その際、Rateボックスに入っている値が事前情報となるため、それらしい値を入れておくとよい。
なお、塩基置換速度は年代推定のみならず、有効集団サイズの推定にも使われるため、
Rateの設定はDemographic reconstructionを行ううえで最も重要。
[Trees]
Treesタブにおいて、BEASTにどのような推定を行わせるかを決定する。
例えば、異なるDemographic scenario間でBayes factor を比較する場合、
Coalescent: Constant size/Expansion/Exponential/LogisticでそれぞれBEASTを走らせることになる。
また、BSP/EBSPの場合、それらを選択する。
特にEBSPでAutosomal nuclearデータを使う場合は、Partitionタブで事前にUnlinkしておき、
設定の際にはPloidy typeからAutosomal Nuclearを選ぶこと。
また、UPGMA starting treeを選択すると、樹形探索の効率が格段に向上する。
（この設定が不具合のもとになる場合もあるので、うまく走らない時はRandom starting treeに切り替える）
[Priors]
MCMCのカギとなる事前情報Priorをここで決定する。
確たる証拠があれば（例えばTransition/Transversion ratio=kappaはどう考えても2.0だ！とか）、
その値を反映させることができる。
が、そんな情報を持ち合わせていない場合はデフォルトの設定で問題ない。
一点、Clock rateの設定がinfinite uniformとなっているので、これを上限値100のUniformに変更する。
（恐らく赤字で示されているはず）
[MCMC]
ここでは、ベイズ推定にかける時間、ランの長さを決定する。
推定のチェーンが長い方が往々にして安定した結果が得られる（真実に近い結果とは限らないことに注意）
デフォルトでは10,000,000ステップのうち1,000回に1回データを取る、となっているが、
可能ならば各値を10倍に変更する。これで、計算にかかる時間も圧倒的に長くなるが、その分収束はよい。
一通り設定が終わったら、右下のGenerate BEAST File を選択し、BEAST .xmlファイルを作成する。
プログラムBEASTを起動し、出来上がったxmlファイルを選択して、走らせる。
（BEAST本体によくあるトラブル）
LnL = -infinity: 事前情報の与え方がおかしいか、設定に問題があるため、尤度が算出できない。たとえば、事前に与えた樹形とtaxa tabで与えたグルーピングに不整合がある、等。
.treeファイルや.logファイルが既に存在する場合も、エラーを吐いて止まる。中途半端に止まってやり直す場合、注意が必要。
（20170204追記)
LogNormal/Random Local等のClockを選択し、化石記録から分岐年代推定、と同時にBSPEBSP、というのも一つの手かもしれない。
Yule-Speciation priorを使って種内の分岐年代を推定するのは、各枝の集団サイズが常に一定であることを仮定している点において現実的でない。
明確な分岐年代のキーがあるならば、Bayesian Skyline proior を適用、MLtree固定、1/2点の年代校正で、こちらのほうが確からしいように思われる。
方法：RAxMLあるいはMrBayes等を用いて最尤樹/ベイズ樹を作成し、Best treeをNewickで手に入れる。
Taxa TabとNewickの不整合が無いよう注意を払いながら、BEAST .xmlファイルを作成。
この際、Operators tabから、Treeに関するOperatorをすべてoffにする。樹形探索でかえって尤度を下げるのを防ぐため。
出来上がった.xmlファイル中にNewick formatでStarting treeを規定する。
"id="startingTree" を含むセクションを探し出し、丸ごと削除、以下の内容に入れ替える。
<newick id="startingTree">
(ここにNewickを規定、文末のセミコロンを忘れずに)
</newick>
上書き保存した.xmlファイルを走らせてみて、確認。
ミトゲノムデータで走らせてみたところ、uniform priorを多用してもうまく走ったので、このほうがよさそう。
MCMCTREEから得られた結果とも整合性が高かったので、実用にはかなうか？