Hi I'm trying ReadDepth to call CNV. Anyone knows how to prepare the input Bed file from the BAM file. There are only short descriptions in the example file. I'm not sure how it works. Is there any tool or script could do this? Anyway, ReadDepth is different from other tools in the input file format. :-)
2 answers
bamToBed would get you part of the way there, but there are two caveats: 1) You need per-chromosome bed files 2) readDepth actually only looks at the first and second columns. You can save yourself some disk by omitting the third column. Use the following command on your sorted bam file instead:
for $chr in $(seq 1 22) X Y;do samtools view myfile.bam $chr:1-999999999 | cut -f 3,4 >$chr.bed;done
If your bam uses the "chr" prefix ("chr1" instead of "1"), use this instead:
for $chr in $(seq 1 22) X Y;do samtools view myfile.bam chr$chr:1-999999999 | cut -f 3,4 >chr$chr.bed;done
Keep in mind that this will use a lot of disk space! Just make sure the files match your annotations. (either use "chr" or don't, but be consistent).
Sidenote: If I had to write this package over again, it would read bam/sam natively, but a lot of this work was done before bam/sam was in widespread use. I'm working on a new followup tool that should be out this summer.
Thanks very much for your reply. Looking forward to your new tools.
Log in to answer this question.
Is it different than normal bed file, if no try bamToBed, if yes how does the input format looks like!!
ReadDepth uses "mapped reads, in one-based bed format, broken into one file for each chromosome" - bamToBed does the job.