Variant calling on custom amplicons with extremely high coverage

memento · 02-11-2013, 12:36 PM

Hi,
I have a hard time obtaining variants from a bam file constructed using amplicon sequencing. When I load bam's on IGV I can see loads of reads (well, there is a limit of how much You can really see but the depth show correct coverage). I have 30.000-1.500.000 reads (identical ones) per each bam file.
I tried to use GATK with specific region bed file (thought that it will speed it up) but it takes forever (already 5h+! for 6 regions of ~250bp).
Tried samtools mpileup (multiple bam's at once) but I cannot find any of known variants (visible under IGV).
Shall I use it one-by-one on each bam (and combine?)?
Thanks!

swbarnes2 · 02-11-2013, 04:30 PM

Sounds like downsampling would make your life easier.

Removing duplicates would be one quick way of accomplishing this. But there are other ways. Picard, for one, can downsample a .bam. You can also do a quick and dirty version yourself, by grepping for reads, say, from one tile.

me_myself_andI · 02-18-2013, 09:38 PM

I would second the downsampling approach.

The removal of duplicates can in theory and depending on your setup introduce some biases. For example if your looking at subpopulations in viral or bacterial sequencing (i.e. not a 'simple' diploid genome) you might end up with only a handful of reads after duplicate removal, and those will not represent the actual 'allele' frequencies.

Similar Threads
Thread	Thread Starter	Forum	Replies	Last Post
Indel detection in high coverage amplicons sequenced by MiSeq	lpalacios	Bioinformatics	2	02-24-2015 07:04 AM
Extremely high Coverage in some regions	JMFA	Bioinformatics	2	01-11-2013 03:42 AM
False negative variant calling in haploids. Call variants using coverage (not stats)?	Genomics101	Bioinformatics	3	06-07-2012 12:29 PM
Variant calling for high-coverage Illumina data	dgmacarthur	Bioinformatics	19	04-08-2011 07:34 AM
New Paper: High Quality SNP Calling Using Illumina Data at Shallow Coverage	nmalhis	Bioinformatics	0	03-01-2010 03:40 PM

02-11-2013, 12:36 PM	#1
memento Member Location: UK Join Date: Nov 2010 Posts: 49	Variant calling on custom amplicons with extremely high coverage Hi, I have a hard time obtaining variants from a bam file constructed using amplicon sequencing. When I load bam's on IGV I can see loads of reads (well, there is a limit of how much You can really see but the depth show correct coverage). I have 30.000-1.500.000 reads (identical ones) per each bam file. I tried to use GATK with specific region bed file (thought that it will speed it up) but it takes forever (already 5h+! for 6 regions of ~250bp). Tried samtools mpileup (multiple bam's at once) but I cannot find any of known variants (visible under IGV). Shall I use it one-by-one on each bam (and combine?)? Thanks!

02-11-2013, 04:30 PM	#2
swbarnes2 Senior Member Location: San Diego Join Date: May 2008 Posts: 912	Sounds like downsampling would make your life easier. Removing duplicates would be one quick way of accomplishing this. But there are other ways. Picard, for one, can downsample a .bam. You can also do a quick and dirty version yourself, by grepping for reads, say, from one tile.

02-18-2013, 09:38 PM	#3
me_myself_andI Member Location: Singapore Join Date: Nov 2010 Posts: 30	I would second the downsampling approach. The removal of duplicates can in theory and depending on your setup introduce some biases. For example if your looking at subpopulations in viral or bacterial sequencing (i.e. not a 'simple' diploid genome) you might end up with only a handful of reads after duplicate removal, and those will not represent the actual 'allele' frequencies.