Hi everyone,
I thought that perhaps others have encountered and previously solved the following problem:
Samtools depth allows you to collect read depth information across a range in the genome using bed files to supply the coordinates.
My problem, is that in my output, the low quality regions or the elimination events appear as gaps rather than as zeros. Is there a way to change the parameters of samtools depth to report for all positions even if it just reports a 0? Or do I just need to write a script to fill in the blanks?
Current output:
Desired output:
**with the .. region filled in as well.
I have managed to write a code that identifies the breakpoints, but I'm not sure how to add additional lines that will run consecutively through the elimination site until it reaches the next set of depth values. Below is the code I've been working on (it's a perl one-liner and pretty raw):
All help is greatly appreciated!
I thought that perhaps others have encountered and previously solved the following problem:
Samtools depth allows you to collect read depth information across a range in the genome using bed files to supply the coordinates.
My problem, is that in my output, the low quality regions or the elimination events appear as gaps rather than as zeros. Is there a way to change the parameters of samtools depth to report for all positions even if it just reports a 0? Or do I just need to write a script to fill in the blanks?
Current output:
Code:
19 80820312 1 0 0 19 80820313 1 0 0 19 80820314 1 0 0 19 80826337 2 0 2 19 80826338 2 0 3 19 80826339 3 0 3 19 80826340 3 0 4
Code:
19 80826312 1 0 0 19 80820313 1 0 0 19 80820314 1 0 0 19 80820315 0 0 0 19 80820316 0 0 0 .. 19 80820335 0 0 0 19 80820336 0 0 0 19 80826337 2 0 2 19 80826338 2 0 3 19 80826339 3 0 3 19 80826340 3 0 4
I have managed to write a code that identifies the breakpoints, but I'm not sure how to add additional lines that will run consecutively through the elimination site until it reaches the next set of depth values. Below is the code I've been working on (it's a perl one-liner and pretty raw):
Code:
perl -e '$last == 0; while (<>){@a = split(/\t/, $_, 5); $n = $a[1]; if ($last == 0) { $last = ($a[1] + 1) ; } elsif ($n == $last){ $last = ($n + 1); next; } elsif ($n > ($last + 1)){ print "$a[0]\t$last\t0\t0\t0\n"; $last = ($n + 1); print;} }' depth1-2-3.bed
Comment