Troubleshooting E. coli Strain Identification: Identifying Transposon Insertion Sites with Paired-End Illumina Data

Question

I have a problem that needs your help.

I am working on a bacterial genomics project working with E. coli. The lab did illumina sequencing of many strains but I have some problems identifying the WT version (MG1655) from the recA, and ruvAB recombination gene mutants. The mutations were done by gene interruption 20+ years ago (I suppose using tn5). While they genetically verified that the genes do not work (also many years ago), they never checked what is the sequence of those mutations, for example, making a PCR and Sanger sequencing it. Which lead to my problem.

I am analyzing the sequencing results of the WT, recA, and ruvAB strains and the results makes more sense if the strains were in a different order: recA, ruvAB, and WT. Therefore, I need to verify that the strains are in a different order. Sequencing was pair-ended and I got about 20 million 50bp reads. Mapping the reads against the reference genome does not show a decrease on the depth of sequencing on any region either of the recA and ruvAB (because the gene got interrupted by a transposon) or at least not enough to consider that the transposon might be inserted on that gene. Which leaves me wondering how can I use the Illumina sequencing reads to verify which mutant does that genome has.

Is it possible to map against the reference genome and ask not to get the reads that map into it, but the ones that don't? That way, the ones I am left with should include the transposon. Then, maybe, assemble them de-novo. Check that the transposon is a Tn5 and what other nucleotides are left at the ends of the tn5 sequence. Finally, check those ends against either the recA or ruvAB gene.

Does this approach makes sense and can be done? If so, how? Does someone has a better idea on how to solve this problem?

Many thanks!

Answer 1

If you're asking about negative selection from mapped reads, yes, that can be done. There are many ways to do it, but the easiest is to use a mapper to create a BAM file containing all reads, then apply a filter on that file to extract the unmapped reads.

An example process would be as follows:

bowtie2-build reference.fasta reference.fasta
bowtie2 -x reference.fasta -U reads.fq.gz | samtools sort > reads_vs_ref.bam
samtools view -b -f 0x04 reads_vs_ref.bam > unmapped_reads_vs_ref.bam