Paired-End breseq Analysis Fails to Detect PCR-Verified Gene Knockout at Full Coverage

[khairnarsnehal940]

I have paired-end clonal whole-genome sequencing data from an E. coli gene knockout strain that I am analyzing using the breseq pipeline. I began with the standard breseq command (breseq -r MG1655.gbk -o breseq_output R1_001.fastq.gz R2_001.fastq.gz) for alignment and mutation calling. The analysis ran successfully, but the expected gene knockout was not detected as a deletion—neither in the missing coverage (MC) nor new junction (JC) sections of the output. When I examined the BAM file in IGV, the deletion region was clearly visible, with only about 5–6 reads aligning there, and the overall coverage across the genome was approximately 300×. However, when I ran the same dataset using the command breseq -l 180 -j 8 -r MG1655.gbk -o l180_breseq_output R1_001.fastq.gz R2_001.fastq.gz, the deletion was successfully identified, but breseq appeared to process only one of the paired-end files (R1). The deletion was confirmed experimentally by PCR, and all reads are of high quality with uniform read lengths. I would like to understand why breseq fails to call this verified deletion at full coverage but detects it when the data are subsampled, and why only one paired-end file seems to be used when the -l option is included. Is this behavior expected for paired-end datasets, or does it indicate a limitation or bug in how breseq handles high-coverage deletion detection?