簡(jiǎn)介
BWA傍睹，即Burrows-Wheeler-Alignment Tool隔盛。BWA 是一種能夠?qū)⒉町惗容^小的序列比對(duì)到一個(gè)較大的參考基因組上的軟件包。它有三個(gè)不同的算法：
- BWA-backtrack: 是用來比對(duì) Illumina 的序列的焰望，reads 長(zhǎng)度最長(zhǎng)能到 100bp骚亿。-
- BWA-SW: 用于比對(duì) long-read ，支持的長(zhǎng)度為 70bp-1Mbp熊赖；同時(shí)支持剪接性比對(duì)来屠。
- BWA-MEM: 推薦使用的算法，支持較長(zhǎng)的read長(zhǎng)度震鹉，同時(shí)支持剪接性比對(duì)（split alignments)俱笛，BWA-MEM是更新的算法，也更快传趾，更準(zhǔn)確迎膜，且 BWA-MEM 對(duì)于 70bp-100bp 的 Illumina 數(shù)據(jù)來說，效果也更好些浆兰。
對(duì)于上述三種算法磕仅，首先需要使用索引命令構(gòu)建參考基因組的索引，用于后面的比對(duì)簸呈。所以榕订，使用BWA整個(gè)比對(duì)過程主要分為兩步，第一步建索引蜕便，第二步使用BWA MEM進(jìn)行比對(duì)劫恒。
bwa安裝
直接使用mamba安裝即可
mamba install bwa
建立參考基因組索引
bwa index [ –p prefix ] [ –a algoType ] <in.db.fasta>
用法說明
-
-pSTR 輸出數(shù)據(jù)庫的前綴；【默認(rèn)和輸入的文件名一致轿腺，輸出的數(shù)據(jù)庫在其輸入文件所在的文件夾两嘴，并以該文件名為前綴∽蹇牵】 -
-a [algoType]構(gòu)建index的算法憔辫，有以下兩個(gè)選項(xiàng)：-
-a is是默認(rèn)的算法，雖然相對(duì)較快仿荆，但是需要較大的內(nèi)存螺垢，當(dāng)構(gòu)建的數(shù)據(jù)庫大于2GB的時(shí)候就不能正常工作了; -
-a bwtsw對(duì)于短的參考序列式不工作的，必須要大于等于10MB, 但能用于較大的基因組數(shù)據(jù)赖歌，比如人的全基因組。
-
示例代碼
bwa index GRCh38.primary_assembly.genome.fa
BWA-MEM 算法
該算法先使用 MEM(maximal exact matches) 進(jìn)行 seeding alignments功茴，再使用 SW(affine-gap Smith-Waterman) 算法進(jìn)行 seeds 的延伸庐冯。BWA–MEM 算法執(zhí)行局部比對(duì)和剪接性】泊可能會(huì)出現(xiàn) query 序列的多個(gè)不同的部位出現(xiàn)各自的最優(yōu)匹配展父，導(dǎo)致 reads 有多個(gè)最佳匹配位點(diǎn)返劲。有些軟件如 Picard’s markDuplicates 跟 mem 的這種剪接性比對(duì)不兼容,在這種情況下，可以使用 –M 選項(xiàng)來將 shorter split hits 標(biāo)記為次優(yōu)栖茉。
對(duì)應(yīng)的子命令為mem, 基本用法如下
bwa mem [options] ref.fa reads.fq [mates.fq]
參數(shù)說明
-
-t INT：線程數(shù)篮绿，默認(rèn)是1，增加線程數(shù)吕漂，會(huì)減少運(yùn)行時(shí)間亲配。 -
-M：將 shorter split hits 標(biāo)記為次優(yōu)，以兼容 Picard’s markDuplicates 軟件惶凝。 -
-p：若無此參數(shù)：輸入文件只有1個(gè)吼虎，則進(jìn)行單端比對(duì)；輸入文件有2個(gè)苍鲜，則作為paired reads進(jìn)行比對(duì)思灰。若有此參數(shù)：則僅以第1個(gè)文件作為輸入(會(huì)忽略第二個(gè)輸入序列文件，把第一個(gè)文件當(dāng)做單端測(cè)序的數(shù)據(jù)進(jìn)行比對(duì))混滔，該文件必須是read1.fq和read2.fa進(jìn)行reads交叉的數(shù)據(jù)洒疚。 -
-R STR： 完整的read group的頭部，可以用 '\t' 作為分隔符坯屿， 在輸出的SAM文件中被解釋為制表符TAB. read group 的ID油湖，會(huì)被添加到輸出文件的每一個(gè)read的頭部。 -
-T INT： 當(dāng)比對(duì)的分值比 INT 小時(shí)愿伴，不輸出該比對(duì)結(jié)果肺魁，這個(gè)參數(shù)只影響輸出的結(jié)果，不影響比對(duì)的過程隔节。 -
-a： 將所有的比對(duì)結(jié)果都輸出鹅经，包括 single-end 和 unpaired paired-end的 reads，但是這些比對(duì)的結(jié)果會(huì)被標(biāo)記為次優(yōu)怎诫。 -
-Y： 對(duì)數(shù)據(jù)進(jìn)行soft clipping, 當(dāng)錯(cuò)配或者gap數(shù)過多比對(duì)不上時(shí)瘾晃，會(huì)對(duì)序列進(jìn)行切除，這里的切除并只是在比對(duì)時(shí)去掉這部分序列幻妓，最終輸出結(jié)果中序列還是存在的蹦误，所以稱為soft clipping。
特別說明
- 如果 mates.fq 缺省肉津，且參數(shù) –p 未設(shè)定强胰，那么 reads.fq 被認(rèn)為是 single-end;
- 如果 mates.fq 存在，且參數(shù) –p 未設(shè)定妹沙，那么 mem 命令會(huì)認(rèn)為 read.fq 和 mates.fq 中的 i-th reads 組成一個(gè)read對(duì) (a read pair)偶洋，這個(gè)模式是常用的 paired-end mode。
- 如果參數(shù) –p 被設(shè)定距糖，那么玄窝，mem 命令會(huì)認(rèn)為 read.fq 中的 第 2i-th 和 第 (2i + 1)-th 的 reads 組成一個(gè) read 對(duì) （a read pair）牵寺，這種方式也被成為交錯(cuò)式的（interleaved paired-end)。 在這種情況下恩脂，即使有 mates.fq帽氓，也會(huì)被忽略。
示例代碼
# single end
bwa mem ref.fa reads.fq > mem-se.sam
# paired end
$ bwa mem ref.fa read1.fq read2.fq > mem-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
參考
http://www.reibang.com/p/19f58a07e6f4
https://blog.csdn.net/weixin_42192188/article/details/132286609
