SEQanswers

Go Back   SEQanswers > Bioinformatics > Bioinformatics



Similar Threads
Thread Thread Starter Forum Replies Last Post
how to calculate base-by-base coverage? shuoguo Bioinformatics 4 02-21-2014 08:06 PM
how to calculate the mean sequencing depth of RNA-seq wangli RNA Sequencing 3 02-19-2014 05:44 PM
calculate the base distribution of a SAM/BAM file weberma Bioinformatics 3 12-05-2012 12:22 AM
calculate a base's distance to the read's 3' end lyz1030 Bioinformatics 2 07-12-2012 06:36 AM
How to calculate proportion of reads with each base at every reference position lindseyjane Bioinformatics 7 01-19-2010 08:43 AM

Reply
 
Thread Tools
Old 10-16-2014, 01:00 AM   #1
Lv Ray
Member
 
Location: GZ,China

Join Date: Jun 2014
Posts: 42
Default Calculate depth of every base position

I have the sorted bam file, and I want to calculate the depth of every position on genome including some positions uncovered by reads(marked as: chr pos 0). I have tried samtools depth, but it seems don't count some position(because I count the row number, and not equal to the total number of genome bases).
Any advice about this would be much appreciated, thank you.
Lv Ray is offline   Reply With Quote
Old 10-16-2014, 01:21 PM   #2
blakeoft
Member
 
Location: Connecticut

Join Date: Oct 2013
Posts: 79
Default

I use bedtools genomecov. There are several output options, but I prefer to run
Code:
bedtools genomecov -bg -ibam your_bam.bam > your_bam.bedgraph
The fields in the output are chr, start, end and then depth. For example, you might see something like
Code:
chr1  554304  554309  5
chr1  554309  554313  6
This means that the depth is 5 from 554304 to 554309 and then the depth is 6 from 554309 to 554313 (I'm pretty sure it's not inclusive on the right endpoint) on chromosome 1.

There are some other options that I haven't actually explored myself, including adding "-d" to the code above that seems to output the depth at every single position, but this would make the file much larger than it would need to be and still contain exactly the same information that you'd have without "-d".

Read more about genomecov here.

Edit: samtools depth might do the same thing as I described above. I don't know since I haven't actually used it. You might want examine its output and see if it's presented in (what I like to call) the RLE way, like genomecov does with just the "-bg" parameter.

Last edited by blakeoft; 10-16-2014 at 01:25 PM.
blakeoft is offline   Reply With Quote
Old 10-17-2014, 02:56 AM   #3
WhatsOEver
Senior Member
 
Location: Germany

Join Date: Apr 2012
Posts: 215
Default

Just a small addition:

You need to call bedtools genomcov with -bga to get regions with 0 coverage. The rest blakeoft explained is correct.
WhatsOEver is offline   Reply With Quote
Old 10-17-2014, 05:35 AM   #4
blakeoft
Member
 
Location: Connecticut

Join Date: Oct 2013
Posts: 79
Default

Oops. I misread what Lv Ray posted. I thought that "marked as: chr pos 0" meant that you wanted to start the indexing at zero rather than one. Yes, WhatsOEver is correct. Thanks for clarifying that by the way.

For the record, this is the line you'll want to execute
Code:
bedtools genomecov -bga -ibam your_bam.bam > your_bam.bedgraph
blakeoft is offline   Reply With Quote
Old 10-17-2014, 09:53 AM   #5
Brian Bushnell
Super Moderator
 
Location: Walnut Creek, CA

Join Date: Jan 2014
Posts: 2,707
Default

BBTools contains a tool called pileup.sh, which can also do this rapidly from an unsorted sam or bam file.

pileup.sh -Xmx16g in=mapped.sam basecov=coverage.txt

It can also generate the coverage while mapping, if desired.
Brian Bushnell is offline   Reply With Quote
Old 10-17-2014, 11:04 PM   #6
Lv Ray
Member
 
Location: GZ,China

Join Date: Jun 2014
Posts: 42
Default

Thank you all. I will try the bedtools, and i will let u know what my result. Thanks again.
Lv Ray is offline   Reply With Quote
Reply

Tags
depth of every base

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:17 AM.


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