Hi all,
I have a bam file and gff file and need to report per base coverage for all intervals (features) from gff file. Instead of having 0 coverage for all bps outside those regions, I simply need them not to be reported. I am familiar with bedtools and examined few tools, such as bedtools intersect, bedtools genomecov and bedtools coverage, but I don't seem to find the exact commands (yet). Number of reads per feature won't help me, I truly need per bp report. Thanks a lot:)
•
written
by
Biomonika (Noolean)
1 answer
You should be able to use bedtools coverage piped to awk. The formatting of the output from coverage with the -d option may or may not be that useful for you, but you can convert it to whatever you want with awk.
Log in to answer this question.
-
Count strand-specific 5' / 3' coverage across whole genome for paired-end fragments
written by acvill •
-
find coverage of specific exomes from a bam file
written by lait
-
Stats on the overlapping length within two bed files
written by beanavarro85 •
-
coverage on exons with bedtools genomecov
written by guillaume.rbt
-
Pick highest quality fastq reads
written by Biomonika (Noolean)
-
KEGG pathway per genes relation
written by bharata1803
-
Intersect multiple bed and keep all fields
written by user230613
-
span of GFF file
written by Biomonika (Noolean)
-
Tophat And Multiple-Mapped Reads
written by Biomonika (Noolean)
-
Coverage In Bam File - Bases And Overall Count
written by Biomonika (Noolean)
Have you considered just using samtools? Namely:
samtools depth -b foo.bed sample.bam? That won't double count overlapping paired end reads and should only print regions specified in the BED file. Remember to convert from 1-based to 0-based coordinates!