Filtering small indel false positives

[BetterPrimate]

I have some Solid exome sequence data including small indel calls. For a single exome (~38mb target) I get ~1000 indels called of which ~150 map to exons of which ~20 are present in dbSNP129.

I don't want to believe that there are ~130 genes blighted with indels in this exome. This seems unlikely because when I compare multiple exomes the same gene names crop up repeated.

How best to sort the wheat from the chaff?

BTW The ABI small indel tool was run with default parameters. i.e. Minimum 2 reads to call an indel. Max (normal reads coverage)/(indel reads coverage) < 12 times, because false +ve more likely at high coverage. etc.

[colindaven]

False positives are also frequent at low coverage. I wouldn't call indels based on two reads. I hope you are planning to validate in the wet lab!

Have you tried more exhaustive / alternate analysis ? I got some indel calls out of Dindel recently with Illumina data, we'll be attempting to validate these soon.

Also, are you really happy with your alignment ? The aligner used might also be the limiting factor here. Perhaps repeat the analysis with a different aligner and downstream indel analysis tool.

[BetterPrimate]

Thanks for the reply Colin. Apparently my question was badly worded. I'm not looking for alternatives to the Abi small indel tool. Initially I'm not assuming there is anything wrong with the data, but am not sure how to interpret it.

[Chipper]

I would look at the number of reads (and their qualities) for the 20 dbSNP indels compared to the others, as well as the position in reads. This should give you an idea if you need to adjust the cut-offs.

[Michael.James.Clark]

I think you're underestimating the number of indels you should be seeing by quite a bit.

I observe 4-6000 indels (depending on read depth) in human exomes using Nimblegen EZ Exome (44Mb target), with 5-700 in RefSeq mRNA codons alone.

That may seem high to you, but we do see a pretty drastic preference against coding indels compared to coding SNPs. About 1% of all SNPs fall in coding regions, but only 0.1% of indels fall in coding regions. Those that do fall in the coding regions have a size distribution that prefers multiples of three to maintain frame as well.

I think you'd be best off re-analyzing indels with Dindel or GATK UnifiedGenotyper. You'll get a higher yield, but you actually should have a higher number than you're quoting now.

Then you can use a functional analysis tool to determine the impact of those indels. For example, Annovar can take in indels and determine if they're causing a frameshift/nonsense/etc. (You can go ahead and do that with the list you have now of course. It's up to you how well you want to analyze this--if you're satisfied with those calls, have at it.)

[arcolombo698]

I think you're underestimating the number of indels you should be seeing by quite a bit.

I observe 4-6000 indels (depending on read depth) in human exomes using Nimblegen EZ Exome (44Mb target), with 5-700 in RefSeq mRNA codons alone.

That may seem high to you, but we do see a pretty drastic preference against coding indels compared to coding SNPs. About 1% of all SNPs fall in coding regions, but only 0.1% of indels fall in coding regions. Those that do fall in the coding regions have a size distribution that prefers multiples of three to maintain frame as well.

I think you'd be best off re-analyzing indels with Dindel or GATK UnifiedGenotyper. You'll get a higher yield, but you actually should have a higher number than you're quoting now.

Then you can use a functional analysis tool to determine the impact of those indels. For example, Annovar can take in indels and determine if they're causing a frameshift/nonsense/etc. (You can go ahead and do that with the list you have now of course. It's up to you how well you want to analyze this--if you're satisfied with those calls, have at it.)

what about using GATK for an RNAseq project? I have been reading a little about GATK and the unifiedgenotyper, and the haplotype caller of GATK.