thanks - much tidier than my way. Janet
On Jun 23, 2011, at 2:08 AM, Martin Morgan wrote: > On 06/22/2011 05:33 PM, Janet Young wrote: >> Hi again, >> >> Thanks, Michael - looking forward to seqselect on GRanges being implemented, >> if possible. Getting a little further through my old code, brings up another >> related request. >> >> Here are some slightly different example ranges: >> myGR<- >> GRanges(seqnames=c("chr1","chr2","chr2"),ranges=IRanges(start=c(25,40,60),end=c(35,45,65)) >> ) >> myRD<- >> RangedData(space=c("chr1","chr2","chr2"),IRanges(start=c(25,40,60),end=c(35,45,65)) >> ) >> >> And I had previously been using this command: >> restrict(ranges(myRD), 30,50) >> I'm doing that because I'm interested in ranges within the same region on >> every single "chromosome". With my real data, they spaces/seqnames are not >> chromosomes, they're promoter regions centered around the transcription >> start site, and I want to take various equivalent sub-portions of the whole >> set of promoters) >> >> I can't get "restrict" to work on the GRanges object at the moment (even if >> I try that same coercion first): >> restrict( myGR , 30,50) >> Error in function (classes, fdef, mtable) : >> unable to find an inherited method for function "restrict", for >> signature "GRanges" >> >> restrict( as(myGR,"RangesList") , 30,50) >> Error in sapply(listData, function(Xi) extends(class(Xi), >> elementTypeX)) : >> error in evaluating the argument 'X' in selecting a method for >> function 'sapply': Error in .CompressedList.list.subscript(X = X, INDEX = >> seq_len(length(X)), : >> invalid output element of class "IRanges" >> >> It would be great if restrict could work directly in the same way as it did >> for RangedData. In the meantime I can probably figure out a way to do it >> with lapply. > > > or > > rng <- restrict(ranges(gr), 30, 50, keep.all.ranges=TRUE) > ranges(gr) <- rng > gr[width(gr) != 0] > > (and similarly for GRangesList). > > Martin > >> thanks, >> >> Janet >> >> >> >> On Jun 22, 2011, at 4:19 PM, Michael Lawrence wrote: >> >>> Hi, >>> >>> There is already a coercion from GRanges to RangesList, i.e., as(gr, >>> "RangesList"). That said, seqselect() should have a method for GRanges, if >>> possible. >>> >>> Michael >>> >>> On Wed, Jun 22, 2011 at 4:15 PM, Janet Young<jayo...@fhcrc.org> wrote: >>> Hi again, >>> >>> A further note - the workaround I first emailed only works when there's one >>> range per chromosome. I have a better workaround now - the coercion looks >>> like this instead: >>> >>> convertGRangesToCompressedIRangesList<- function (myGRobject) { >>> myGR_convert<- split(myGRobject, seqnames(myGRobject) ) >>> names(myGR_convert)<- unlist(lapply(myGR_convert, function(x) { >>> as.character(seqnames(x))[1] })) >>> >>> myGR_convert<- lapply(myGR_convert, ranges ) >>> class(myGR_convert)<- "CompressedIRangesList" >>> myGR_convert >>> } >>> >>> Janet >>> >>> >>> >>> Begin forwarded message: >>> >>>> From: Janet Young<jayo...@fhcrc.org> >>>> Date: June 22, 2011 4:01:03 PM PDT >>>> To: Bioc-sig-sequencing@r-project.org >>>> Subject: seqselect using GRanges object on multiple chromosomes >>>> >>>> Hi there, >>>> >>>> I'm updating some of my older code - I was previously storing regions of >>>> interest as RangedData but now I'm switching to GRanges. I'm running into >>>> a little trouble with seqselect - I've found a workaround but wanted to >>>> suggest extending seqselect so it can work with GRanges objects directly. >>>> >>>> I have some scores for each base-pair I've stored as a SimpleRleList >>>> object. I want to use GRanges object with seqselect to pull out scores >>>> from my regions of interest, but to make that work I first have to do an >>>> odd (and slightly wrong-looking, to me) coercion of my GRanges object to a >>>> CompressedIRangesList. >>>> >>>> I think the code below explains all (?). >>>> >>>> thanks very much, >>>> >>>> Janet >>>> >>>> >>>> >>>> >>>> >>>>> library(GenomicRanges) >>>> Loading required package: IRanges >>>> >>>> Attaching package: 'IRanges' >>>> >>>> The following object(s) are masked from 'package:base': >>>> >>>> cbind, eval, intersect, Map, mapply, order, paste, pmax, pmax.int, pmin, >>>> pmin.int, rbind, rep.int, setdiff, table, union >>>> >>>>> >>>>> ### make some scores objects, for single chromosomes, or across several >>>>> chrs >>>>> tempscores<- Rle(1:20) >>>>> tempscores2<- Rle(101:120) >>>>> allscores<- RleList(chr1=tempscores,chr2=tempscores2) # yields a >>>>> SimpleRleList >>>>> >>>>> ## make some ranges objects >>>>> myIR<- IRanges(start=3,end=5) >>>>> myGR<- GRanges(seqnames=c("chr1","chr2"),ranges=IRanges(start=3,end=5)) >>>>> myRD<- RangedData(space=c("chr1","chr2"),IRanges(start=c(3,3),end=c(5,5)) >>>>> ) >>>>> >>>>> >>>>> ### test seqselect: >>>>> seqselect(tempscores,myIR) #works >>>> 'integer' Rle of length 3 with 3 runs >>>> Lengths: 1 1 1 >>>> Values : 3 4 5 >>>>> seqselect(allscores,myGR) # doesn't work >>>> Error in seqselect(allscores, myGR) : unrecognized 'start' type >>>>> seqselect(allscores,myRD) # doesn't work >>>> Error in seqselect(allscores, myRD) : unrecognized 'start' type >>>>> >>>>> seqselect(allscores,ranges(myRD)) # works. ranges(myRD) is a >>>>> CompressedIRangesList >>>> SimpleRleList of length 2 >>>> $chr1 >>>> 'integer' Rle of length 3 with 3 runs >>>> Lengths: 1 1 1 >>>> Values : 3 4 5 >>>> >>>> $chr2 >>>> 'integer' Rle of length 3 with 3 runs >>>> Lengths: 1 1 1 >>>> Values : 103 104 105 >>>> >>>>> >>>>> seqselect(allscores,ranges(myGR)) # doesn't work >>>> Error in .bracket.Index(start, length(x), names(x), asRanges = TRUE) : >>>> range index out of bounds >>>>> seqselect(allscores,split(myGR)) # doesn't work >>>> Error in seqselect(allscores, split(myGR)) : unrecognized 'start' type >>>>> >>>>> #### coerce myGR to something that looks a bit like a >>>>> CompressedIRangesList >>>>> myGR_convert<- split(myGR) >>>>> names(myGR_convert)<- names(seqlengths(myGR_convert)) >>>>> myGR_convert<- lapply(myGR_convert, ranges ) >>>>> class(myGR_convert)<- "CompressedIRangesList" >>>>> >>>>> >>>>> seqselect(allscores,myGR_convert) # works >>>> SimpleRleList of length 2 >>>> $chr1 >>>> 'integer' Rle of length 3 with 3 runs >>>> Lengths: 1 1 1 >>>> Values : 3 4 5 >>>> >>>> $chr2 >>>> 'integer' Rle of length 3 with 3 runs >>>> Lengths: 1 1 1 >>>> Values : 103 104 105 >>>> >>>>> >>>>> sessionInfo() >>>> R version 2.13.0 (2011-04-13) >>>> Platform: i386-apple-darwin9.8.0/i386 (32-bit) >>>> >>>> locale: >>>> [1] en_US.UTF-8/en_US.UTF-8/C/C/en_US.UTF-8/en_US.UTF-8 >>>> >>>> attached base packages: >>>> [1] stats graphics grDevices utils datasets methods base >>>> >>>> other attached packages: >>>> [1] GenomicRanges_1.4.6 IRanges_1.10.4 >>>> >>>> >>> >>> _______________________________________________ >>> Bioc-sig-sequencing mailing list >>> Bioc-sig-sequencing@r-project.org >>> https://stat.ethz.ch/mailman/listinfo/bioc-sig-sequencing >>> >> >> _______________________________________________ >> Bioc-sig-sequencing mailing list >> Bioc-sig-sequencing@r-project.org >> https://stat.ethz.ch/mailman/listinfo/bioc-sig-sequencing > > > -- > Computational Biology > Fred Hutchinson Cancer Research Center > 1100 Fairview Ave. N. PO Box 19024 Seattle, WA 98109 > > Location: M1-B861 > Telephone: 206 667-2793 _______________________________________________ Bioc-sig-sequencing mailing list Bioc-sig-sequencing@r-project.org https://stat.ethz.ch/mailman/listinfo/bioc-sig-sequencing