## ---- message=FALSE, warning=FALSE, echo=FALSE, fig.width=8.5, fig.height=6.5---- ## Define colors colors <- c("#e19995", "#adaf64", "#4fbe9b", "#6eb3d9", "#d098d7") ## Create artificial GRanges library(GenomicRanges) set.seed(5) pool <- GRanges(seqnames = "chr1", ranges = IRanges(start = sample(1:800, 120, replace = TRUE), width = sample(25:200, 120, replace = TRUE)), color = sample(1:5, 120, replace = TRUE)) focal <- GRanges(seqnames = "chr1", ranges = IRanges(start = sample(1:800, 16, replace = TRUE), width = sample(25:200, 16, replace = TRUE)), color = sample(1:5, 16, replace = TRUE)) ## Add width to metadata pool$length <- width(pool) focal$length <- width(focal) ## Match ranges library(nullranges) set.seed(123) x <- matchRanges(focal = focal, pool = pool, covar = ~color + length, method = 'n', replace = TRUE) ## Visualize sets library(plotgardener) library(grid) set.seed(123) pageCreate(width = 8.5, height = 6.5, showGuides = FALSE, xgrid = 0, ygrid = 0) ## Define common parameters p <- pgParams(chrom = "chr1", chromstart = 1, chromend = 1000) ## Pool set poolSet <- plotRanges(data = pool, params = p, x = 1, y = 1, width = 2.5, height = 2.5, fill = colors, colorby = colorby("color")) annoGenomeLabel(plot = poolSet, x = 1, y = 3.55) plotText(label = "Pool Set", x = 2.25, y = 0.9, just = c("center", "bottom"), fontcolor = "#33A02C", fontface = "bold", fontfamily = 'mono') ## Focal set focalSet <- plotRanges(data = focal, params = p, x = 5, y = 1, width = 2.5, height = 1, fill = colors, colorby = colorby("color")) annoGenomeLabel(plot = focalSet, x = 5, y = 2.05) plotText(label = "Focal Set", x = 6.25, y = 0.9, just = c("center", "bottom"), fontcolor = "#1F78B4", fontface = "bold", fontfamily = 'mono') ## Matched set matchedSet <- plotRanges(data = matched(x), params = p, x = 5, y = 2.5, width = 2.5, height = 1, fill = colors, colorby = colorby("color")) annoGenomeLabel(plot = matchedSet, x = 5, y = 3.55) plotText(label = "Matched Set", x = 6.25, y = 2.75, just = c("center", "bottom"), fontcolor = "#A6CEE3", fontface = "bold", fontfamily = 'mono') ## Arrow and matchRanges label plotSegments(x0 = 3.5, y0 = 3, x1 = 5, y1 = 3, arrow = arrow(type = "closed", length = unit(0.1, "inches")), fill = "black", lwd = 2) plotText(label = "matchRanges()", fontfamily = 'mono', x = 4.25, y = 2.9, just = c("center", "bottom")) ## Matching plots library(ggplot2) smallText <- theme(legend.title = element_text(size=8), legend.text=element_text(size=8), title = element_text(size=8), axis.title.x = element_text(size=8), axis.title.y = element_text(size=8)) plot1 <- plotPropensity(x, sets=c('f','m','p')) + smallText + theme(legend.key.size = unit(0.5, 'lines'), title = element_blank()) plot2 <- plotCovariate(x=x, covar=covariates(x)[1], sets=c('f','m','p')) + smallText + theme(legend.text = element_blank(), legend.position = 'none') plot3 <- plotCovariate(x=x, covar=covariates(x)[2], sets=c('f','m','p'))+ smallText + theme(legend.key.size = unit(0.5, 'lines')) ## Propensity scores plotText(label = "plotPropensity()", x = 1.10, y = 4.24, just = c("left", "bottom"), fontface = "bold", fontfamily = 'mono') plotText(label = "~color + length", x = 1.25, y = 4.5, just = c("left", "bottom"), fontsize = 10, fontfamily = "mono") plotGG(plot = plot1, x = 1, y = 4.5, width = 2.5, height = 1.5, just = c("left", "top")) ## Covariate balance plotText(label = "plotCovariate()", x = 3.75, y = 4.24, just = c("left", "bottom"), fontface = "bold", fontfamily = "mono") plotText(label = covariates(x), x = c(4, 5.9), y = 4.5, just = c("left", "bottom"), fontsize = 10, fontfamily = "mono") plotGG(plot = plot2, x = 3.50, y = 4.5, width = 1.8, height = 1.5, just = c("left", "top")) plotGG(plot = plot3, x = 5.30, y = 4.5, width = 2.75, height = 1.5, just = c("left", "top")) ## ---- message=FALSE, warning=FALSE-------------------------------------------- library(nullrangesData) ## Load example data bins <- hg19_10kb_bins() bins ## ----------------------------------------------------------------------------- library(nullranges) ## Match ranges set.seed(123) mgr <- matchRanges(focal = bins[bins$looped], pool = bins[!bins$looped], covar = ~dnaseSignal + n_dnase_sites) mgr ## ---- message=FALSE, warning=FALSE-------------------------------------------- library(GenomicRanges) library(plyranges) library(ggplot2) ## Summarize ctcfSignal by n_ctcf_sites mgr %>% group_by(n_ctcf_sites) %>% summarize(ctcfSignal = mean(ctcfSignal)) %>% as.data.frame() %>% ggplot(aes(x = n_ctcf_sites, y = ctcfSignal)) + geom_line() + geom_point(shape = 21, stroke = 1, fill = 'white') + theme_minimal() + theme(panel.border = element_rect(color = 'black', fill = NA)) ## ----------------------------------------------------------------------------- overview(mgr) ## ---- message=FALSE----------------------------------------------------------- plotPropensity(mgr, sets = c('f', 'p', 'm'), type = 'ridges') ## ---- message=FALSE, warning=FALSE, fig.height=6, fig.width=5----------------- library(patchwork) plots <- lapply(covariates(mgr), plotCovariate, x=mgr, sets = c('f', 'm', 'p')) Reduce('/', plots) ## ----fig.width=4.5, fig.height=5---------------------------------------------- ## Percent of bins with CTCF g1 <- (sum(focal(mgr)$n_ctcf_sites >= 1) / length(focal(mgr))) * 100 g2 <- (sum(pool(mgr)$n_ctcf_sites >= 1) / length(pool(mgr))) * 100 g3 <- (sum(mgr$n_ctcf_sites >= 1) / length(mgr)) * 100 ## Visualize barplot(height = c(g1, g2, g3), names = c('looped\n(focal)', 'unlooped\n(pool)', 'unlooped\n(matched)'), ylab = "CTCF occupied bins (%)", col = c('#1F78B4', '#33A02C', '#A6CEE3'), main = 'CTCF occupancy', border = NA, las = 1) ## ----------------------------------------------------------------------------- sessionInfo()