How to collapse samtools depth generated bed file to regions with same coverage
2
0
Entering edit mode
6.9 years ago
BalderNGS • 0

I am creating per base depth files samtools depth. I want to collapse to a proper BED file for regions with similar depth:

Chr_01  140     10
Chr_01  141     10
Chr_01  142     10
Chr_01  143     30
Chr_01  144     30
Chr_01  145     30
Chr_01  146     30
Chr_01  147     20
Chr_01  148     20
Chr_01  149     20
Chr_01  150     20
Chr_01  151     20
Chr_01  152     10
Chr_01  153     10
Chr_01  154     10

Which should become (with trivial chosen names):

Chr_01  140     142    normal 10
Chr_01  143     146    high 30
Chr_01  147     151    medium  20
Chr_01  152     154    normal 10

Is this possible with bedtools or other existing tools?

samtools bedtools • 3.5k views
ADD COMMENT
1
Entering edit mode
6.9 years ago
$ cut -f 3 input.txt | sort| uniq | while read X; do awk -v X=$X '($3==X) { printf("%s\t%d\t%d\n",$1,$2,int($2)+1);}' input.txt | sort -k1,1 -k2,2n | bedtools merge -i - | sed "s/\$/\t${X}/" ; done

Chr_01  140 143 10
Chr_01  152 155 10
Chr_01  147 152 20
Chr_01  143 147 30

as it's BED, there is a +1 in the 'end' position.

ADD COMMENT
0
Entering edit mode

Thanks Pierre, I knew a single commandline would be out there. Thanks for the +1 on the end of the 'feature'. i missed that

ADD REPLY
0
Entering edit mode
6.9 years ago
ATpoint 85k

Tao Liu (one of the authors of MACS) has a script for this in his Git. Just double check if his script assumes 1-based or 0-based coordinates. You might need to add +1 to the coordinates of the output. Anyway, usage is ./script.pl samtools_depth.txt > your.output

ADD COMMENT
0
Entering edit mode

Thanks, I will checkout which of the solutions is faster on large files.

ADD REPLY

Login before adding your answer.

Traffic: 1953 users visited in the last hour
Help About
FAQ
Access RSS
API
Stats

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.

Powered by the version 2.3.6