2019,11,15 paftools.js call のコマンドエラー修正
2019 11/26 vatiantコールコマンドのエラー修正
2021 7/21 コマンド例追加
Minimap2には、PAFフォーマットのアライメントを処理する(java)スクリプトpaftools.jsが付属している。paftoolsを使うことで、 アセンブリをリファレンスゲノムにアラインメントしてバリアントをコールしたり、PAF/SAMからBEDなどのフォーマットに変換したりできる。
インストール
#minimap2インストールで導入できる。(link)
mamba install -c bioconda -y minimap2
#またはソースからビルドする(Github参照)
> paftools.js
# paftools.js
Usage: paftools.js <command> [arguments]
Commands:
view convert PAF to BLAST-like (for eyeballing) or MAF
splice2bed convert spliced alignment in PAF/SAM to BED12
sam2paf convert SAM to PAF
delta2paf convert MUMmer's delta to PAF
gff2bed convert GTF/GFF3 to BED12
stat collect basic mapping information in PAF/SAM
asmstat collect basic assembly information
asmgene evaluate gene completeness (EXPERIMENTAL)
liftover simplistic liftOver
call call variants from asm-to-ref alignment with the cs tag
bedcov compute the number of bases covered
version print paftools.js version
mapeval evaluate mapping accuracy using mason2/PBSIM-simulated FASTQ
mason2fq convert mason2-simulated SAM to FASTQ
pbsim2fq convert PBSIM-simulated MAF to FASTQ
junceval evaluate splice junction consistency with known annotations
ov-eval evaluate read overlap sensitivity using read-to-ref mapping
テストラン
git clone https://github.com/lh3/minimap2.git
cd minimap2/test/
1、paftools.js view - アラインメント表示
Usage: paftools.js view [options] <in.paf>
Options:
-f STR output format: aln (BLAST-like), maf or lastz-cigar [aln]
-l INT line length in BLAST-like output [80]
minimap2 --cs MT-human.fa MT-orang.fa | paftools.js view -
#実際にはlessなどに渡したり保存("> alignment")して見る
minimap2 --cs MT-human.fa MT-orang.fa | paftools.js view - |less
- --cs output the cs tag; STR is 'short' (if absent) or 'long' [none]
2、paftools.js stat - アラインメントstatistics
Usage: paftools.js stat [-l gapOutLen] <in.sam>|<in.paf>
minimap2 -c MT-human.fa MT-orang.fa | paftools.js stat -
- -c output CIGAR in PAF
出力
Number of mapped sequences: 1
Number of primary alignments: 1
Number of secondary alignments: 0
Number of primary alignments with >65535 CIGAR operations: 0
Number of bases in mapped sequences: 16499
Number of mapped bases: 16025
Number of insertions in [0,50): 51
Number of insertions in [50,100): 0
Number of insertions in [100,300): 0
Number of insertions in [300,400): 0
Number of insertions in [400,1000): 0
Number of insertions in [1000,inf): 0
Number of deletions in [0,50): 50
Number of deletions in [50,100): 0
Number of deletions in [100,300): 0
Number of deletions in [300,400): 0
Number of deletions in [400,1000): 0
Number of deletions in [1000,inf): 0
3、paftools.js call - バリアントコール
Usage: sort -k6,6 -k8,8n <with-cs.paf> | paftools.js call [options] -
Options:
-l INT min alignment length to compute coverage [10000]
-L INT min alignment length to call variants [50000]
-q INT min mapping quality [5]
-g INT short/long gap threshold (for statistics only) [50]
-f FILE reference sequences (enabling VCF output) [null]
-s NAME sample name in VCF header [sample]
#sortコマンドに渡しcoordinate sortしてからバリアントコールする
minimap2 -c --cs test/MT-human.fa test/MT-orang.fa \
| sort -k6,6 -k8,8n \
| paftools.js call -L 15000 -
#haploid assemblyのバリアントコール。de novoアセンブリ配列とリファレンスを比較し、結果を保存。
minimap2 -cx asm5 -t8 --cs ref.fa asm.fa \
| sort -k6,6 -k8,8n - \
| paftools.js call - > variant
#vcf output (ここではmap-ontプリセット)(100bp以上の違いのみ)
minimap2 -cx map-ont -t8 --cs ref.fa targets.fa | sort -k6,6 -k8,8n - | paftools.js call -f ref.fa -L 100 -l 100 - > variant.vcf
- -L min alignment length to call variants [50000]
- -x asm5|asm10|asm20 asm-to-ref mapping, for ~0.1/1/5% sequence divergence
- -f FILE reference sequences (enabling VCF output) [null]
出力最後尾
V MT_human 16528 16529 1 60 t c MT_orang 15983 15984 +
V MT_human 16531 16531 1 60 - g MT_orang 15986 15987 +
V MT_human 16533 16534 1 60 a c MT_orang 15989 15990 +
15993 reference bases covered by exactly one contig
2150 substitutions; ts/tv = 3.498
36 1bp deletions
34 1bp insertions
9 2bp deletions
5 2bp insertions
5 [3,50) deletions
12 [3,50) insertions
0 [50,1000) deletions
0 [50,1000) insertions
0 >=1000 deletions
0 >=1000 insertions
4、paftools.js liftover - リフトオーバー
Usage: paftools.js liftover [options] <aln.paf> <query.bed>
Options:
-q INT min mapping quality [5]
-l INT min alignment length [50000]
-d FLOAT max sequence divergence (>=1 to disable) [1]
クエリのbedファイルを指定する。
minimap2 -c MT-human.fa MT-orang.fa \
| paftools.js liftover -l 10000 - <(echo -e "MT_orang\t2000\t5000")
- -l min alignment length [50000]
5、paftools.js junceval - アノテーションとスプライスジャンクションを比較
Usage: paftools.js junceval [options] <gene.gtf> <aln.sam>
Options:
-l INT tolerance of junction positions (0 for exact) [0]
-p print overlapping introns
-e print erroreous overlapping introns
-c only consider alignments to /^(chr)?([0-9]+|X|Y)$/
gtfとマッピングのsam(pafではない)を指定する。
paftools.js junceval -e anno.gtf splice.sam > out.txt
6、paftools.js splice2bed - GTF/GFF3をBED12フォーマットに変換
Usage: paftools.js splice2bed [options] <in.paf>|<in.sam>
Options:
-m keep multiple mappings (SAM flag 0x100)
gtf/gffを指定する。
paftools.js splice2bed anno.gtf > anno.bed
7、paftools.js sam2paf - SAMをBED12フォーマットに変換
Usage: paftools.js sam2paf [options] <in.sam>
Options:
-p convert primary or supplementary alignments only
-L output the cs tag in the long form
paftools.js sam2paf input.sam > output.bed
引用
Minimap2: pairwise alignment for nucleotide sequences.
Li H
Bioinformatics. 2018 May 10.
関連