[BioC] IRanges : Strange behavior subsetting an IntervalTree with an indexing variable from within a function

[Steve Lianoglou]

Argh!

Wrapping into the for loop actually didn't fix anything -- the  
symptoms of the problem have matured in a way that make even less  
sense to me now.

It just seems as if the value of my indexing variable (``hits``) is  
stuck on the value that was first used when trying to index into my  
IntervalTree -- every time I use ``hits`` it always returns the same  
bunch of Intervals that were returned the first time my interval tree  
was indexed, even though its (hits) current value is different.

For instance, in the 6th iteration of the loop:

Browse[1]> hits
  [1] 3025 3026 3027 3028 3094 3095 3096 3097 3098 3099 3100 3101 3102  
3103

Browse[1]> ir[hits]
IntervalTree instance:
     start   end width
[1] 34576 34581     6
[2] 34608 34613     6
[3] 34635 34640     6
[4] 34888 34893     6

Clearly that's wrong there should be 13 intervals returned.

Those values being returned are from the "hits" vector that was  
calculated on the first pass of the loop, which was 158:161. This was  
calc'd 5 iterations ago, though, and should be long since gone!.

Browse[1]> ir[158:161]
IntervalTree instance:
     start   end width
[1] 34576 34581     6
[2] 34608 34613     6
[3] 34635 34640     6
[4] 34888 34893     6

Let's make ``hits`` something completely inane:

Browse[1]> hits <- 0
Browse[1]> ir[hits]
IntervalTree instance:
     start   end width
[1] 34576 34581     6
[2] 34608 34613     6
[3] 34635 34640     6
[4] 34888 34893     6

Meh ... this is only happening when indexing my IntervalTree, though:

Browse[1]> m <- matrix(1:100, 10,10)
Browse[1]> m[hits]
integer(0)


Still -- indexing with a NEW var, one that was never defined before,  
fails again on the IntervalTree:

Browse[1]> yy <- 1:10
Browse[1]> ir[yy]
Error in eval(expr, envir, enclos) : object 'yy' not found

Browse[1]> m[yy]
[1]  1  2  3  4  5  6  7  8  9 10

I'm stuck in bizarro land, and don't know how to get out ... anybody  
have a map?

-steve

On Sep 9, 2009, at 1:03 PM, Steve Lianoglou wrote:

> Hmm ... continuing from the previous email, switching the code out  
> of lapply and into a for loop seems to fix.
>
> So, instead of:
>
> dist <- lapply(range.list, function(rl) {
> hits <- subjectHits(overlap(ir, rl))
> browser() # The debug calls below start from here
> if (length(hits) > 1) {
>   diff(start(ir[hits]))
> } else {
>   NA
> }
> })
>
> do this:
>
> dist <- list()
> for (idx in seq(range.list)) {
>  hits <- subjectHits(overlap(ir, range.list[[idx]]))
>  dist[[idx]] <- if (length(hits) > 1) {
>    diff(start(ir[hits]))
>  } else {
>    NA
>  }
> }
>
> That works. Weird, no?
>
> Of course, I forgot to mention sessionInfo() -- sorry.
>
> R> > sessionInfo()
> R version 2.10.0 Under development (unstable) (2009-09-07 r49613)
> x86_64-apple-darwin9.8.0
>
> 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] RSQLite_0.7-2      DBI_0.2-4          ShortRead_1.3.33    
> lattice_0.17-25    BSgenome_1.13.11   doMC_1.1.1          
> multicore_0.1-3
> [8] foreach_1.2.1      codetools_0.2-2    iterators_1.0.2     
> Biostrings_2.13.36 IRanges_1.3.69     ARE.utils_0.1.0
>
> loaded via a namespace (and not attached):
> [1] Biobase_2.5.5 grid_2.10.0   hwriter_1.1   tools_2.10.0
>
> -steve
>
> --
> Steve Lianoglou
> Graduate Student: Computational Systems Biology
>  |  Memorial Sloan-Kettering Cancer Center
>  |  Weill Medical College of Cornell University
> Contact Info: http://cbio.mskcc.org/~lianos/contact
>