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:
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...
Thanks in advance
•
written
by
rotemkat •
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.
Was not aware of that - Thanks!! (and yes, i'm familiar with samtools, but it's so slow and i have many bam files...)
Log in to answer this question.
I don't understand what is your problem. All the data you want are here.
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.txtto get positions with depth more than 2. To get per base/per position coverage, you can tryigvtoolswith window size 1. You can also try https://github.com/sstadick/perbase tool in addition to igvtools. Perbase uses 1 based index by default.Thanks for your help!
the link is actually broken!
fixed it