2019 7/28 condaインストール追記
PBSIMはPacbioリードのシミュレーションを行うツール。ユーザーの持っているPacbioデータをもとにリードの長さやクオリティをシミュレートすることもできるため、実際の解析に適用しやすい。
インストール
GitHub - pfaucon/PBSIM-PacBio-Simulator: This is an updated mirror of the original PacBio Read Simulatorからソースコードをダウンロードしてビルドする。
git clone https://github.com/pfaucon/PBSIM-PacBio-Simulator.git
autoreconf -i
./configure
make
#bioconda (link)
conda install -c bioconda -y pbsim
ビルドが終わると、src/にpbsimができる。これにパスを通す。Pacbioのモデルデータとして/dataのファイルを参照する。
> pbsim
USAGE: pbsim [options] <reference>
<reference> FASTA format file.
[general options]
--prefix prefix of output files (sd).
--data-type data type. CLR or CCS (CLR).
--depth depth of coverage (CLR: 20.0, CCS: 50.0).
--length-min minimum length (100).
--length-max maximum length (CLR: 25000, CCS: 2500).
--accuracy-min minimum accuracy.
(CLR: 0.75, CCS: fixed as 0.75).
this option can be used only in case of CLR.
--accuracy-max maximum accuracy.
(CLR: 1.00, CCS: fixed as 1.00).
this option can be used only in case of CLR.
--difference-ratio ratio of differences. substitution:insertion:deletion.
each value up to 1000 (CLR: 10:60:30, CCS:6:21:73).
--seed for a pseudorandom number generator (Unix time).
[options of sampling-based simulation]
--sample-fastq FASTQ format file to sample.
--sample-profile-id sample-fastq (filtered) profile ID.
when using --sample-fastq, profile is stored.
'sample_profile_<ID>.fastq', and
'sample_profile_<ID>.stats' are created.
when not using --sample-fastq, profile is re-used.
Note that when profile is used, --length-min,max,
--accuracy-min,max would be the same as the profile.
[options of model-based simulation].
--model_qc model of quality code.
--length-mean mean of length model (CLR: 3000.0, CCS:450.0).
--length-sd standard deviation of length model.
(CLR: 2300.0, CCS: 170.0).
--accuracy-mean mean of accuracy model.
(CLR: 0.78, CCS: fixed as 0.98).
this option can be used only in case of CLR.
--accuracy-sd standard deviation of accuracy model.
(CLR: 0.02, CCS: fixed as 0.02).
this option can be used only in case of CLR.
root@3dff20ac5227:/data#
実行方法
リードのシミュレーション
pbsim --data-type CCS --depth 20 \
--model_qc PBSIM-PacBio-Simulator/data/model_qc_ccs ref.fa
--depth カバレッジ指定
--data-type CCSを指定。他にCLRがある。
--model_qc
エラー率とリード長はオーサーらの設定したモデルが反映される。CLR、CCSは大崎さんのブログを参照してください。
終わると、fasta一つに付き、3つのファイルが出力される。一つはシミュレートされたfastqファイルで、欲しいのはこれである。他2つはリファレンスから作られたエラー導入前のfasta配列とリファレンス、アライメント結果のファイルである。
生成されたfastqをもとの配列にアライメントしてみる。
bwa mem -x pacbio -t 12 ref.fa sd_0001.fastq > sd_0001.sam
sd_0001.fastqがPBSIMでジェネレートされたfastqファイルである。-x pacbio をつけ、
Pacbioのアライメントに最適化している（"-k17 -W40 -r10 -A1 -B1 -O1 -E1 -L0"をつけるのと同じになる）。
sd_0001.samをsamtoolsでbamにしてソートし、IGVで開く。下はCLRのシミュレートリード。
数kbpのリードがアライメントされているが、indel、ミスマッチはやや多い。
CCSもテストする。
pbsim --data-type CCS --depth 20 --model_qc PBSIM-PacBio-Simulator-master/data/model_qc_ccs ref.fa
bamにして、CLRのリードとエラーを比べてみる。
上がCCS、下が先ほどのCLRのリードである。一目見てCCSのエラーが減ってるのが分かる。
引用
PBSIM: PacBio reads simulator--toward accurate genome assembly
Ono Y, Asai K, Hamada M
Bioinformatics. 2013 Jan 1;29(1):119-21.
