[BioC] Turning a GRanges Metadata Column into Rle List

Hervé Pagès hpages at fhcrc.org
Wed Jan 9 02:10:47 CET 2013


Hi Dario,

Alternatively you could call coverage() a 2nd time with no weights to
find the regions with no coverage, and set them to NA:

   cvg <- coverage(gr, weight="scores")
   cvg[coverage(gr) == 0] <- NA

H.


On 01/07/2013 04:00 PM, Michael Lawrence wrote:
> I was thinking you could just do (assuming score is always >= 0, and that
> the ranges do not overlap, which seems to be the case from the initial
> code):
>
> gr$score <- gr$score + 1L
> cov <- coverage(gr, weight  = "score")
> cov[cov == 0L] <- NA
> cov <- cov - 1L
>
> Note that it should not be necessary to use endoapply or any other apply
> function, as it is all implicit in the List API.
>
> Not tested; just an illustration.
>
> Michael
>
>
> On Mon, Jan 7, 2013 at 2:55 PM, Aaron Statham <a.statham at garvan.org.au>wrote:
>
>> Here's a quick workaround (would be more elegant if weights in coverage()
>> could be NA, but that would probably break other functionality)
>>
>> g <- GRanges(c("chr1", "chr2"), IRanges(c(10, 50), c(15, 55)), scores
>> = c(20, 10))
>> seqlengths(g) <- c(100, 100)
>>
>> mask=-1 #placeholder for NA
>> #Define the ranges we want to be NA
>> gInv <- setdiff(GRanges(seqnames(g), IRanges(1, seqlengths(g))), g)
>> values(gInv)$scores <- mask
>> g <- c(g, gInv)
>>
>> endoapply(coverage(g, weight="scores"), function(x) {
>>      #replace mask with NA
>>      runValue(x)[runValue(x==mask)] <- NA
>>      x
>> })
>>
>> Aaron
>>
>>
>> On 8 January 2013 03:59, Michael Lawrence <lawrence.michael at gene.com>wrote:
>>
>>> One hack would be to add 1 to all of the scores, replace the zeros in the
>>> coverage() result with NAs and subtract 1. I've hit something like this a
>>> couple times in the past. The coverage method for GRanges could gain a
>>> default value argument.
>>>
>>> Michael
>>>
>>>
>>>
>>> On Sun, Jan 6, 2013 at 11:00 PM, Dario Strbenac <D.Strbenac at garvan.org.au
>>>> wrote:
>>>
>>>>>    coverage(g, weight=values(g)$scores)
>>>>
>>>> That is close to what I was after, although the ranges which aren't
>>>> present in the GRanges object effectively have missing scores and need to
>>>> be NA because 0 is a valid score for the ranges which are present. So, if I
>>>> had a view defined that only partially overlapped with a range in the
>>>> example, I would want any summary to give NA, rather than taking 0 into the
>>>> calculation.
>>>>
>>>> _______________________________________________
>>>> Bioconductor mailing list
>>>> Bioconductor at r-project.org
>>>> https://stat.ethz.ch/mailman/listinfo/bioconductor
>>>> Search the archives:
>>>> http://news.gmane.org/gmane.science.biology.informatics.conductor
>>>>
>>>
>>>
>>
>>
>> --
>> Aaron Statham
>> Postgraduate Scholar, Cancer Epigenetics
>> Garvan Institute of Medical Research   Tel: (02) 9295 8393
>> 384 Victoria St Darlinghurst 2010   Fax: (02) 9295 8316
>> NSW Australia         email: a.statham at garvan.org.au
>>
>
> 	[[alternative HTML version deleted]]
>
> _______________________________________________
> Bioconductor mailing list
> Bioconductor at r-project.org
> https://stat.ethz.ch/mailman/listinfo/bioconductor
> Search the archives: http://news.gmane.org/gmane.science.biology.informatics.conductor
>

-- 
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



More information about the Bioconductor mailing list