pindel not picking up on a lab-made gene knockout deletion

Question

[pdm37]

I have used pindel (version 0.2.4t) to look for indels in two strains that had a gene knocked out and replaced with an antibiotic resistance gene. The two mutants grew despite the knockout and we were looking for other mutations that may have led to this phenotype. Pindel didn't identify the knockout as either a deletion or an insertion, or any of the other outputs for that matter. Should it have? I have identify another deletion in another strain with pindel (not a deletion/insertion combo as in my two strains, just a regular deletion), so I'm pretty sure I'm running pindel correctly. Any thoughts?

I ran it with default settings, using wildtype as reference (fully assembled, closed genome).

Thank you.

[Chris Miller]

What's are the sizes of the deletion and insertion? Large events, especially, are hard to detect with pindel (and many other tools), because read mapping is often poor.

[pdm37]

the deletion is 199 bp, the inserted sequence is 251bp. I can see the deletion in the bam file and the insertion comes up in the unmapped reads

Answer 1

[Brian Bushnell]

I suggest you look at the mapped reads in IGV and inspect the region manually. You will either see:

Normal coverage;
Reduced coverage;
No coverage;
or Reads aligned with deletions.

The ability to identify deletions depends a lot on the aligner.

[genya35]

I'm experiencing a similar problem with detecting 19 pb indels with Pindel in (EGFR exon 19 E746_A750del) in Ion Torrent amplicon base assay. I've tested with several patient samples and the score for that deletion is 10-16 for some reason. Freebayes variant caller picks them up without a problem. When I look in IGV the reads are clean and abundant. I'm using a target approach and run Pindel on small bams with chr specific file as reference. It seems to work well with other regions. Is there anything that can be done to improve the performance? Thanks

[Brian Bushnell]

To be honest... I can't think of a good solution here other than not using Pindel. I encourage you to send this as feedback to the Pindel developers, and use Freebayes in the interim.

[pdm37]

thank you for the suggestion! The alignment has a gap in that region, or rather all the bp in the reads that span it are crossed out (due to mismatch, which is due to the insertion, I would presume). I used BWA mem as the aligner.

[Brian Bushnell]

You might try mapping with BBMap and calling variants with BBMap's callvariants.sh. BBMap is very tolerant of long indels.

[genya35]

@Brian Please explain how to run bbmap.sh for BBMap to index ta reference. I'm trying to run it from windows command line. I have Java 8 and successfully tested: C:> java -cp C:\temp\bbmap\current jgi.AssemblyStats2 in=C:\temp\bbmap\resources\phix174_ill.ref.fa.gz worked fine.

When I run C:\temp\bbmap>bbmap.sh ref=hg19.fasta. it tries to open a window and doesn't run it. I must be running it wrong. but what is the correct way? Please help it's driving me crazy and I really want to try align RNA-seq fastq. Thanks

[Brian Bushnell]

The shell script "bbmap.sh" won't work in Windows; it's for Linux (or MacOS). The command should be:

java -cp C:\temp\bbmap\current align2.BBMap -Xmx26g in=reads.fq ref=hg19.fasta out=mapped.sam maxindel=400k

...where "reads.fq" and "hg19.fasta" will need the full path if the are not in the current directory. BBMap needs around 24GB RAM or so for HG19, by the way. If you have paired reads in 2 files the syntax is "in1=read1.fq in2=read2.fq".