################################################### ### chunk number 1: setup ################################################### options(width=50) set.seed(0) library(IRanges) 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) } ################################################### ### chunk number 2: sequences-example ################################################### s <- c(rep(0, 40), rep(c(1:10, seq(10, 1, -2), rep(1,6), seq(1,5,2), rep(5, 5), 4:1), rep(1:3, 11)), rep(0, 50)) ################################################### ### chunk number 3: sequences-example-plot ################################################### plot(s, type="l") ################################################### ### chunk number 4: sequences-diff ################################################### sum(diff(s) == 0) ################################################### ### chunk number 5: sequences-rle ################################################### sRle <- Rle(s) sRle ################################################### ### chunk number 6: rle-basic ################################################### sRle > 0 | rev(sRle) > 0 ################################################### ### chunk number 7: rle-summary ################################################### sum(sRle > 0) ################################################### ### chunk number 8: rle-cor ################################################### cor(sRle, rev(sRle)) ################################################### ### chunk number 9: rle-split ################################################### sRleList <- split(sRle, rep(c("chr1", "chr2"), each = length(sRle)/2)) sRleList[1] ################################################### ### chunk number 10: ranges-construct ################################################### ir <- IRanges(c(1, 8, 14, 15, 19, 34, 40), width=c(12, 6, 6, 15, 6, 2, 7)) ################################################### ### chunk number 11: ranges-plot ################################################### plotRanges(ir) ################################################### ### chunk number 12: ranges-start ################################################### start(ir) ################################################### ### chunk number 13: ranges-end ################################################### end(ir) ################################################### ### chunk number 14: ranges-width ################################################### width(ir) ################################################### ### chunk number 15: ranges-subset ################################################### ir[1:5] ################################################### ### chunk number 16: ranges-split ################################################### rl <- split(ir, c(rep("chr1", 2), rep("chr2", 3), "chr1", "chr2")) rl[1] ################################################### ### chunk number 17: ranges-shift ################################################### shift(ir, 1) ################################################### ### chunk number 18: ranges-shift-plot ################################################### plotRanges(shift(ir, 1)) ################################################### ### chunk number 19: ranges-resize ################################################### ir15 <- resize(ir, 15) # forward ir15 <- resize(ir, 15, start=FALSE) # backward ir15[1:3] ################################################### ### chunk number 20: ranges-restrict ################################################### restrict(ir15, 1) ################################################### ### chunk number 21: ranges-reduce ################################################### reduce(ir) ################################################### ### chunk number 22: ranges-reduce-plot ################################################### plotRanges(reduce(ir)) ################################################### ### chunk number 23: ranges-gaps ################################################### gaps(ir) ################################################### ### chunk number 24: ranges-gaps-plot ################################################### plotRanges(gaps(ir), ir) ################################################### ### chunk number 25: ranges-disjoin ################################################### disjoin(ir) ################################################### ### chunk number 26: ranges-disjoin-plot ################################################### plotRanges(disjoin(ir)) ################################################### ### chunk number 27: ranges-overlap ################################################### ol <- findOverlaps(ir, reduce(ir)) as.matrix(ol) ################################################### ### chunk number 28: ranges-coverage ################################################### cov <- coverage(ir) ################################################### ### chunk number 29: ranges-coverage-plot ################################################### 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) ################################################### ### chunk number 30: ranges-coverage ################################################### covL <- coverage(rl) covL ################################################### ### chunk number 31: views-sequence-plot ################################################### plot(s, type="l") abline(h = 3, col="red") ################################################### ### chunk number 32: views-rle ################################################### Views(sRle, as(sRle > 3, "IRanges")) ################################################### ### chunk number 33: views-slice ################################################### sViews <- slice(sRle, 4) ################################################### ### chunk number 34: views-slice-list ################################################### sViewsList <- slice(sRleList, 4) sViewsList[1] ################################################### ### chunk number 35: views-viewSums ################################################### viewSums(sViews) viewSums(sViewsList) ################################################### ### chunk number 36: views-viewMaxs ################################################### viewMaxs(sViews) viewMaxs(sViewsList) ################################################### ### chunk number 37: naive-RangedData-1 ################################################### chr <- c("chr1", "chr2", "chr1") strand <- c("+", "+", "-") start <- c(3L, 4L, 1L) end <- c(7L, 5L, 3L) score <- c(1L, 3L, 2L) naiveTable <- data.frame(chr, strand, score, start, end) naiveTable ################################################### ### chunk number 38: naive-RangedData-2 ################################################### getRange <- function(x) range(x[c("start", "end")]) by(naiveTable, naiveTable[["chr"]], getRange) ################################################### ### chunk number 39: RangedData-construction ################################################### rdTable <- RangedData(IRanges(start, end), strand, score, space = chr) ################################################### ### chunk number 40: RangedData-display ################################################### rdTable ################################################### ### chunk number 41: rd-ranges ################################################### ranges(rdTable)[1] ################################################### ### chunk number 42: rd-values ################################################### values(rdTable)[1] ################################################### ### chunk number 43: attribute-accessors ################################################### start(rdTable) ################################################### ### chunk number 44: access-column ################################################### rdTable$strand ################################################### ### chunk number 45: subset-features ################################################### ## get the first 2 features first10 <- rdTable[1:2,] ## get features on positive strand pos <- rdTable[rdTable$strand == "+",] ## get first chromosome chr1Table <- rdTable[1] ## with only the score column scoreTable <- rdTable[,"score"] ################################################### ### chunk number 46: rle-to-rd ################################################### head(as(sRleList, "RangedData"), 3) ################################################### ### chunk number 47: rle-views-rd ################################################### height <- unlist(viewMaxs(sViewsList)) RangedData(sViewsList, height) ################################################### ### chunk number 48: rd-to-rle ################################################### coverage(rdTable, weight = "score")[1]