targeted bed file average coverage
I am trying to modify this very useful post by Michael James Clark which he wrote a perl script to calculate average coverage in a bam.
coverage.pl
($num,$den)=(0,0);
while ($cov=<STDIN>) {
$num=$num+$cov;
$den++;
}
$cov=$num/$den;
print "Mean Coverage = $cov\n";What I am trying to do is use that code to calculate average coverage in a bam for targeted regions in a bed file.
/path/to/samtools view -b in.bam <genomic region> | /path/to/samtools mpileup - | awk '{print $4}' | perl ~/coverage.plThank you :).
0 answers
No answers yet.
Log in to answer this question.
More posts like this
-
Request to validate my RNAseq workflow
written by umn_bist •
-
picard error version 1.139 on ubuntu 14.04
written by bioguy24
-
Piping samtools commands
written by jyu429
-
subset fastq file
written by cara78 •
-
Trouble Passing Parameters To Bowtie
written by Nick
-
Problem Running Htseq
written by lkmklsmn
-
Blast Xml Document Iterations
written by pld
-
How To Calculate The Coverage Of A Specific Interval Based On Data From A Bam File
written by rehma.ar •
-
Problem With Mpileup
written by virtuouslycan •
-
Why Import Command In Samtools Is Not Shown In The Short Help Of Samtools?
written by Pablo Marin-Garcia
If you looking for alternative, both
bedtoolsand bedmap frombedOpscan do it standalone.Thank you :)
Here is a nice blog post describing cumulative coverage across bait regions using bedtools. Very elegant.
Thank you, I will try bedtools. Does samtools have a similar feature?