diff --git a/r/R/genomeAnnotations.R b/r/R/genomeAnnotations.R
index 158a6378..33620ee7 100644
--- a/r/R/genomeAnnotations.R
+++ b/r/R/genomeAnnotations.R
@@ -328,6 +328,16 @@ canonical_gene_symbol <- function(query, gene_mapping = human_gene_mapping) {
 #' @param extend_bp Bases to extend region upstream and downstream of gene. If length 1, extension is 
 #'     symmetric. If length 2, provide upstream extension then downstream extension as positive distances.
 #' @return List of chr, start, end positions for use with trackplot functions.
+#' @examples
+#' ## Prep data
+#' genes <- read_gencode_transcripts(
+#'   file.path(tempdir(), "references"), release = "42",
+#'   annotation_set = "basic",
+#'   features = "transcript"
+#' )
+#' 
+#' ## Get gene region
+#' gene_region(genes, "CD19", extend_bp = 1e5)
 #' @export
 gene_region <- function(genes, gene_symbol, extend_bp = c(1e4, 1e4), gene_mapping = human_gene_mapping) {
   genes <- normalize_ranges(genes, metadata_cols = c("strand", "gene_name"))
diff --git a/r/R/plots.R b/r/R/plots.R
index 5da3aaeb..e37ad784 100644
--- a/r/R/plots.R
+++ b/r/R/plots.R
@@ -108,6 +108,13 @@ continuous_palette <- function(name) {
 #' @param return_data If true, return data from just before plotting rather than a plot.
 #' @param apply_styling If false, return a plot without pretty styling applied
 #' @return ggplot2 plot object
+#' @examples
+#' ## Prep data
+#' mat <- get_demo_mat(filter_qc = FALSE, subset = FALSE)
+#' reads_per_cell <- colSums(mat)
+#' 
+#' # Render knee plot
+#' plot_read_count_knee(reads_per_cell, cutoff = 1e3)
 #' @export
 plot_read_count_knee <- function(read_counts, cutoff = NULL, return_data = FALSE, apply_styling = TRUE) {
   # Keep ~1k cells per order of magnitude to speed up plotting
@@ -178,6 +185,20 @@ plot_read_count_knee <- function(read_counts, cutoff = NULL, return_data = FALSE
 #' @param min_tss Minimum TSS Enrichment cutoff
 #' @param bins Number of bins for density calculation
 #' @inheritParams plot_embedding
+#' @examples
+#' ## Prep data
+#' frags <- get_demo_frags(filter_qc = FALSE, subset = FALSE)
+#' genes <- read_gencode_transcripts(
+#'   file.path(tempdir(), "references"), release = "42",
+#'   annotation_set = "basic",
+#'   features = "transcript"
+#' )
+#' blacklist <- read_encode_blacklist(file.path(tempdir(), "references"), genome="hg38")
+#' atac_qc <- qc_scATAC(frags, genes, blacklist)
+#' 
+#' 
+#' ## Render tss enrichment vs fragment plot
+#' plot_tss_scatter(atac_qc, min_frags = 1000, min_tss = 10)
 #' @export
 plot_tss_scatter <- function(atac_qc, min_frags = NULL, min_tss = NULL, bins = 100, apply_styling = TRUE) {
   assert_is(atac_qc, "data.frame")
@@ -255,6 +276,9 @@ plot_tss_scatter <- function(atac_qc, min_frags = NULL, min_tss = NULL, bins = 1
 #' @param max_length Maximum length to show on the plot
 #' @inheritParams plot_embedding
 #' @return Numeric vector where index i contans the number of length-i fragments
+#' @examples
+#' frags <- get_demo_frags(filter_qc = FALSE, subset = FALSE)
+#' plot_fragment_length(frags)
 #' @export
 plot_fragment_length <- function(fragments, max_length = 500, return_data = FALSE, apply_styling = TRUE) {
   assert_is(fragments, "IterableFragments")
@@ -297,6 +321,17 @@ plot_fragment_length <- function(fragments, max_length = 500, return_data = FALS
 #' @param genes Coordinate ranges for genes (must include strand)
 #' @param smooth Number of bases to smooth over (rolling average)
 #' @seealso `footprint()`, `plot_tf_footprint()`
+#' @examples
+#' ## Prep data
+#' frags <- get_demo_frags(filter_qc = FALSE, subset = FALSE)
+#' genes <- read_gencode_transcripts(
+#'   file.path(tempdir(), "references"), release = "42",
+#'   annotation_set = "basic",
+#'   features = "transcript"
+#' )
+#' 
+#' ## Plot tss profile
+#' plot_tss_profile(frags, genes)
 #' @export
 plot_tss_profile <- function(fragments, genes, cell_groups = rlang::rep_along(cellNames(fragments), "all"),
                              flank = 2000L, smooth = 0L, zero_based_coords = !is(genes, "GRanges"),
@@ -485,6 +520,48 @@ collect_features <- function(source, features = NULL, gene_mapping = human_gene_
 #' @return By default, returns a ggplot2 object with all the requested features plotted
 #' in a grid. If `return_data` or `return_plot_list` is called, the return value will
 #' match that argument.
+#' @examples
+## Prep data
+#' set.seed(123)
+#' mat <- get_demo_mat()
+#' ## Normalize matrix
+#' mat_norm <- log1p(multiply_cols(mat, 1/colSums(mat)) * 10000) %>% write_matrix_memory(compress = FALSE)
+#' ## Get variable genes
+#' stats <- matrix_stats(mat, row_stats = "variance")
+#' variable_genes <- order(stats$row_stats["variance",], decreasing=TRUE) %>% 
+#'   head(1000) %>% 
+#'   sort()
+#' # Z score normalize genes
+#' mat_norm <- mat[variable_genes, ]
+#' gene_means <- stats$row_stats['mean', variable_genes]
+#' gene_vars <- stats$row_stats['variance', variable_genes]
+#' mat_norm <- (mat_norm - gene_means) / gene_vars
+#' ## Save matrix to memory
+#' mat_norm <- mat_norm %>% write_matrix_memory(compress = FALSE)
+#' ## Run SVD
+#' svd <- BPCells::svds(mat_norm, k = 10)
+#' pca <- multiply_cols(svd$v, svd$d)
+#' ## Get UMAP
+#' umap <- uwot::umap(pca)
+#' ## Get clusters
+#' clusts <- knn_hnsw(pca, ef = 500) %>%
+#'   knn_to_snn_graph() %>%
+#'   cluster_graph_louvain()
+#' 
+#' 
+#' ## Plot embedding
+#' print(length(clusts))
+#' 
+#' plot_embedding(clusts, umap)
+#' 
+#' 
+#' ### Can also plot by features
+#' #plot_embedding(
+#' #  source = mat,
+#' #  umap,
+#' #  features = c("MS4A1", "CD3E"),
+#' #)
+#' 
 #' @export
 plot_embedding <- function(source, embedding, features = NULL,
                            quantile_range = c(0.01, 0.99),
@@ -741,6 +818,17 @@ rotate_x_labels <- function(degrees = 45) {
 #' @inheritParams plot_embedding
 #' @param gene_mapping An optional vector for gene name matching with match_gene_symbol().
 #' @param colors Color scale for plot
+#' @examples
+#' ## Prep data
+#' mat <- get_demo_mat()
+#' cell_types <- paste("Group", rep(1:3, length.out = length(colnames(mat))))
+#' 
+#' ## Plot dot
+#' plot <- plot_dot(mat, c("MS4A1", "CD3E"), cell_types)
+#' 
+#' BPCells:::render_plot_from_storage(
+#'   plot, width = 4, height = 5
+#' )
 #' @export
 plot_dot <- function(source, features, groups, group_order = NULL, gene_mapping = human_gene_mapping,
                      colors = c("lightgrey", "#4682B4"),
diff --git a/r/R/trackplots.R b/r/R/trackplots.R
index 998bfecd..5f04cbb9 100644
--- a/r/R/trackplots.R
+++ b/r/R/trackplots.R
@@ -252,6 +252,22 @@ trackplot_normalize_ranges_with_metadata <- function(data, metadata) {
   return(data)
 }
 
+#' Render a plot with intermediate disk storage step
+#' 
+#' Take a plotting object and save in temp storage, so it can be outputted with exact dimensions.
+#' Primarily used to allow for adjusting plot dimensions within function reference examples.
+#' @param plot (ggplot) ggplot output from a plotting function
+#' @param width (numeric) width of rendered plot
+#' @param height (numeric) height of rendered plot
+#' @keywords internal
+render_plot_from_storage <- function(plot, width, height) {
+  assert_is(plot, "ggplot")
+  image_path <- tempfile(fileext = ".png")
+  ggplot2::ggsave(image_path, plot, width = width, height = height)
+  img <- png::readPNG(image_path)
+  grid::grid.raster(img)
+}
+
 #' Combine track plots
 #' 
 #' Combines multiple track plots of the same region into a single grid.
@@ -268,6 +284,39 @@ trackplot_normalize_ranges_with_metadata <- function(data, metadata) {
 #'    the text label, y-axis, and plot body. The relative height of each row is given
 #'    by heights. A shared title and x-axis are put at the top.
 #' @seealso `trackplot_coverage()`, `trackplot_gene()`, `trackplot_loop()`, `trackplot_scalebar()`
+#' @examples
+#' ## Prep data
+#' frags <- get_demo_frags()
+#' 
+#' ## Use genes and blacklist to determine proper number of reads per cell
+#' genes <- read_gencode_transcripts(
+#'   file.path(tempdir(), "references"), release = "42",
+#'   annotation_set = "basic",
+#'   features = "transcript"
+#' )
+#' blacklist <- read_encode_blacklist(file.path(tempdir(), "references"), genome="hg38")
+#' read_counts <- qc_scATAC(frags, genes, blacklist)$nFrags
+#' region <- "chr4:3034877-4034877"
+#' cell_types <- paste("Group", rep(1:3, length.out = length(cellNames(frags))))
+#' transcripts <- read_gencode_transcripts(
+#'   file.path(tempdir(), "references"), release = "42",
+#'   annotation_set = "basic"
+#' )
+#' region <- "chr4:3034877-4034877"
+#' 
+#' 
+#' ## Get all trackplots and scalebars to combine
+#' plot_scalebar <- trackplot_scalebar(region)
+#' plot_gene <- trackplot_gene(transcripts, region)
+#' plot_coverage <- trackplot_coverage(frags, region, groups = cell_types, cell_read_counts = read_counts)
+#' 
+#' 
+#' ## Combine trackplots and render
+#' ## Also remove colors from gene track
+#' plot <- trackplot_combine(
+#'     list(plot_scalebar, plot_coverage, plot_gene + ggplot2::guides(color = "none"))
+#' )
+#' BPCells:::render_plot_from_storage(plot, width = 6, height = 4)
 #' @export
 trackplot_combine <- function(tracks, side_plot = NULL, title = NULL, side_plot_width = 0.3) {
   for (plot in tracks) {
@@ -407,6 +456,26 @@ trackplot_combine <- function(tracks, side_plot = NULL, title = NULL, side_plot_
 #' specify the labels for each track. If `return_data` or `return_plot_list` is
 #' `TRUE`, the return value will be modified accordingly.
 #' @seealso `trackplot_combine()`, `trackplot_gene()`, `trackplot_loop()`, `trackplot_scalebar()`
+#' @examples
+## Prep data
+#' frags <- get_demo_frags()
+#' 
+#' ## Use genes and blacklist to determine proper number of reads per cell
+#' genes <- read_gencode_transcripts(
+#'   file.path(tempdir(), "references"), release = "42",
+#'   annotation_set = "basic",
+#'   features = "transcript"
+#' )
+#' blacklist <- read_encode_blacklist(file.path(tempdir(), "references"), genome="hg38")
+#' read_counts <- qc_scATAC(frags, genes, blacklist)$nFrags
+#' region <- "chr4:3034877-4034877"
+#' cell_types <- paste("Group", rep(1:3, length.out = length(cellNames(frags))))
+#' 
+#' 
+#' BPCells:::render_plot_from_storage(
+#'   trackplot_coverage(frags, region, groups = cell_types, cell_read_counts = read_counts),
+#'   width = 6, height = 3
+#' )
 #' @export
 trackplot_coverage <- function(fragments, region, groups,
                            cell_read_counts,
@@ -510,6 +579,18 @@ trackplot_coverage <- function(fragments, region, groups,
 #' @param label_size size for transcript labels in units of mm
 #' @return Plot of gene locations
 #' @seealso `trackplot_combine()`, `trackplot_coverage()`, `trackplot_loop()`, `trackplot_scalebar()`
+#' @examples
+#' ## Prep data
+#' transcripts <- read_gencode_transcripts(
+#'   file.path(tempdir(), "references"), release = "42",
+#'   annotation_set = "basic", features = "transcript"
+#' )
+#' region <- "chr4:3034877-4034877"
+#' 
+#' 
+#' ## Plot gene trackplot
+#' plot <- trackplot_gene(transcripts, region)
+#' BPCells:::render_plot_from_storage(plot, width = 6, height = 1)
 #' @export
 trackplot_gene <- function(transcripts, region, exon_size = 2.5, gene_size = 0.5, label_size = 11*.8/ggplot2::.pt, track_label="Genes", return_data = FALSE) {
   region <- normalize_ranges(region)
@@ -607,6 +688,28 @@ trackplot_gene <- function(transcripts, region, exon_size = 2.5, gene_size = 0.5
 #' @param show_strand If TRUE, show strand direction as arrows
 #' @return Plot of genomic loci if return_data is FALSE, otherwise returns the data frame used to generate the plot
 #' @seealso `trackplot_combine()`, `trackplot_coverage()`, `trackplot_loop()`, `trackplot_scalebar()`, `trackplot_gene()`
+#' @examples
+#' ## Prep data
+#' ## Peaks generated from demo frags, as input into `call_peaks_tile()`
+#' peaks <- tibble::tibble(
+#'   chr = factor(rep("chr4", 16)),
+#'   start = c(3041400, 3041733, 3037400, 3041933, 3040466, 3041200, 
+#'             3038200, 3038000, 3040266, 3037733, 3040800, 3042133, 
+#'             3038466, 3037200, 3043333, 3040066),
+#'   end = c(3041600, 3041933, 3037600, 3042133, 3040666, 3041400, 
+#'           3038400, 3038200, 3040466, 3037933, 3041000, 3042333, 
+#'           3038666, 3037400, 3043533, 3040266),
+#'   enrichment = c(46.4, 43.5, 28.4, 27.3, 17.3, 11.7, 
+#'                  10.5, 7.95, 7.22, 6.86, 6.32, 6.14, 
+#'                  5.96, 5.06, 4.51, 3.43)
+#' )
+#' region <- "chr4:3034877-3044877"
+#' 
+#' ## Plot peaks
+#' BPCells:::render_plot_from_storage(
+#'   trackplot_genome_annotation(peaks, region, color_by = "enrichment"),
+#'   width = 6, height = 1
+#' )
 #' @export
 trackplot_genome_annotation <- function(loci, region, color_by = NULL, colors = NULL, label_by = NULL, label_size = 11*.8/ggplot2::.pt, show_strand=FALSE,
                                         annotation_size = 2.5, track_label="Peaks", return_data = FALSE) {
@@ -729,6 +832,19 @@ trackplot_genome_annotation <- function(loci, region, color_by = NULL, colors =
 #' 
 #' @return Plot of loops connecting genomic coordinates
 #' @seealso `trackplot_combine()`, `trackplot_coverage()`, `trackplot_gene()`, `trackplot_scalebar()`, `trackplot_genome_annotation()`
+#' @examples
+#' peaks <- c(3054877, 3334877, 3534877, 3634877, 3734877)
+#' loops <- tibble::tibble(
+#'   chr = "chr4",
+#'   start = peaks[c(1,1,2,3)],
+#'   end = peaks[c(2,3,4,5)],
+#'   score = c(4,1,3,2)
+#' )
+#' region <- "chr4:3034877-4034877"
+#' 
+#' ## Plot loops
+#' plot <- trackplot_loop(loops, region, color_by = "score")
+#' BPCells:::render_plot_from_storage(plot, width = 6, height = 1.5)
 #' @export
 trackplot_loop <- function(loops, region, color_by=NULL, colors=NULL, allow_truncated=TRUE, curvature=0.75, track_label="Links", return_data = FALSE) {
   region <- normalize_ranges(region)
@@ -826,6 +942,11 @@ trackplot_loop <- function(loops, region, color_by=NULL, colors=NULL, allow_trun
 #' 
 #' @return Plot with coordinates and scalebar for plotted genomic region
 #' @seealso `trackplot_combine()`, `trackplot_coverage()`, `trackplot_gene()`, `trackplot_loop()`
+#' @examples
+#' region <- "chr4:3034877-3044877"
+#' BPCells:::render_plot_from_storage(
+#'   trackplot_scalebar(region), width = 6, height = 1
+#' )
 #' @export
 trackplot_scalebar <- function(region, font_pt=11) {
   region <- normalize_ranges(region)
diff --git a/r/man/gene_region.Rd b/r/man/gene_region.Rd
index 81082017..4855da39 100644
--- a/r/man/gene_region.Rd
+++ b/r/man/gene_region.Rd
@@ -35,3 +35,14 @@ Conveniently look up the region of a gene by gene symbol. The value returned by
 can be used as the \code{region} argument for trackplot functions such as
 \code{trackplot_coverage()} or \code{trackplot_gene()}
 }
+\examples{
+## Prep data
+genes <- read_gencode_transcripts(
+  file.path(tempdir(), "references"), release = "42",
+  annotation_set = "basic",
+  features = "transcript"
+)
+
+## Get gene region
+gene_region(genes, "CD19", extend_bp = 1e5)
+}
diff --git a/r/man/plot_dot.Rd b/r/man/plot_dot.Rd
index 61bd55fe..e92da530 100644
--- a/r/man/plot_dot.Rd
+++ b/r/man/plot_dot.Rd
@@ -38,3 +38,15 @@ map to zero)}
 Plot feature levels per group or cluster as a grid of dots.
 Dots are colored by z-score normalized average expression, and sized by percent non-zero.
 }
+\examples{
+## Prep data
+mat <- get_demo_mat()
+cell_types <- paste("Group", rep(1:3, length.out = length(colnames(mat))))
+
+## Plot dot
+plot <- plot_dot(mat, c("MS4A1", "CD3E"), cell_types)
+
+BPCells:::render_plot_from_storage(
+  plot, width = 4, height = 5
+)
+}
diff --git a/r/man/plot_embedding.Rd b/r/man/plot_embedding.Rd
index 10d4dc9c..281c49ae 100644
--- a/r/man/plot_embedding.Rd
+++ b/r/man/plot_embedding.Rd
@@ -96,3 +96,45 @@ are repeatedly multiplied by the smoothing matrix and re-scaled so the average
 value stays the same.
 }
 }
+\examples{
+set.seed(123)
+mat <- get_demo_mat()
+## Normalize matrix
+mat_norm <- log1p(multiply_cols(mat, 1/colSums(mat)) * 10000) \%>\% write_matrix_memory(compress = FALSE)
+## Get variable genes
+stats <- matrix_stats(mat, row_stats = "variance")
+variable_genes <- order(stats$row_stats["variance",], decreasing=TRUE) \%>\% 
+  head(1000) \%>\% 
+  sort()
+# Z score normalize genes
+mat_norm <- mat[variable_genes, ]
+gene_means <- stats$row_stats['mean', variable_genes]
+gene_vars <- stats$row_stats['variance', variable_genes]
+mat_norm <- (mat_norm - gene_means) / gene_vars
+## Save matrix to memory
+mat_norm <- mat_norm \%>\% write_matrix_memory(compress = FALSE)
+## Run SVD
+svd <- BPCells::svds(mat_norm, k = 10)
+pca <- multiply_cols(svd$v, svd$d)
+## Get UMAP
+umap <- uwot::umap(pca)
+## Get clusters
+clusts <- knn_hnsw(pca, ef = 500) \%>\%
+  knn_to_snn_graph() \%>\%
+  cluster_graph_louvain()
+
+
+## Plot embedding
+print(length(clusts))
+
+plot_embedding(clusts, umap)
+
+
+### Can also plot by features
+#plot_embedding(
+#  source = mat,
+#  umap,
+#  features = c("MS4A1", "CD3E"),
+#)
+
+}
diff --git a/r/man/plot_fragment_length.Rd b/r/man/plot_fragment_length.Rd
index c86f0ca4..0931b560 100644
--- a/r/man/plot_fragment_length.Rd
+++ b/r/man/plot_fragment_length.Rd
@@ -29,3 +29,7 @@ x-axis, and proportion of fragments on the y-axis. Typical plots will show
 10-basepair periodicity, as well as humps spaced at multiples of a nucleosome
 width (about 150bp).
 }
+\examples{
+frags <- get_demo_frags(filter_qc = FALSE, subset = FALSE)
+plot_fragment_length(frags)
+}
diff --git a/r/man/plot_read_count_knee.Rd b/r/man/plot_read_count_knee.Rd
index 3fe07c7b..b4ffbaad 100644
--- a/r/man/plot_read_count_knee.Rd
+++ b/r/man/plot_read_count_knee.Rd
@@ -29,3 +29,11 @@ Plots read count rank vs. number of reads on a log-log scale.
 \details{
 Performs logarithmic downsampling to reduce the number of points plotted
 }
+\examples{
+## Prep data
+mat <- get_demo_mat(filter_qc = FALSE, subset = FALSE)
+reads_per_cell <- colSums(mat)
+
+# Render knee plot
+plot_read_count_knee(reads_per_cell, cutoff = 1e3)
+}
diff --git a/r/man/plot_tss_profile.Rd b/r/man/plot_tss_profile.Rd
index e2f99919..08a5334f 100644
--- a/r/man/plot_tss_profile.Rd
+++ b/r/man/plot_tss_profile.Rd
@@ -40,6 +40,18 @@ Plot the enrichmment of insertions relative to transcription start sites (TSS).
 Typically, this plot shows strong enrichment of insertions near a TSS, and a
 small bump downstream around 220bp downstream of the TSS for the +1 nucleosome.
 }
+\examples{
+## Prep data
+frags <- get_demo_frags(filter_qc = FALSE, subset = FALSE)
+genes <- read_gencode_transcripts(
+  file.path(tempdir(), "references"), release = "42",
+  annotation_set = "basic",
+  features = "transcript"
+)
+
+## Plot tss profile
+plot_tss_profile(frags, genes)
+}
 \seealso{
 \code{footprint()}, \code{plot_tf_footprint()}
 }
diff --git a/r/man/plot_tss_scatter.Rd b/r/man/plot_tss_scatter.Rd
index 94d3ddb0..9c75ffd4 100644
--- a/r/man/plot_tss_scatter.Rd
+++ b/r/man/plot_tss_scatter.Rd
@@ -28,3 +28,18 @@ Density scatter plot with log10(fragment_count) on the x-axis and TSS
 enrichment on the y-axis. This plot is most useful to select which cell barcodes
 in an experiment correspond to high-quality cells
 }
+\examples{
+## Prep data
+frags <- get_demo_frags(filter_qc = FALSE, subset = FALSE)
+genes <- read_gencode_transcripts(
+  file.path(tempdir(), "references"), release = "42",
+  annotation_set = "basic",
+  features = "transcript"
+)
+blacklist <- read_encode_blacklist(file.path(tempdir(), "references"), genome="hg38")
+atac_qc <- qc_scATAC(frags, genes, blacklist)
+
+
+## Render tss enrichment vs fragment plot
+plot_tss_scatter(atac_qc, min_frags = 1000, min_tss = 10)
+}
diff --git a/r/man/render_plot_from_storage.Rd b/r/man/render_plot_from_storage.Rd
new file mode 100644
index 00000000..bfcea1f6
--- /dev/null
+++ b/r/man/render_plot_from_storage.Rd
@@ -0,0 +1,20 @@
+% Generated by roxygen2: do not edit by hand
+% Please edit documentation in R/trackplots.R
+\name{render_plot_from_storage}
+\alias{render_plot_from_storage}
+\title{Render a plot with intermediate disk storage step}
+\usage{
+render_plot_from_storage(plot, width, height)
+}
+\arguments{
+\item{plot}{(ggplot) ggplot output from a plotting function}
+
+\item{width}{(numeric) width of rendered plot}
+
+\item{height}{(numeric) height of rendered plot}
+}
+\description{
+Take a plotting object and save in temp storage, so it can be outputted with exact dimensions.
+Primarily used to allow for adjusting plot dimensions within function reference examples.
+}
+\keyword{internal}
diff --git a/r/man/trackplot_combine.Rd b/r/man/trackplot_combine.Rd
index 6e0bc283..a18d0b34 100644
--- a/r/man/trackplot_combine.Rd
+++ b/r/man/trackplot_combine.Rd
@@ -32,6 +32,40 @@ by heights. A shared title and x-axis are put at the top.
 Combines multiple track plots of the same region into a single grid.
 Uses the \code{patchwork} package to perform the alignment.
 }
+\examples{
+## Prep data
+frags <- get_demo_frags()
+
+## Use genes and blacklist to determine proper number of reads per cell
+genes <- read_gencode_transcripts(
+  file.path(tempdir(), "references"), release = "42",
+  annotation_set = "basic",
+  features = "transcript"
+)
+blacklist <- read_encode_blacklist(file.path(tempdir(), "references"), genome="hg38")
+read_counts <- qc_scATAC(frags, genes, blacklist)$nFrags
+region <- "chr4:3034877-4034877"
+cell_types <- paste("Group", rep(1:3, length.out = length(cellNames(frags))))
+transcripts <- read_gencode_transcripts(
+  file.path(tempdir(), "references"), release = "42",
+  annotation_set = "basic"
+)
+region <- "chr4:3034877-4034877"
+
+
+## Get all trackplots and scalebars to combine
+plot_scalebar <- trackplot_scalebar(region)
+plot_gene <- trackplot_gene(transcripts, region)
+plot_coverage <- trackplot_coverage(frags, region, groups = cell_types, cell_read_counts = read_counts)
+
+
+## Combine trackplots and render
+## Also remove colors from gene track
+plot <- trackplot_combine(
+    list(plot_scalebar, plot_coverage, plot_gene + ggplot2::guides(color = "none"))
+)
+BPCells:::render_plot_from_storage(plot, width = 6, height = 4)
+}
 \seealso{
 \code{trackplot_coverage()}, \code{trackplot_gene()}, \code{trackplot_loop()}, \code{trackplot_scalebar()}
 }
diff --git a/r/man/trackplot_coverage.Rd b/r/man/trackplot_coverage.Rd
index c6225c9f..9abf0fe1 100644
--- a/r/man/trackplot_coverage.Rd
+++ b/r/man/trackplot_coverage.Rd
@@ -58,6 +58,26 @@ specify the labels for each track. If \code{return_data} or \code{return_plot_li
 Plot a pseudobulk genome track, showing the number of fragment insertions
 across a region for each cell type or group.
 }
+\examples{
+frags <- get_demo_frags()
+
+## Use genes and blacklist to determine proper number of reads per cell
+genes <- read_gencode_transcripts(
+  file.path(tempdir(), "references"), release = "42",
+  annotation_set = "basic",
+  features = "transcript"
+)
+blacklist <- read_encode_blacklist(file.path(tempdir(), "references"), genome="hg38")
+read_counts <- qc_scATAC(frags, genes, blacklist)$nFrags
+region <- "chr4:3034877-4034877"
+cell_types <- paste("Group", rep(1:3, length.out = length(cellNames(frags))))
+
+
+BPCells:::render_plot_from_storage(
+  trackplot_coverage(frags, region, groups = cell_types, cell_read_counts = read_counts),
+  width = 6, height = 3
+)
+}
 \seealso{
 \code{trackplot_combine()}, \code{trackplot_gene()}, \code{trackplot_loop()}, \code{trackplot_scalebar()}
 }
diff --git a/r/man/trackplot_gene.Rd b/r/man/trackplot_gene.Rd
index 7991e14a..96285492 100644
--- a/r/man/trackplot_gene.Rd
+++ b/r/man/trackplot_gene.Rd
@@ -47,6 +47,19 @@ Plot of gene locations
 \description{
 Plot transcript models
 }
+\examples{
+## Prep data
+transcripts <- read_gencode_transcripts(
+  file.path(tempdir(), "references"), release = "42",
+  annotation_set = "basic", features = "transcript"
+)
+region <- "chr4:3034877-4034877"
+
+
+## Plot gene trackplot
+plot <- trackplot_gene(transcripts, region)
+BPCells:::render_plot_from_storage(plot, width = 6, height = 1)
+}
 \seealso{
 \code{trackplot_combine()}, \code{trackplot_coverage()}, \code{trackplot_loop()}, \code{trackplot_scalebar()}
 }
diff --git a/r/man/trackplot_genome_annotation.Rd b/r/man/trackplot_genome_annotation.Rd
index 71af5964..f5bf7129 100644
--- a/r/man/trackplot_genome_annotation.Rd
+++ b/r/man/trackplot_genome_annotation.Rd
@@ -47,6 +47,29 @@ Plot of genomic loci if return_data is FALSE, otherwise returns the data frame u
 \description{
 Plot range-based annotation tracks (e.g. peaks)
 }
+\examples{
+## Prep data
+## Peaks generated from demo frags, as input into `call_peaks_tile()`
+peaks <- tibble::tibble(
+  chr = factor(rep("chr4", 16)),
+  start = c(3041400, 3041733, 3037400, 3041933, 3040466, 3041200, 
+            3038200, 3038000, 3040266, 3037733, 3040800, 3042133, 
+            3038466, 3037200, 3043333, 3040066),
+  end = c(3041600, 3041933, 3037600, 3042133, 3040666, 3041400, 
+          3038400, 3038200, 3040466, 3037933, 3041000, 3042333, 
+          3038666, 3037400, 3043533, 3040266),
+  enrichment = c(46.4, 43.5, 28.4, 27.3, 17.3, 11.7, 
+                 10.5, 7.95, 7.22, 6.86, 6.32, 6.14, 
+                 5.96, 5.06, 4.51, 3.43)
+)
+region <- "chr4:3034877-3044877"
+
+## Plot peaks
+BPCells:::render_plot_from_storage(
+  trackplot_genome_annotation(peaks, region, color_by = "enrichment"),
+  width = 6, height = 1
+)
+}
 \seealso{
 \code{trackplot_combine()}, \code{trackplot_coverage()}, \code{trackplot_loop()}, \code{trackplot_scalebar()}, \code{trackplot_gene()}
 }
diff --git a/r/man/trackplot_loop.Rd b/r/man/trackplot_loop.Rd
index 188517f9..5ebbb94b 100644
--- a/r/man/trackplot_loop.Rd
+++ b/r/man/trackplot_loop.Rd
@@ -41,6 +41,20 @@ Plot of loops connecting genomic coordinates
 \description{
 Plot loops
 }
+\examples{
+peaks <- c(3054877, 3334877, 3534877, 3634877, 3734877)
+loops <- tibble::tibble(
+  chr = "chr4",
+  start = peaks[c(1,1,2,3)],
+  end = peaks[c(2,3,4,5)],
+  score = c(4,1,3,2)
+)
+region <- "chr4:3034877-4034877"
+
+## Plot loops
+plot <- trackplot_loop(loops, region, color_by = "score")
+BPCells:::render_plot_from_storage(plot, width = 6, height = 1.5)
+}
 \seealso{
 \code{trackplot_combine()}, \code{trackplot_coverage()}, \code{trackplot_gene()}, \code{trackplot_scalebar()}, \code{trackplot_genome_annotation()}
 }
diff --git a/r/man/trackplot_scalebar.Rd b/r/man/trackplot_scalebar.Rd
index 9de670db..f52f663e 100644
--- a/r/man/trackplot_scalebar.Rd
+++ b/r/man/trackplot_scalebar.Rd
@@ -18,6 +18,12 @@ Plot with coordinates and scalebar for plotted genomic region
 \description{
 Plots a human-readable scale bar and coordinates of the region being plotted
 }
+\examples{
+region <- "chr4:3034877-3044877"
+BPCells:::render_plot_from_storage(
+  trackplot_scalebar(region), width = 6, height = 1
+)
+}
 \seealso{
 \code{trackplot_combine()}, \code{trackplot_coverage()}, \code{trackplot_gene()}, \code{trackplot_loop()}
 }