depth per position using mosdepth
1
0
Entering edit mode
4 months ago
rotemkat ▴ 10

Hi everyone, I'm trying to use "mosdepth" in order to get the depth in each position in a bam file i have (of rna-seq data). I tried to run -

mosdepth  output_file bam_file

and one of the files i get is:

output_file.per_base.bed.gz

When I open this file i see that actually the depth i get is not always per position, but most of the times it is for some preceding positions:

1. 1 18500 18505 8
2. 1 18505 18514 7
3. 1 18514 18517 6
4. 1 18517 18519 5

and so on...

Because i need the positions with depth >= 2, i'm not sure how i can use this output.

btw, i tried to use --thresholds 2, but it didn't work...

mosdepth RNA-Seq depth • 220 views
0
Entering edit mode

Because i need the positions with depth >= 2, i'm not sure how i can use this output.

I don't understand what is your problem. All the data you want are here.

0
Entering edit mode

Please read about bed format. Start position (2nd column) is inclusive and end position (3rd column) is non-inclusive. you can use awk \$4>=2 input.txt to get positions with depth more than 2. To get per base/per position coverage, you can try igvtools with window size 1. You can also try https://github.com/sstadick/perbase tool in addition to igvtools. Perbase uses 1 based index by default.

0
Entering edit mode

0
Entering edit mode

0
Entering edit mode

fixed it

0
Entering edit mode
4 months ago
ATpoint 49k

That is a bedGraph, it summarizes positions with the same depth into one interval, meaning that all the nucleotides at 1 18500 18505 have a coverage of 8. You can simply filter for column 4 to be greater than your threshold, e.g. with awk. A "real" per-base output with one entry per nucleotide could come from samtools depth.

0
Entering edit mode

Was not aware of that - Thanks!! (and yes, i'm familiar with samtools, but it's so slow and i have many bam files...)