SEQanswers

Go Back   SEQanswers > Bioinformatics > Bioinformatics



Similar Threads
Thread Thread Starter Forum Replies Last Post
Local realignment around indels (454 technology) desmo Bioinformatics 2 05-17-2012 06:31 AM
GATK: VCF file for Local realignment around indels jorge Bioinformatics 2 10-11-2011 12:15 AM
Interpreting Local Realignment results from GATK alma Bioinformatics 0 07-07-2011 03:46 AM
GATK recalibration score and local realignment: what is the first step? m_elena_bioinfo Bioinformatics 0 01-24-2011 03:59 AM
Local realignment using GATK and smra seq_GA Bioinformatics 28 01-17-2011 09:10 AM

Reply
 
Thread Tools
Old 07-21-2011, 11:15 AM   #1
Hkins552
Member
 
Location: Baltimore MD

Join Date: Jun 2011
Posts: 18
Post Samtools mpileup creates extra large file after local realignment

Hi All,
Typically when running mpileup, my variants.raw.vcf files run around 4-6G (whole exome resequencing). I recently began using local realignment and recalibration of quality score with GATK and now the .vcf files from mpileup are 25-45G. What would be causing this?
Hkins552 is offline   Reply With Quote
Old 07-21-2011, 11:27 AM   #2
Richard Finney
Senior Member
 
Location: bethesda

Join Date: Feb 2009
Posts: 700
Default

Show the actual command you use.

Check your bcftools flags. Do you want all locations? Or just variants?

Last edited by Richard Finney; 07-21-2011 at 11:39 AM. Reason: typo fixed
Richard Finney is offline   Reply With Quote
Old 07-21-2011, 11:35 AM   #3
Hkins552
Member
 
Location: Baltimore MD

Join Date: Jun 2011
Posts: 18
Default

$ samtools mpileup -uf /hg19.fasta input.sorted.rmdup.reordered.realigned.recalibrated.bam > input.variants.raw

I would like to call just variants, not all loci.

I do not pipe directly to bcftools like the manual, but later use this command on the file from above:

$ bcftools view -bvcg input.variants.raw > input.variants.raw.bcf

$ bcftools view input.variants.raw.bcf | vcfutils.pl varFilter -d 3 -D 1000 -G 20 > input.variants.flt.vcf

Last edited by Hkins552; 07-21-2011 at 11:43 AM. Reason: added info
Hkins552 is offline   Reply With Quote
Old 07-21-2011, 12:17 PM   #4
oiiio
Senior Member
 
Location: USA

Join Date: Jan 2011
Posts: 105
Default

Have you already taken a look at the lines of your smaller and larger file? Are there differences?

Maybe worth mentioning is using the unix 'comm' command, it will compare the two files for you.
oiiio is offline   Reply With Quote
Old 07-21-2011, 12:20 PM   #5
swbarnes2
Senior Member
 
Location: San Diego

Join Date: May 2008
Posts: 912
Default

I do exome capture, and I use bedtools to filter my .bams against the capture probe .bed file.

That might help. It should winnow out some false aligning.
swbarnes2 is offline   Reply With Quote
Old 07-21-2011, 12:37 PM   #6
Michael.James.Clark
Senior Member
 
Location: Palo Alto

Join Date: Apr 2009
Posts: 213
Default

How many lines are there in input.variants.raw.bcf? In input.variants.flt.vcf?

Using the GATK pipeline, my exome-seq VCF files are on the order of 50-80,000 variants (lines) and a size of around 10-20Mb (depending on platform used for enrichment).

Quote:
Originally Posted by swbarnes2 View Post
I do exome capture, and I use bedtools to filter my .bams against the capture probe .bed file.

That might help. It should winnow out some false aligning.
If you do this, I advise you to use a modified capture probe .bed file where you've added 50 or 100 bases (or however much) to the end of each target region. You'll drop a lot of good data if you cut off right at the boundaries of the target intervals.
__________________
Mendelian Disorder: A blogshare of random useful information for general public consumption. [Blog]
Breakway: A Program to Identify Structural Variations in Genomic Data [Website] [Forum Post]
Projects: U87MG whole genome sequence [Website] [Paper]
Michael.James.Clark is offline   Reply With Quote
Old 07-21-2011, 01:03 PM   #7
swbarnes2
Senior Member
 
Location: San Diego

Join Date: May 2008
Posts: 912
Default

Quote:
Originally Posted by Michael.James.Clark View Post
If you do this, I advise you to use a modified capture probe .bed file where you've added 50 or 100 bases (or however much) to the end of each target region. You'll drop a lot of good data if you cut off right at the boundaries of the target intervals.
I'm pretty sure BEDTools includes reads that hang off the edge of your target, so you can still call SNPs that are just off target. But yes, I usually align to padded targets, to be sure, though I generally count coverage against unpadded target.

We use agilent capturing, and from the non-random sizes of the targets, I think the targets must be padded as well, with respect to the exons.
swbarnes2 is offline   Reply With Quote
Old 07-21-2011, 02:20 PM   #8
Michael.James.Clark
Senior Member
 
Location: Palo Alto

Join Date: Apr 2009
Posts: 213
Default

Quote:
Originally Posted by swbarnes2 View Post
I'm pretty sure BEDTools includes reads that hang off the edge of your target, so you can still call SNPs that are just off target. But yes, I usually align to padded targets, to be sure, though I generally count coverage against unpadded target.
Yeah, I think that's fair for coverage to be sure, but for variant calling I typically run GATK with an intervalsList containing the targets +50bp on each end.

Quote:
We use agilent capturing, and from the non-random sizes of the targets, I think the targets must be padded as well, with respect to the exons.
With respect to the exons that's typically true (with regards to RefSeq at least). The Agilent baits tend to extend outside the exons a bit. Still, you do end up pulling down significant excess in flanking regions as expected, so you gain even more. This is why I typically see a lot more variants than expected for exome alone in one of these experiments.
__________________
Mendelian Disorder: A blogshare of random useful information for general public consumption. [Blog]
Breakway: A Program to Identify Structural Variations in Genomic Data [Website] [Forum Post]
Projects: U87MG whole genome sequence [Website] [Paper]
Michael.James.Clark is offline   Reply With Quote
Reply

Tags
gatk, mpileup, realignment, recalibration

Thread Tools

Posting Rules
You may not post new threads
You may not post replies
You may not post attachments
You may not edit your posts

BB code is On
Smilies are On
[IMG] code is On
HTML code is Off




All times are GMT -8. The time now is 12:03 AM.


Powered by vBulletin® Version 3.8.9
Copyright ©2000 - 2020, vBulletin Solutions, Inc.
Single Sign On provided by vBSSO