[Bioc-devel] Why be default 'D' is not dropped from coverage?

[Leonardo Collado Torres]

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 <hpages at 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: hpages at fhcrc.org
> Phone:  (206) 667-5791
> Fax:    (206) 667-1319