生信生信工具与原理

【比对软件】BWA使用

2020-03-18  本文已影响0人  生物信息与育种

BWA简介

BWA(Burrow-Wheeler Aligner),是一款将DNA序列mapping到参考基因组上的软件。有三个比对算法:BWA-backtrack,BWA-SW和BWA-MEM。详情可以看看李恒的解释。

BWA使用

1. 构建索引

reads比对前,需要对fasta构建FM-index。

index Usage:
      bwa index [ –p prefix ] [ –a algoType ] <in.db.fasta>
OPTIONS: 
      -p STR   输出数据库的前缀(默认和输入文件名一致)
      -a [is|bwtsw]   构建index的算法: is 是默认算法,相对较快,但需要较大内存,不适合超过2G的基因组。 bwtsw 对于基因组大于2G的,如人类基因组。

示例:

bwa index ref.fa
bwa index -a is ref.fa             #对小基因组建立index,速度快,内存消耗大
bwa index ref.fa -p genome
bwa index -a bwtsw ref.fa         #对于大基因组建立FM-Index

2. 比对

1). mem比对
mem Usage: bwa mem [options] ref.fa reads.fq [mates.fq]
最常用的比对方法,mem 使用的 MEMs(maximal exact matches) 进行seedling alignments, 再使用 SW(affine-gap Smith-Waterman 算法)进行seeds延伸,mem 进行局部比对和剪接性,因此,对于一条序列的不同区域可能会产生多种最优匹配结果, 这对于long reads 来说尤为重要。 有些软件如 Picard’s markDuplicates 跟 mem 的这种剪接性比对不兼容,在这种情况下,可以使用 –M 选项来将 shorter split hits 标记为次优。
常用参数:

-t   线程数,默认1。
-M   将 shorter split hits 标记为次优,以兼容 Picard’s markDuplicates 软件。
-p   若无此参数:输入文件只有1个,则进行单端比对;若输入文件有2个,则作为paired reads进行比对。若加入此参数:则仅以第1个文件作为输入(输入的文件若有2个,则忽略之),该文件必须是read1.fq和read2.fa进行reads交叉的数据。
-R   STR 完整的read标头,可以用 '\t' 作为分隔符, 在输出的SAM文件中被解释为制表符TAB. read group 的ID,会被添加到输出文件的每一个read的头部。

-T   INT   当比对的分值比 INT 小时,不输出该比对结果,这个参数只影响输出的结果,不影响比对的过程。-a 将所有的比对结果都输出,包括 single-end 和 unpaired paired-end的 reads,但是这些比对的结果会被标记为次优。

示例:

bwa mem ref.fa reads.fq > aln-se.sam
bwa mem ref.fa read1.fq read2.fq > aln-pe.sam
bwa mem -t 4 -M -R "\@RG\tID:{library}\tLB:{library}\tPL:Illumina\tPU:{sample}\tSM:{sample}\" ref.fa read1.fastq read2.fastq > mem-pe.sam 2> ./mem-pe.log

2). align/samse/sampe比对
对于single-read:

bwa aln [options] ref.fa read.fq > aln_sa.sai
bwa samse [options] ref.fa aln_sa.sai read.fq > aln-se.sam

对于pair-reads:

bwa aln [options] ref.fa read1.fq > aln_sa1.sai
bwa aln [options] ref.fa read2.fq > aln_sa2.sai
bwa sampe [options] ref.fa aln_sa1.sai aln_sa2.sai read1.fq read2.fq > aln-pe.sam

经典的 bwa 先使用 aln 命令将单独的 reads 比对到参考序列,再使用 samse 或 sampe 生成 sam 文件。用法:

bwa aln ref.fa short_read.fq > aln_sa.sai
bwa samse ref.fa aln_sa.sai short_read.fq > aln-se.sam
bwa sampe ref.fa aln_sa1.sai aln_sa2.sai read1.fq read2.fq > aln-pe.sam
bwa bwasw ref.fa long_read.fq > aln.sam

3). bwasw比对

bwa bwasw genome long_read.fq > aln.sam
bwa bwasw genome read1.fq read2.fq > aln-pe.sam

3. 人类示例

# 建索引
bwa index -a bwtsw -p hg19 hg19.fa 1>hg19.bwa_index.log 2>&1
# 比对
bwa mem -t 5 -M -R @RG\tID:KPGP-00001_L1\tSM:KPGP-00001_L1\tLB:WGS\tPL:Illumina ~/reference/index/bwa/hg19  KPGP-00001_L1_R1.fq.gz KPGP-00001_L1_R2.fq.gz 1>KPGP-00001_L1.sam 2>KPGP-00001_L1.bwa.align.log

Ref:https://www.jianshu.com/p/3b86615d647b
https://www.jianshu.com/p/1f6899d0fb71
https://www.bioinfo-scrounger.com/archives/181/

上一篇下一篇

猜你喜欢

热点阅读