### R code from vignette source 'IRangesOverview.Rnw' ### Encoding: UTF-8 ################################################### ### code chunk number 1: options ################################################### options(width=72) ################################################### ### code chunk number 2: biocLite (eval = FALSE) ################################################### ## source("http://bioconductor.org/biocLite.R") ## biocLite("IRanges") ################################################### ### code chunk number 3: initialize ################################################### library(IRanges) ################################################### ### code chunk number 4: initialize ################################################### set.seed(0) lambda <- c(rep(0.001, 4500), seq(0.001, 10, length = 500), seq(10, 0.001, length = 500)) xVector <- Rle(rpois(1e7, lambda)) yVector <- Rle(rpois(1e7, lambda[c(251:length(lambda), 1:250)])) ################################################### ### code chunk number 5: basic-ops ################################################### length(xVector) xVector[1] zVector <- c(xVector, yVector) ################################################### ### code chunk number 6: seq-extraction ################################################### xSnippet <- xVector[IRanges(4751, 4760)] xSnippet head(xSnippet) tail(xSnippet) rev(xSnippet) rep(xSnippet, 2) subset(xSnippet, xSnippet >= 5L) ################################################### ### code chunk number 7: seq-combine ################################################### c(xSnippet, rev(xSnippet)) append(xSnippet, xSnippet, after = 3) ################################################### ### code chunk number 8: aggregate ################################################### xSnippet aggregate(xSnippet, start = 1:8, width = 3, FUN = median) ################################################### ### code chunk number 9: shiftApply-cor ################################################### cor(xVector, yVector) shifts <- seq(235, 265, by=3) corrs <- shiftApply(shifts, yVector, xVector, FUN = cor) ################################################### ### code chunk number 10: figshiftcorrs ################################################### plot(shifts, corrs) ################################################### ### code chunk number 11: Rle-construction ################################################### xRle <- Rle(xVector) yRle <- Rle(yVector) xRle yRle ################################################### ### code chunk number 12: Rle-vector-compare ################################################### as.vector(object.size(xRle) / object.size(xVector)) identical(as.vector(xRle), xVector) ################################################### ### code chunk number 13: Rle-accessors ################################################### head(runValue(xRle)) head(runLength(xRle)) ################################################### ### code chunk number 14: Rle-ops ################################################### xRle > 0 xRle + yRle xRle > 0 | yRle > 0 ################################################### ### code chunk number 15: Rle-summary ################################################### range(xRle) sum(xRle > 0 | yRle > 0) ################################################### ### code chunk number 16: Rle-math ################################################### log1p(xRle) ################################################### ### code chunk number 17: Rle-cor ################################################### cor(xRle, yRle) shiftApply(249:251, yRle, xRle, FUN = function(x, y) var(x, y) / (sd(x) * sd(y))) ################################################### ### code chunk number 18: list-intro ################################################### getClassDef("RleList") ################################################### ### code chunk number 19: list-construct ################################################### args(IntegerList) cIntList1 <- IntegerList(x = xVector, y = yVector) cIntList1 sIntList2 <- IntegerList(x = xVector, y = yVector, compress = FALSE) sIntList2 ## sparse integer list xExploded <- lapply(xVector[1:5000], function(x) seq_len(x)) cIntList2 <- IntegerList(xExploded) sIntList2 <- IntegerList(xExploded, compress = FALSE) object.size(cIntList2) object.size(sIntList2) ################################################### ### code chunk number 20: list-length ################################################### length(cIntList2) Rle(elementNROWS(cIntList2)) ################################################### ### code chunk number 21: list-lapply ################################################### system.time(sapply(xExploded, mean)) system.time(sapply(sIntList2, mean)) system.time(sapply(cIntList2, mean)) identical(sapply(xExploded, mean), sapply(sIntList2, mean)) identical(sapply(xExploded, mean), sapply(cIntList2, mean)) ################################################### ### code chunk number 22: list-groupgenerics ################################################### xRleList <- RleList(xRle, 2L * rev(xRle)) yRleList <- RleList(yRle, 2L * rev(yRle)) xRleList > 0 xRleList + yRleList sum(xRleList > 0 | yRleList > 0) ################################################### ### code chunk number 23: list-endoapply ################################################### safe.max <- function(x) { if(length(x)) max(x) else integer(0) } endoapply(sIntList2, safe.max) endoapply(cIntList2, safe.max) endoapply(sIntList2, safe.max)[[1]] ################################################### ### code chunk number 24: iranges-constructor ################################################### ir1 <- IRanges(start = 1:10, width = 10:1) ir2 <- IRanges(start = 1:10, end = 11) ir3 <- IRanges(end = 11, width = 10:1) identical(ir1, ir2) & identical(ir2, ir3) ir <- IRanges(c(1, 8, 14, 15, 19, 34, 40), width = c(12, 6, 6, 15, 6, 2, 7)) ################################################### ### code chunk number 25: iranges-start ################################################### start(ir) ################################################### ### code chunk number 26: iranges-end ################################################### end(ir) ################################################### ### code chunk number 27: iranges-width ################################################### width(ir) ################################################### ### code chunk number 28: iranges-subset-numeric ################################################### ir[1:4] ################################################### ### code chunk number 29: iranges-subset-logical ################################################### ir[start(ir) <= 15] ################################################### ### code chunk number 30: ranges-extraction ################################################### ir[[1]] ################################################### ### code chunk number 31: plotRanges ################################################### plotRanges <- function(x, xlim = x, main = deparse(substitute(x)), col = "black", sep = 0.5, ...) { height <- 1 if (is(xlim, "Ranges")) xlim <- c(min(start(xlim)), max(end(xlim))) bins <- disjointBins(IRanges(start(x), end(x) + 1)) plot.new() plot.window(xlim, c(0, max(bins)*(height + sep))) ybottom <- bins * (sep + height) - height rect(start(x)-0.5, ybottom, end(x)+0.5, ybottom + height, col = col, ...) title(main) axis(1) } ################################################### ### code chunk number 32: ir-plotRanges ################################################### plotRanges(ir) ################################################### ### code chunk number 33: ranges-reduce ################################################### reduce(ir) plotRanges(reduce(ir)) ################################################### ### code chunk number 34: rangeslist-contructor ################################################### rl <- RangesList(ir, rev(ir)) ################################################### ### code chunk number 35: rangeslist-start ################################################### start(rl) ################################################### ### code chunk number 36: bracket-ranges ################################################### irextract <- IRanges(start = c(4501, 4901) , width = 100) xRle[irextract] ################################################### ### code chunk number 37: overlap-ranges ################################################### ol <- findOverlaps(ir, reduce(ir)) as.matrix(ol) ################################################### ### code chunk number 38: ranges-coverage ################################################### cov <- coverage(ir) plotRanges(ir) cov <- as.vector(cov) mat <- cbind(seq_along(cov)-0.5, cov) d <- diff(cov) != 0 mat <- rbind(cbind(mat[d,1]+1, mat[d,2]), mat) mat <- mat[order(mat[,1]),] lines(mat, col="red", lwd=4) axis(2) ################################################### ### code chunk number 39: ranges-shift ################################################### shift(ir, 10) ################################################### ### code chunk number 40: ranges-narrow ################################################### narrow(ir, start=1:5, width=2) ################################################### ### code chunk number 41: ranges-restrict ################################################### restrict(ir, start=2, end=3) ################################################### ### code chunk number 42: ranges-threebands ################################################### threebands(ir, start=1:5, width=2) ################################################### ### code chunk number 43: ranges-plus ################################################### ir + seq_len(length(ir)) ################################################### ### code chunk number 44: ranges-asterisk ################################################### ir * -2 # double the width ################################################### ### code chunk number 45: ranges-disjoin ################################################### disjoin(ir) plotRanges(disjoin(ir)) ################################################### ### code chunk number 46: ranges-disjointBins ################################################### disjointBins(ir) ################################################### ### code chunk number 47: ranges-reflect ################################################### reflect(ir, IRanges(start(ir), width=width(ir)*2)) ################################################### ### code chunk number 48: ranges-flank ################################################### flank(ir, width = seq_len(length(ir))) ################################################### ### code chunk number 49: ranges-gaps ################################################### gaps(ir, start=1, end=50) plotRanges(gaps(ir, start=1, end=50), c(1,50)) ################################################### ### code chunk number 50: ranges-pgap ################################################### ################################################### ### code chunk number 51: ranges-union ################################################### ################################################### ### code chunk number 52: ranges-punion ################################################### ################################################### ### code chunk number 53: ranges-intersect ################################################### ################################################### ### code chunk number 54: ranges-pintersect ################################################### ################################################### ### code chunk number 55: ranges-setdiff ################################################### ################################################### ### code chunk number 56: ranges-psetdiff ################################################### ################################################### ### code chunk number 57: Views-constructors ################################################### xViews <- Views(xRle, xRle >= 1) xViews <- slice(xRle, 1) xViewsList <- slice(xRleList, 1) ################################################### ### code chunk number 58: views-looping ################################################### head(viewSums(xViews)) viewSums(xViewsList) head(viewMaxs(xViews)) viewMaxs(xViewsList) ################################################### ### code chunk number 59: sessionInfo ################################################### toLatex(sessionInfo())