bwa序列比对

BWA主要用于把基因组测序数据比对到参考基因组上，由三种类似算法组成：BWA-backtrack，BWA-SW和BWA-MEM。首推BWA-MEM。

1. 建立索引

bwa index hg38.fa

注：索引构建时间比较长，个人电脑上一般一两个小时，依据配置而定。在同一目录下生成以hg38.fa为前缀的5个文件，输入不能是.gz压缩文件。索引建立好了，bwa比对时不再需要原始的hg38.fa 文件。

2. 用默认参数比对

# SE （single end）
bwa mem hg38.fa test.fq.gz > test_se.sam
bwa mem hg38.fa test.fq.gz -o test_se.sam
# hg38.fa 为索引前缀

# PE（paired end）
bwa mem ref.fa test_read1.fq test_read2.fq > test_pe.sam

3. bwa 和samtools 联用

bwa mem hg38.fa test_read1.fq test_read2.fq|samtools view -bS |samtools sort >test_pe.sorted.bam

4. 常用参数比对

bwa mem -t 10 -K 1000000000 -Y -M -R "@RG\tID:${ID}\tPL:Illumina\tLB:${library}\tSM:${sample}" hg38.fa read1.fq.gz read2.fq.gz >mem-pe.sam 2>./mem-pe.log
  
# ${ID}，${library},${sample}分别为测序id，建库类型以及样品号的变量值。

参数说明：

-t 整数，线程数，默认是1，增加线程数，会减少运行时间。

-K 每批处理的碱基数     

-Y 对数据进行soft clipping, 当错配或者gap数过多比对不上时，会对序列进行切除，这里的切除并只是在比对时去掉这部分序列，最终输出结果中序列还是存在的，所以称为soft clipping。 
-M  将 shorter split hits 标记为次优，以兼容 Picard’s markDuplicates 软件。

-R 字符串  read group header line such as '@RG\tID:foo\tSM:bar' ，默认没有。

2>  表示标准错误重定向到文件。

参考

Burrows-Wheeler Aligner