pyRAD

2016年03月17日 更新

pyRADとは

RAD-seqやそれに類似するデータ解析用のpython2のソフトウェア.
RAD, ddRAD, PE-ddRAD, GBS, PE-GBS, EzRAD, PR-EzRAD, 2B-RAD, nextRADなど.

インストール

これまでpythonを使ったことのない人が新たにMacに導入することを想定します.
このために作ったPYRADディレクトリ内で扱うことを想定します.pyRAD以外にも使える部分がかなりあるので,保存するディレクトリをいじるなどは各自必要に応じて試してみてください.

まずは,pythonの環境を整えて必要なライブラリをインストール.
pip pythonインストーラ
% sudo python2.7 -m ensurepip

virtualenv 仮想環境構築ライブラリ
% sudo pip2.7 install virtualenv

virtualenv実行 仮想環境で独立にpythonを操作します.
% virtualenv -p python2.7 --no-site-packages py27
% source py27/bin/activate

numpy,scipy 数値演算・科学演算ライブラリ
% pip install numpy
% pip install scipy
virtualenvは必須ではないですが,python2と3は互換性がないので,できるだけ区別したほうが良いという僕の思想です.

続いて,アラインメントソフトMUSCLEをインストール.
ここからMac用のMUSCLEをDL.解凍します.

最後にVSEARCHをインストール.
% git clone https://github.com/torognes/vsearch
% cd vsearch
% ./autogen.sh
% ./configure
% make
% sudo make install
./autogen.shで
./autogen.sh: line 2: autoreconf: command not found
と出る場合は,automakeautoconfを入れてください.

ようやく,pyRADです.
pyRADはgit-hubからDLできます.
% git clone https://github.com/dereneaton/pyrad.git
これで準備完了.

使い方

pyradを使うときは常に仮想環境で行います.もし(py27)のような環境名がターミナル上に表示されていなければ仮想環境を実行します.
% source py27/bin/activate
ちなみに、virtualenvをオフにするには
% deactivate
バーコードファイルがある場合は,サンプル名と配列をタブ区切りで並べた形式にします.
sample1    ACAGG
sample2    ATTCA
sample3    CGGCATA
sample4    AAGAACA
pyRAD用パラメータファイルを作ります.
% python pyrad/pyrad/pyRAD.py -n
ファイル内のパラメータの値を編集します.詳細は下に記述.ファイル内や本家のチュートリアルも参照.1から14は必須,それ以降はオプションです.
準備ができたら実行します.
pyRADの全プロセスを実行する.
% python pyrad/pyrad/pyRAD.py -p params.txt

pyRADの一部のステップのみ実行する.
% python pyrad/pyrad/pyRAD.py -p params.txt -s 1
% python pyrad/pyrad/pyRAD.py -p params.txt -s 234567
これで,時間をおけば様々な形式のアウトプットが出てきます.

パラメータ詳細

以下の1-14は必須のパラメータ.
## 1. 作業ディレクトリ.
## 2. バーコードごとに分離されていないfastq.既に分離されていれば空欄で良い.
## 3. バーコードファイル.既に分離されていれば空欄で良い.
## 4. vsearchのコマンド(とパス).
## 5. muscleのコマンド(とパス).
## 6. 制限酵素の認識配列.複数用いた場合はコンマで区切る.幅広い制限酵素,レアな制限酵素の順.
## 7. 並列に実行するジョブの数.Macならコア数と同じにすると速い.サンプル数より多くしてもほとんど意味はない.
## 8. 最小カバレッジ.これよりリード数の少ない座位は捨てて質を保つ.典型的には5とか10とか.
## 9. リード中にいくつまでNを許すか.
## 10. サンプル内,サンプル間でのクラスタリングにおける類似性の閾値.
## 11. データのタイプ.rad, pairddrad, gbsなどがあります.
## 12. 最小サンプル数.外群を除いてこのサンプル数以上の座位のみ用いる.集団遺伝学では用いたサンプル数の75%とすることが多い.
## 13. 座位内で多型部位を共有する最大個体数.同じ部位で同じヘテロを多くのサンプルが共有する場合,ヘテロではなくパラログかもしれない.
## 14. アウトプットファイルの頭につける名前.お好きにどうぞ.

以下の15-37はオプション.必要に応じて指定するが,なくても動く. ## 15. サブセットがある場合に指定する."SUB1"と書くと,SUB1から始まるサンプルがサブセットとなる.
## 16. 外群の指定.outg1,outg2のように書く.12行目の最小サンプル数はこの外群を除く.
## 17. 取り除きたいサンプルがある場合,ここで指定する.sample3,sample4のように書く.
## 18. バーコードで分離されたfastqファイル.パラメータ2と3が空欄の場合,必須.
## 19.
## 20.
## 21.
## 22.
## 23.
## 24.
## 25.
## 26.
## 27.
## 28. 乱数のシードの指定.
## 29.
## 30. 出力のフォーマット.pnasvutmkg*が使える.全て出力する*がお勧め.
## 31.
## 32.
## 33.
## 34.
## 35.
## 36.
## 37. vsearchで使うジョブあたりスレッド数.

ステップ詳細

Step1. バーコードファイルに従ってfastqの配列を分離.fastq/を作って出力.
Step2.
Step3.
Step4.
Step5.
Step6.
Step7.

その他注意点

パスやコマンドの指定は,絶対パスでも相対パスでも可.
*はワイルドカードとして同じディレクトリ内の複数ファイルを指定.
paired-endの場合,ペアとなるfastqファイルの名前が_R1_と_R2_以外同じとする.シングルエンドならどんなファイル名でも良い.
その他パラメータやステップ,アウトプットの詳細は暇があれば追記するかもしれません.

出力のvcfファイルをADMIXTUREで解析する場合、PLINKによってvcfファイルを変換する必要があります.pyRADで出力されるvcfファイルの#CHROM POS ID ...の行にはタブだけでなくスペースを含む部分があり,それが原因でサンプル名が見つからないエラーが出ます.とりあえずスペースを消去してタブのみにすれば解決します.

参考文献

pyRADのチュートリアル
ensurepipについて
virtualenvについて