max_depth of 8,000 too low for targeted sequencings.

[charles-plessy]

Dear developers,

in Debian we apply the following patch to raise the maximal depth to 1,000,000 when running samtools depth, because the original limit of 8,000 gets too easily reached in targeted sequencing applications.

Index: samtools/bam2depth.c
===================================================================
--- samtools.orig/bam2depth.c
+++ samtools/bam2depth.c
@@ -106,6 +106,7 @@

        // the core multi-pileup loop
        mplp = bam_mplp_init(n, read_bam, (void**)data); // initialization
+       bam_mplp_set_maxcnt(mplp,1000000); // set maxdepth to 1M
        n_plp = calloc(n, sizeof(int)); // n_plp[i] is the number of covering reads from the i-th BAM
        plp = calloc(n, sizeof(void*)); // plp[i] points to the array of covering reads (internal in mplp)
        while (bam_mplp_auto(mplp, &tid, &pos, n_plp, plp) > 0) { // come to the next covered position

Would you consider applying the patch or adding a run-time option ?

Have a nice day,

Charles Plessy
Debian Med packaging team
https://www.debian.org/devel/debian-med
Tsurumi, Kanagawa, Japan

[charles-plessy]

Oops sorry, here is a version that really applies.

--- a/bam2depth.c
+++ b/bam2depth.c
@@ -154,6 +154,7 @@

     // the core multi-pileup loop
     mplp = bam_mplp_init(n, read_bam, (void**)data); // initialization
+    bam_mplp_set_maxcnt(mplp,1000000); // set maxdepth to 1M
     n_plp = calloc(n, sizeof(int)); // n_plp[i] is the number of covering reads from the i-th BAM
     plp = calloc(n, sizeof(bam_pileup1_t*)); // plp[i] points to the array of covering reads (internal in mplp)
     while (bam_mplp_auto(mplp, &tid, &pos, n_plp, plp) > 0) { // come to the next covered position

[SamStudio8]

@charles-plessy It appears this patch has already been should have been added to samtools as per this message from @lh3 on the samtools-devel mailing list in 2012; http://sourceforge.net/p/samtools/mailman/message/29589210/

In future could you perhaps make pull requests so as to trigger the Travis tests and allow for diffing?

[SamStudio8]

@charles-plessy Sorry I managed to misread your patch file! The line is clearly still missing from bam2depth.c but was acknowledged for inclusion on the samtools-devel mailing list in 2012. Would you still be able to provide a pull request to trigger testing? I believe this would make it easier for the package maintainers to merge any such patch also.

[charles-plessy]

Thanks for the quick answer. See pull request #285. I did not know for Travis… nice feature !

[dcalderon]

Hello, I'm working with a file that also reaches the 8000 coverage limit too easily. Has this patch been added? Could you perhaps point me to the updated code that I could download and set up?

[SamStudio8]

@dcalderon The patch can be found at #285, it is still awaiting approval.
For the time being you can merge the patch to your own develop branch:

git remote add samtools_debian https://github.com/Debian/samtools.git
git fetch samtools_debian
git checkout develop
git merge samtools_debian/coverage-depth-cap
make
make install

[dcalderon]

Works like a charm. Thank you @SamStudio8 !

For anyone else reading this, remember to set the max depth to a larger value otherwise it defaults to 8000. So something like:

samtools mpileup -d 1000000 ...

[peterjc]

Can this issue be closed now since the -d option is present in official samtools releases?

[charles-plessy]

Indeed, -d obsoletes this pull request. I am closing it as you suggested.