In this vignette, we demonstrate the unsegmented block bootstrap functionality implemented in nullranges. “Unsegmented” refers to the fact that this implementation does not consider segmentation of the genome for sampling of blocks, see the segmented block bootstrap vignette for the alternative implementation.

Timing on DHS peaks

First we use the DNase hypersensitivity peaks in A549 downloaded from AnnotationHub, and pre-processed as described in the nullrangesData package.

library(nullrangesData)
dhs <- DHSA549Hg38()
library(nullranges)

The following chunk of code evaluates various types of bootstrap/permutation schemes, first within chromosome, and then across chromosome (the default). The default type is bootstrap, and the default for withinChrom is FALSE (bootstrapping with blocks moving across chromosomes).

set.seed(5) # reproducibility
library(microbenchmark)
blockLength <- 5e5
microbenchmark(
  list=alist(
    p_within=bootRanges(dhs, blockLength=blockLength,
                        type="permute", withinChrom=TRUE),
    b_within=bootRanges(dhs, blockLength=blockLength,
                        type="bootstrap", withinChrom=TRUE),
    p_across=bootRanges(dhs, blockLength=blockLength,
                        type="permute", withinChrom=FALSE),
    b_across=bootRanges(dhs, blockLength=blockLength,
                        type="bootstrap", withinChrom=FALSE)
  ), times=10)
## Unit: milliseconds
##      expr      min        lq      mean    median        uq      max neval cld
##  p_within 980.5564 1023.8027 1076.8785 1068.1294 1140.0632 1172.357    10   b
##  b_within 905.4653  928.8574 1066.3188  944.3148  958.8301 2081.844    10   b
##  p_across 230.3390  252.2741  273.5545  268.1419  306.2476  313.606    10  a 
##  b_across 282.5991  301.1342  438.2609  337.2092  342.1123 1445.015    10  a

Visualize on synthetic data

We create some synthetic ranges in order to visualize the different options of the unsegmented bootstrap implemented in nullranges.

library(GenomicRanges)
seq_nms <- rep(c("chr1","chr2","chr3"),c(4,5,2))
gr <- GRanges(seqnames=seq_nms,
              IRanges(start=c(1,101,121,201,
                              101,201,216,231,401,
                              1,101),
                      width=c(20, 5, 5, 30,
                              20, 5, 5, 5, 30,
                              80, 40)),
              seqlengths=c(chr1=300,chr2=450,chr3=200),
              chr=factor(seq_nms))

The following function uses functionality from plotgardener to plot the ranges. Note in the plotting helper function that chr will be used to color ranges by chromosome of origin.

suppressPackageStartupMessages(library(plotgardener))
plotGRanges <- function(gr) {
  pageCreate(width = 5, height = 2, xgrid = 0,
                ygrid = 0, showGuides = FALSE)
  for (i in seq_along(seqlevels(gr))) {
    chrom <- seqlevels(gr)[i]
    chromend <- seqlengths(gr)[[chrom]]
    suppressMessages({
      p <- pgParams(chromstart = 0, chromend = chromend,
                    x = 0.5, width = 4*chromend/500, height = 0.5,
                    at = seq(0, chromend, 50),
                    fill = colorby("chr", palette=palette.colors))
      prngs <- plotRanges(data = gr, params = p,
                          chrom = chrom,
                          y = 0.25 + (i-1)*.7,
                          just = c("left", "bottom"))
      annoGenomeLabel(plot = prngs, params = p, y = 0.30 + (i-1)*.7)
    })
  }
}
plotGRanges(gr)

Within chromosome

Visualizing two permutations of blocks within chromosome:

for (i in 1:2) {
  gr_prime <- bootRanges(gr, blockLength=100, type="permute", withinChrom=TRUE)
  plotGRanges(gr_prime)
}

Visualizing two bootstraps within chromosome:

for (i in 1:2) {
  gr_prime <- bootRanges(gr, blockLength=100, withinChrom=TRUE)
  plotGRanges(gr_prime)
}

Across chromosome

Visualizing two permutations of blocks across chromosome. Here we use larger blocks than previously.

for (i in 1:2) {
  gr_prime <- bootRanges(gr, blockLength=200, type="permute", withinChrom=FALSE)
  plotGRanges(gr_prime)
}

Visualizing two bootstraps across chromosome:

for (i in 1:2) {
  gr_prime <- bootRanges(gr, blockLength=200, withinChrom=FALSE)
  plotGRanges(gr_prime)
}

Session information

