SEQanswers

Go Back   SEQanswers > Bioinformatics > Bioinformatics



Similar Threads
Thread Thread Starter Forum Replies Last Post
Extract genotype consensus from BAM files at specific locations gr8guns Bioinformatics 8 06-19-2015 09:04 AM
How to extract a regional consensus sequence of reads from BAM file jeffseq Bioinformatics 8 01-29-2015 01:05 PM
How to extract unmapped reference sequence from BAM/SAM file.? Pinal Bioinformatics 2 03-06-2014 01:40 PM
fill gaps in genome by transcripts neurospora Bioinformatics 2 08-28-2012 01:47 AM
Want to generate more data to fill in gaps Salmon 454 Pyrosequencing 2 03-19-2009 08:33 AM

Reply
 
Thread Tools
Old 03-26-2014, 10:15 AM   #1
Danko1
Junior Member
 
Location: Germany

Join Date: Mar 2014
Posts: 6
Default Extract consensus from BAM + fill gaps from reference

Hi everyone!

I am new to this forum and have the following question:

How can I extract the consensus sequence from a BAM alignment AND fill positions where there might be gaps in the alignment with the reference sequence from the alignment?

I know that there is a way using samtools mpileup to extract a consensus sequence from a BAM file but the "fill from reference" option is very important to me.

Any ideas what tools I could use to easily extract the consensus this way?

Thanks for your help!
Danko1 is offline   Reply With Quote
Old 03-26-2014, 10:20 AM   #2
dpryan
Devon Ryan
 
Location: Freiburg, Germany

Join Date: Jul 2011
Posts: 3,480
Default

Quoting from one of the samtools pages:

Code:
samtools mpileup -uf ref.fa aln.bam | bcftools view -cg - | vcfutils.pl vcf2fq > cns.fq
dpryan is offline   Reply With Quote
Old 03-26-2014, 10:25 AM   #3
Danko1
Junior Member
 
Location: Germany

Join Date: Mar 2014
Posts: 6
Default

Thanks for your quick answer dpryan!

I know this samtools code, but this won't allow to fill gaps in the alignment with the reference sequence of the alignment at that position so that I end up with a consensus fasta file without any gaps, right?

How can I extract the consensus sequence from an alignment where positions where there are gaps are filled with the corresponding nucleotide(s) from the reference sequence of the alignment?
Danko1 is offline   Reply With Quote
Old 03-26-2014, 10:32 AM   #4
dpryan
Devon Ryan
 
Location: Freiburg, Germany

Join Date: Jul 2011
Posts: 3,480
Default

That actually should produce a consensus file. At least the output from bcftools will contain the SNP calls for every base and the reference sequence (regardless of whether the base is covered or not), so if that's not working correctly then presumably there's been a change in vcfutils.pl. To confirm this, what happens if you:

Code:
samtools mpileup -uf ref.fa aln.bam | bcftools view -cg - >foo.vcf
and then look at one of the gaps? It should be covered in the VCF (if not, then I likely mistyped something). The only trick then is converting the reference given the VCF file, which which there are other options (GATK has a tool for that).
dpryan is offline   Reply With Quote
Old 03-27-2014, 11:59 AM   #5
Danko1
Junior Member
 
Location: Germany

Join Date: Mar 2014
Posts: 6
Default

So I tested the proposed code to generate a consensus from a BAM alignment:

samtools mpileup -uf ref.fa aln.bam | bcftools view -cg - | vcfutils.pl vcf2fq > cns.fq

It works fine, despite the fact that I end up with 'N' at every position where my reads don't map to the reference in the alignment.

How can I replace those 'N's' by the original nucleotides from the reference at this position? I am glad for any help!

@dpryan: Did you mean the FastaAlternateReferenceMaker by GATK? Have you used this tool before?
Danko1 is offline   Reply With Quote
Old 03-28-2014, 08:13 AM   #6
Danko1
Junior Member
 
Location: Germany

Join Date: Mar 2014
Posts: 6
Default

Figured it out! The GATK tool did the trick! Thanks for your help!
Danko1 is offline   Reply With Quote
Reply

Tags
alignment, bioinformatics, consensus extraction, fill from reference

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 03:24 PM.


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