This is a test version of Biostars. For the public version, visit https://www.biostars.org.
Coveragebed, Depth/Breadth Of Coverage

I'm using coverageBed to calculate the depth and breadth of coverage, but I'm not sure I'm doing this right. I want to calculate the two values for each human chromosome.

For example, I've created a bed file with 1 chromosome. When I input my BAM file and the BED file, I get the following output:

chr1    0       249250621       103718897       224950839       249250621       0.9025086

I know the first 3 fields are from my chr BED file, the 4th field is the # of reads, 5th is # of bases covered, 6th is length of chromosome (redundant to field 3), and the last column is the fraction of bases covered (5th field/6th field).

So the 7th/last field gives the breadth of coverage, but I don't see a depth of coverage value. How do I get a depth of coverage?

bedtools coverage
• 12,541 views • link
updated by Istvan Albert Istvan Albert • written by Newbionf • Newbionf

I see you are looking for coverage on the entire human chrom1. The genomeCoverageBed tool is better suited to "genome scale" coverage statistics. coverageBed is designed for arbitrary BED intervals.

• 0 views • link
written by Aaronquinlan Aaronquinlan

Hello Newbionf. I just read this thread and am interested in how do you create a BED file per chromosome. Thanks. GP.

• 1 views • link
written by GPR GPR

1 answer

You will most likely need to pass another flag to customize the output. See the BEDTools manual::

One should use the “-hist” option to create, for each interval in B, a histogram of coverage of the features in A across B.

• 1 views • link
written by Istvan Albert Istvan Albert

I've done this for a 10kb region of chr1 on human and get the following output:

$ more chr1.bed chr1 0 249250621

[golharam@vail 8b-coverage]$ more S.chr1.10kb.coverage chr1 0 10000 0 9989 10000 0.9989000 chr1 0 10000 1 8 10000 0.0008000 chr1 0 10000 2 1 10000 0.0001000 chr1 0 10000 5 1 10000 0.0001000 chr1 0 10000 41 1 10000 0.0001000 all 0 9989 10000 0.9989000 all 1 8 10000 0.0008000 all 2 1 10000 0.0001000 ...

• 0 views • link
written by Newbionf • Newbionf

I've done this for a 10kb region on chr1 and get the followng output:

$ more S.chr1.10kb.coverage chr1 0 10000 0 9989 10000 0.9989000 chr1 0 10000 1 8 10000 0.0008000 chr1 0 10000 2 1 10000 0.0001000 chr1 0 10000 5 1 10000 0.0001000 chr1 0 10000 41 1 10000 0.0001000 all 0 9989 10000 0.9989000 all 1 8 10000 0.0008000 all 2 1 10000 0.0001000 all 5 1 10000 0.0001000 all 41 1 10000 0.0001000

• 0 views • link
written by Newbionf • Newbionf

Log in to answer this question.

More posts like this