sessionInfo()
## R version 4.2.1 (2022-06-23 ucrt)
## Platform: x86_64-w64-mingw32/x64 (64-bit)
## Running under: Windows Server x64 (build 20348)
## 
## Matrix products: default
## 
## locale:
## [1] LC_COLLATE=C                          
## [2] LC_CTYPE=English_United States.utf8   
## [3] LC_MONETARY=English_United States.utf8
## [4] LC_NUMERIC=C                          
## [5] LC_TIME=English_United States.utf8    
## 
## attached base packages:
## [1] grid      stats4    stats     graphics  grDevices utils     datasets 
## [8] methods   base     
## 
## other attached packages:
##  [1] microbenchmark_1.4.9        purrr_0.3.5                
##  [3] ggridges_0.5.4              tidyr_1.2.1                
##  [5] EnsDb.Hsapiens.v86_2.99.0   ensembldb_2.22.0           
##  [7] AnnotationFilter_1.22.0     GenomicFeatures_1.50.0     
##  [9] AnnotationDbi_1.60.0        patchwork_1.1.2            
## [11] plyranges_1.18.0            nullrangesData_1.3.0       
## [13] ExperimentHub_2.6.0         AnnotationHub_3.6.0        
## [15] BiocFileCache_2.6.0         dbplyr_2.2.1               
## [17] ggplot2_3.3.6               plotgardener_1.4.0         
## [19] nullranges_1.4.0            InteractionSet_1.26.0      
## [21] SummarizedExperiment_1.28.0 Biobase_2.58.0             
## [23] MatrixGenerics_1.10.0       matrixStats_0.62.0         
## [25] GenomicRanges_1.50.0        GenomeInfoDb_1.34.0        
## [27] IRanges_2.32.0              S4Vectors_0.36.0           
## [29] BiocGenerics_0.44.0        
## 
## loaded via a namespace (and not attached):
##   [1] RcppHMM_1.2.2                 lazyeval_0.2.2               
##   [3] splines_4.2.1                 BiocParallel_1.32.0          
##   [5] TH.data_1.1-1                 digest_0.6.30                
##   [7] yulab.utils_0.0.5             htmltools_0.5.3              
##   [9] fansi_1.0.3                   magrittr_2.0.3               
##  [11] memoise_2.0.1                 ks_1.13.5                    
##  [13] Biostrings_2.66.0             sandwich_3.0-2               
##  [15] prettyunits_1.1.1             jpeg_0.1-9                   
##  [17] colorspace_2.0-3              blob_1.2.3                   
##  [19] rappdirs_0.3.3                xfun_0.34                    
##  [21] dplyr_1.0.10                  crayon_1.5.2                 
##  [23] RCurl_1.98-1.9                jsonlite_1.8.3               
##  [25] survival_3.4-0                zoo_1.8-11                   
##  [27] glue_1.6.2                    gtable_0.3.1                 
##  [29] zlibbioc_1.44.0               XVector_0.38.0               
##  [31] strawr_0.0.9                  DelayedArray_0.24.0          
##  [33] scales_1.2.1                  mvtnorm_1.1-3                
##  [35] DBI_1.1.3                     Rcpp_1.0.9                   
##  [37] xtable_1.8-4                  progress_1.2.2               
##  [39] gridGraphics_0.5-1            bit_4.0.4                    
##  [41] mclust_6.0.0                  httr_1.4.4                   
##  [43] RColorBrewer_1.1-3            speedglm_0.3-4               
##  [45] ellipsis_0.3.2                pkgconfig_2.0.3              
##  [47] XML_3.99-0.12                 farver_2.1.1                 
##  [49] sass_0.4.2                    utf8_1.2.2                   
##  [51] DNAcopy_1.72.0                ggplotify_0.1.0              
##  [53] tidyselect_1.2.0              labeling_0.4.2               
##  [55] rlang_1.0.6                   later_1.3.0                  
##  [57] munsell_0.5.0                 BiocVersion_3.16.0           
##  [59] tools_4.2.1                   cachem_1.0.6                 
##  [61] cli_3.4.1                     generics_0.1.3               
##  [63] RSQLite_2.2.18                evaluate_0.17                
##  [65] stringr_1.4.1                 fastmap_1.1.0                
##  [67] yaml_2.3.6                    knitr_1.40                   
##  [69] bit64_4.0.5                   KEGGREST_1.38.0              
##  [71] mime_0.12                     pracma_2.4.2                 
##  [73] xml2_1.3.3                    biomaRt_2.54.0               
##  [75] compiler_4.2.1                filelock_1.0.2               
##  [77] curl_4.3.3                    png_0.1-7                    
##  [79] interactiveDisplayBase_1.36.0 tibble_3.1.8                 
##  [81] bslib_0.4.0                   stringi_1.7.8                
##  [83] highr_0.9                     lattice_0.20-45              
##  [85] ProtGenerics_1.30.0           Matrix_1.5-1                 
##  [87] vctrs_0.5.0                   pillar_1.8.1                 
##  [89] lifecycle_1.0.3               BiocManager_1.30.19          
##  [91] jquerylib_0.1.4               data.table_1.14.4            
##  [93] bitops_1.0-7                  httpuv_1.6.6                 
##  [95] rtracklayer_1.58.0            R6_2.5.1                     
##  [97] BiocIO_1.8.0                  promises_1.2.0.1             
##  [99] KernSmooth_2.23-20            codetools_0.2-18             
## [101] MASS_7.3-58.1                 assertthat_0.2.1             
## [103] rjson_0.2.21                  withr_2.5.0                  
## [105] GenomicAlignments_1.34.0      Rsamtools_2.14.0             
## [107] multcomp_1.4-20               GenomeInfoDbData_1.2.9       
## [109] parallel_4.2.1                hms_1.1.2                    
## [111] rmarkdown_2.17                shiny_1.7.3                  
## [113] restfulr_0.0.15