Error with BWA Mem input having multiple fastq files using cat and process substitution

Question

Running BWA Mem on a number of paired end fastq files using process substitution on the inputs results in this error:

[mem_sam_pe] paired reads have different names: "HWI-D00360:5:H814YADXX:1:2111:20284:33068", "HWI-D00360:5:H814YADXX:1:1213:17580:34847"

The command:

bwa mem -M ref.fa "<(cat U0a_CGATGT_L001_R1_004.fastq U0a_CGATGT_L001_R1_005.fastq U0a_CGATGT_L002_R1_001.fastq U0a_CGATGT_L002_R1_002.fastq U0a_CGATGT_L002_R1_003.fastq U0a_CGATGT_L001_R1_001.fastq U0a_CGATGT_L001_R1_002.fastq U0a_CGATGT_L002_R1_004.fastq U0a_CGATGT_L001_R1_003.fastq U0a_CGATGT_L002_R1_005.fastq)" "<(cat U0a_CGATGT_L001_R2_003.fastq U0a_CGATGT_L002_R2_005.fastq U0a_CGATGT_L001_R2_004.fastq U0a_CGATGT_L001_R2_005.fastq U0a_CGATGT_L002_R2_001.fastq U0a_CGATGT_L002_R2_002.fastq U0a_CGATGT_L001_R2_001.fastq U0a_CGATGT_L002_R2_003.fastq U0a_CGATGT_L001_R2_002.fastq U0a_CGATGT_L002_R2_004.fastq)" > out.sam

However the read (there were several) that threw this error does seem to be properly paired:

grep 'HWI-D00360:5:H814YADXX:1:2111:20284:33068' *.fastq
U0a_CGATGT_L001_R1_004.fastq:@HWI-D00360:5:H814YADXX:1:2111:20284:33068 1:N:0:CGATGT
U0a_CGATGT_L001_R2_004.fastq:@HWI-D00360:5:H814YADXX:1:2111:20284:33068 2:N:0:CGATGT

I can see that the concatenation of the input files is happening in a different order, so this is the next thing I'm going to try to fix, but has anyone else run into this?

(I probably just answered my own question but I'll post it anyway, that may or may not be the issue)

Answer 1

You are correct that the concatenation of input read files will confuse BWA. BWA depends on the input reads being in the same order in R1 and R2 files for paired-end alignment.

From bwa docs: If mates.fq is present, this command assumes the i-th read in reads.fq and the i-th read in mates.fq constitute a read pair. BWA will not search through FASTQ files trying to find matching mate IDs.

You will notice that your initial error message is reporting two different read IDs that are improperly paired. Presumably the mispaired read corresponds to the first read in the out-of-order file.

Answer 2

AS both the OP and @MaximilianPress suggested, the paired-end reads should be provided as 2 separate arguments to bwa. This is true whether the reads are in files or streamed. For the sake of completeness, I provide the corrected command here.

bwa mem -M ref_file.fa <(cat A.read1.fq B.read1.fq ...) <(cat A.read2.fq B.read2.fq ...)

The order of the files is important: if the order is A, B for read1, then it should be A, B and not B, A for read2.

EDIT:

(Thanks to @JohnMarshall for pointing this out.)

Avoid double quotes around process substitution. Quote the file names (if needed). The quotes around process substitution, as apparently intended, work only on some systems/versions¹, and seem like an undocumented feature to me.

DON'T:

bwa mem -M ref_file.fa "<(cat A.read1.fq B.read1.fq ...)" "<(cat A.read2.fq B.read2.fq ...)"

---

¹ Oddly enough, process substitution inside double quotes works at least in the case I tested before posting: bwa mem -M test1.fa "<(cat test_assembler_tiny_1.1.fastq)" "<(cat test_assembler_tiny_1.2.fastq)" This contradicts the answer here: https://unix.stackexchange.com/q/486019 . I used bwa 0.7.17-r1188, with both GNU bash, version 3.2.57(1)-release (x86_64-apple-darwin17) and zsh 5.7.1 (x86_64-apple-darwin17.7.0), on macOS High Sierra 10.13.6. From OP's error msg, it seems it worked for @glyph too.