2010年12月20日月曜日

何か作ってみた2

下の考えから、RAxML最新版の使ってブートストラップ解析するとき各サンプルにおいて通常のML探索と同様に複数系統樹を使って解析ができるようなスクリプトを書いたわけであるが、ぶっちゃけSeqboot使ってブートストラップデータを出力してそれをいちいち自分でコマンドを指定して解析するのを自動化しただけだけどね。

まあ、複数遺伝子でそれをやるときにちょっと操作がめんどくさかったけど。

ちなみに、複数遺伝子解析(塩基配列とアミノ酸配列の混合データ含む)では、一度配列をばらし、それぞれの遺伝子についてSeqbootでブートストラップデータを作成し、それを改めてがっちゃんこする方法をとった。

これは、結合データをそのままブートストラップリサンプリングしたら、各サンプルにおいて保持されるべき本来の遺伝子の配列長が変化してしまうからであり、多分基本的にはこれでいいと思う。ただ、RAxMLで複数遺伝子使ったブートストラップするときに具体的にどういう風にやってるか詳しく知らないので、細かなプロセスが違うかも知れんが。そこらへんはseqbootの方のオプション指定でなんとかなるかね。


また、RAxMLの方でブートストラップデータを作るコマンドはあるが、これで出力されるデータがあまりにもむちゃくちゃ(Aが100以上続くとか)なもので正直信用できなかったので、今回はseqbootを使うことにした。



また、このスクリプトでは各サンプルにおけるML探索を並列処理(RAxMLの並列版を使うだけ)がすること、またサンプルの解析自体も並列に行う(forkを使って例えば2サンプルずつ解析を進める)こと、これらを組み合わせて行う(前者に2スレッド、後者に4プロセス当てはめて8コア使うとか)ことができるようにした。

あと、これは単一遺伝子解析でも、複数遺伝子解析(Concatenate model or Separate model)でもできる。


まあ、こんくらいのスクリプトだったら探せば他の所にもありそうですが(笑)


とりあえず、塩基単一、アミノ酸単一、塩基+塩基、アミノ酸+アミノ酸、塩基+アミノ酸それぞれのデータを使って動作確認をし、ちゃんと動いたのでWikiの方にあげときます。需要があったらどうぞ。


全く使われないのもアレなので、何かいいテストデータ集めて色々検証してみますか。
ちょっとネタになるか考え中。。。

が、もう眠いのでもう寝る。疲れたぜ・・・

0 件のコメント: