Thank you Hervé! I failed to realized that this was included in the docs. Anyhow, we've been discussing internally what default value to use. We also noted that the `genomecov` tool from BED in it's split mode ignores the D's: http://bedtools.readthedocs.org/en/latest/content/tools/genomecov.html
Cheers, Leo On Mon, Oct 20, 2014 at 5:45 PM, Hervé Pagès <hpa...@fhcrc.org> wrote: > Hi Leo, > > > On 10/18/2014 10:50 AM, Leonardo Collado Torres wrote: >> >> Hi, >> >> A collaborator of mine is working on a new software and we while we were >> doing a sanity check to compare the base-level coverage from BAM files >> and bigWig files generated from his software we realized that by default >> bases corresponding to a 'D' on the CIGAR string get counted when >> reading coverage from BAM files. >> >> That is: >> >> > getMethod('coverage', 'GAlignments') >> Method Definition: >> >> function (x, shift = 0L, width = NULL, weight = 1L, ...) >> { >> .local <- function (x, shift = 0L, width = NULL, weight = 1L, >> method = c("auto", "sort", "hash"), drop.D.ranges = FALSE) >> coverage(grglist(x, drop.D.ranges = drop.D.ranges), shift = shift, >> width = width, weight = weight, method = method) >> .local(x, shift, width, weight, ...) >> } >> <environment: namespace:GenomicAlignments> >> >> Signatures: >> x >> target "GAlignments" >> defined "GAlignments" >> > packageVersion('GenomicAlignments') >> [1] ‘1.2.0’ >> >> >> You can see that this is default elsewhere, for example: >> >> > extractAlignmentRangesOnReference >> function (cigar, pos = 1L, drop.D.ranges = FALSE, f = NULL) >> { >> if (!isTRUEorFALSE(drop.D.ranges)) >> stop("'drop.D.ranges' must be TRUE or FALSE") >> if (drop.D.ranges) { >> ops <- c("M", "=", "X", "I") >> } >> else { >> ops <- c("M", "=", "X", "I", "D") >> } >> cigarRangesAlongReferenceSpace(cigar, flag = NULL, pos = pos, >> f = f, ops = ops, drop.empty.ranges = FALSE, reduce.ranges = >> TRUE) >> } >> <environment: namespace:GenomicAlignments> >> >> >> What is the rationale for setting `drop.D.ranges` by default to FALSE? > > > Because last time I checked (but that was 3 or 4 years ago), that's what > Samtools pileup was doing. > > Cheers, > H. > >> >> Thanks, >> Leo >> >> > > -- > Hervé Pagès > > Program in Computational Biology > Division of Public Health Sciences > Fred Hutchinson Cancer Research Center > 1100 Fairview Ave. N, M1-B514 > P.O. Box 19024 > Seattle, WA 98109-1024 > > E-mail: hpa...@fhcrc.org > Phone: (206) 667-5791 > Fax: (206) 667-1319 _______________________________________________ Bioc-devel@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/bioc-devel