使用 BWA 创建用于下游分析的输入 BAM 文件

Question

我有大肠杆菌的 Illumina 双端读段，用于创建草图程序集 (SPAdes)。我现在的任务是创建我将与 Pilon 一起使用的输入 BAM 文件——用于改进草图组件。

我决定使用此处的文档来使用 BWA：http://bio-bwa.sourceforge.net/bwa.shtml#3

计划创建参考基因组索引，创建比对，然后转换为 BAM 文件。

这是我用来索引参考的命令：

bwa index -p bwa_indices/B055 temp/spades/scaffolds.fasta

此命令输出以下文件：B055.amb B055.ann B055.bwt B055.pac B055.sa

下一步是生成比对文件——为此我使用了以下命令：

bwa aln -t 20 temp/spades/scaffolds.fasta temp/spades/corrected/B055_S5_R1_filtered_1P.fastq.00.0_0.cor.fastq.gz > bwa_indices/B055_R1.sai
#bwa aln -t 20 temp/spades/scaffolds.fasta temp/spades/corrected/B055_S5_R1_filtered_2P.fastq.00.0_0.cor.fastq.gz > bwa_indices/B055_R2.sai

在运行第一个命令之后，我收到以下输出：

[bwa_aln] 17bp reads: max_diff = 2
[bwa_aln] 38bp reads: max_diff = 3
[bwa_aln] 64bp reads: max_diff = 4
[bwa_aln] 93bp reads: max_diff = 5
[bwa_aln] 124bp reads: max_diff = 6
[bwa_aln] 157bp reads: max_diff = 7
[bwa_aln] 190bp reads: max_diff = 8
[bwa_aln] 225bp reads: max_diff = 9
[bwa_aln] fail to locate the index

最后一行让我有点烦恼。有一个输出文件 (B055_R1.sai)，但它是空的。

我可以清楚地看到在我的对齐命令中，没有引用任何以前创建的索引文件，但是当我查看文档时 (http://bio-bwa.sourceforge.net/bwa.shtm)，我看不到任何选项引用任何索引文件。谷歌搜索让我找到了一个网站，该网站说我需要将我的参考 fasta 文件放在与索引文件相同的目录中，我什至将我的程序集 fasta 文件草稿的名称从 scaffolds.fasta 更改为 B055.fasta——但无济于事。我还解压缩了 fastq.gz 文件并将扩展名从 fastq 更改为 fq——所有这些都没有成功。这些可能仍然是问题，但在我看来，在上次 bwa aln 调用中引用索引文件是最紧迫的问题。

谁能给我指明正确的方向？我用的BWA Version: 0.7.5a-r405 (我也安装了最新版本(Version: 0.7.12-r1039) 没有任何改进), CentOS 6.7, 34核内存充足

提前致谢。

Answer 1

基于 this discussion 我认为你应该使用 bwa mem 来执行对齐。

首先，生成您的索引文件：

# with this the index file will be in the same dir of your reference
bwa index path/to/your/index/scaffolds.fasta

然后进行比对：

bwa mem \
    path/to/your/index/scaffolds.fasta \
    /path/to/R1.fastq.gz \
    /path/to/R2.fastq.gz | samtools view -1 -bS - > youBamFile.bam

您可以将代码放在 bash 脚本中，或者像这样在同一行中键入所有内容：

bwa mem path/to/your/index/scaffolds.fasta /path/to/R1.fastq.gz /path/to/R2.fastq.gz | samtools view -1 -bS - > youBamFile.bam

解释：

bwa mem：比对双端读段并将其通过管道传输到 samtools 程序。
samtools view -1 -bS: 将您的 sam 文件排序并压缩为 bam 格式。

我没有测试命令，但可能会起作用。

Answer 2

根据其他论坛中某人的一些建议，我更改了我的文件的名称，以便它们在所有方面都保持一致。

mkdir -p bwa_indices
bwa index -p B055 -a is B055.fa
bwa aln -t 20 B055.fa ../temp/spades/corrected/B055_S5_R1_filtered_1P.fq > B055_R1.sai

但是，我仍然收到错误消息。我认为这是 outdated/incorrect 文档的问题。

文档 (http://bio-bwa.sourceforge.net/bwa.shtml#3) 有以下对齐方式（注意 in.db.fasta）：

aln     bwa aln [-n maxDiff] [-o maxGapO] [-e maxGapE] [-d nDelTail] [-i nIndelEnd] [-k maxSeedDiff] [-l seedLen] [-t nThrds] [-cRN] [-M misMsc] [-O gapOsc] [-E gapEsc] [-q trimQual] <in.db.fasta> <in.query.fq> > <out.sai>

我一直在使用以下内容（我尝试了 .fa 和 .fasta 扩展名）：

bwa aln -t 20 B055.fa B055_R1_1P.fq  > B055_R1.sai

我删除了 .fa 扩展名，它运行。

bwa aln -t 20 B055 B055_R1_1P.fq  > B055_R1.sai

使用 BWA 创建用于下游分析的输入 BAM 文件

Using BWA to create input BAM files for downstream analysis

bioinformatics

alignment

genome