feat: geom_aggrcoverage

Author: js2264

DESCRIPTION

@@ -35,6 +35,7 @@ Imports:
     rtracklayer,
     methods,
     tidyr, 
+    ggplot2,
     dplyr,
     fansi, 
     pillar, 
@@ -46,7 +47,6 @@ Imports:
 Suggests: 
     tidySummarizedExperiment,
     plyranges,
-    ggplot2,
     TxDb.Hsapiens.UCSC.hg19.knownGene,
     AnnotationHub,
     GenomicFeatures, 

NAMESPACE

@@ -11,6 +11,7 @@ export(aggregate)
 export(as_tibble)
 export(coarsen)
 export(expand)
+export(geom_aggrcoverage)
 export(show)
 exportMethods(CoverageExperiment)
 exportMethods(aggregate)
@@ -19,6 +20,7 @@ exportMethods(show)
 import(GenomeInfoDb)
 import(GenomicRanges)
 import(SummarizedExperiment)
+import(ggplot2)
 import(methods)
 importFrom(BiocIO,resource)
 importFrom(BiocParallel,bplapply)

R/geoms.R

@@ -0,0 +1,87 @@
+#' geom_aggrcoverage
+#'
+#' #' @description
+#' 
+#' `geom_aggrcoverage()`
+#' 
+#' @name geom_aggrcoverage
+#' @rdname geom_aggrcoverage
+#' 
+#' @return A `ggplot` object`
+#'
+#' @import ggplot2
+#'
+#' @examples
+#' library(rtracklayer)
+#' library(plyranges)
+#' library(ggplot2)
+#' TSSs_bed <- system.file("extdata", "TSSs.bed", package = "tidyCoverage")
+#' features <- list(
+#'     TSS_fwd = import(TSSs_bed) |> filter(strand == '+'), 
+#'     TSS_rev = import(TSSs_bed) |> filter(strand == '-')
+#' )
+#' tracks <- list(
+#'     RNA_fwd = system.file("extdata", "RNA.fwd.bw", package = "tidyCoverage"),
+#'     RNA_rev = system.file("extdata", "RNA.rev.bw", package = "tidyCoverage")
+#' ) |> map(import, as = 'Rle')
+#' df <- CoverageExperiment(tracks, features, width = 5000, ignore.strand = FALSE) |> 
+#'  aggregate() |> 
+#'  as_tibble()
+NULL 
+
+GeomAggrCoverage <- ggplot2::ggproto("GeomAggrCoverage", ggplot2::Geom,
+    setup_params = function(data, params) {
+        params$ci <- params$ci
+        params
+    },
+    extra_params = c("na.rm"),
+    required_aes = c("x", "y", "ymin", "ymax"), 
+    default_aes = ggplot2::aes(
+        colour = "black", 
+        linewidth = 1, 
+        linetype = 1, 
+        alpha = 0.4
+    ),
+    draw_group = function(data, params, coord, ci = TRUE, ...) {
+        forLine <- transform(data, alpha = 1)
+        forRibbon <- transform(
+            data, 
+            alpha = data$alpha, 
+            fill = data$colour, 
+            colour = NA
+        )
+        grid::gList(
+            if (ci) ggplot2::GeomRibbon$draw_panel(forRibbon, params, coord, ...),
+            ggplot2::GeomLine$draw_panel(forLine, params, coord, ...)
+        )
+    }, 
+    draw_key = ggplot2::draw_key_smooth
+)
+
+#' @rdname geom_aggrcoverage
+#' @export
+
+geom_aggrcoverage <- function(
+    mapping = NULL, 
+    data = NULL, 
+    ..., 
+    ci = TRUE, 
+    na.rm = FALSE, 
+    show.legend = NA, 
+    inherit.aes = TRUE
+) {
+    m <- ggplot2::aes(x = coord, y = mean, ymin = ci_low, ymax = ci_high, group = interaction(.sample, .feature))
+    if (!is.null(mapping)) m <- utils::modifyList(m, mapping)
+
+    ggplot2::layer(
+        data = data, 
+        mapping = m,  
+        stat = "identity", 
+        geom = GeomAggrCoverage, 
+        position = "identity", 
+        show.legend = show.legend, 
+        inherit.aes = inherit.aes,
+        params = list(na.rm = na.rm, ci = ci, ...)
+    )
+}
+

man/geom_aggrcoverage.Rd

@@ -403,39 +403,27 @@ aggregated coverage signal of multiple tracks over multiple sets of genomic rang
 ```{r}
 AC |> 
     as_tibble() |> 
-    ggplot(aes(x = coord, y = mean, group = interaction(features, track), col = track)) + 
-    geom_line()
+    ggplot() + 
+    geom_aggrcoverage()
 ```
 
-This plot is way too busy. Let's first split into facets using `features`:
+Oopsie, a little busy here. Let's color by tracks and split facets by `features`:
 
 ```{r}
 AC |> 
     as_tibble() |> 
-    ggplot(aes(x = coord, y = mean, col = track)) + 
-    geom_line() + 
+    ggplot(aes(col = track)) + 
+    geom_aggrcoverage() + 
     facet_grid(features ~ .)
 ```
 
-Better, but what about adding confidence interval for each coverage track: 
+Nearly there, few cosmethic changes and we are done! 
 
 ```{r}
 AC |> 
     as_tibble() |> 
-    ggplot(aes(x = coord, y = mean)) + 
-    geom_ribbon(aes(ymin = ci_low, ymax = ci_high, fill = track), alpha = 0.2) + 
-    geom_line(aes(col = track)) + 
-    facet_grid(features ~ .)
-```
-
-Nearly there, few cosmethic changes and we're done! 
-
-```{r}
-AC |> 
-    as_tibble() |> 
-    ggplot(aes(x = coord, y = mean)) + 
-    geom_ribbon(aes(ymin = ci_low, ymax = ci_high, fill = track), alpha = 0.2) + 
-    geom_line(aes(col = track)) + 
+    ggplot(aes(col = track, linetype = track %in% c('RNA_fwd', 'RNA_rev'))) + 
+    geom_aggrcoverage() + 
     facet_grid(features ~ .) + 
     labs(x = 'Distance from genomic feature', y = 'Mean coverage (± 95% conf. intervale)') + 
     theme_bw() + 
@@ -470,9 +458,8 @@ if the `tidySummarizedExperiment` package id loaded.
 
 ```{r}
 AC |> 
-    ggplot(aes(x = coord, y = mean)) + 
-    geom_line(aes(col = track)) + 
-    facet_grid(track ~ .) + 
+    ggplot() + 
+    geom_aggrcoverage() + 
     labs(x = 'Distance from locus of convergent transcription', y = 'Scc1 coverage') + 
     theme_bw() + 
     theme(legend.position = 'top')
@@ -487,10 +474,8 @@ AC |>
 CoverageExperiment(tracks, features, width = 5000, scale = TRUE, center = TRUE) |> 
     filter(track == 'RNA_fwd') |> 
     aggregate(bin = 20) |> 
-    ggplot(aes(x = coord, y = mean)) + 
-    geom_ribbon(aes(ymin = ci_low, ymax = ci_high, fill = features), alpha = 0.2) + 
-    geom_line(aes(col = features)) + 
-    facet_grid(features ~ .) + 
+    ggplot(col = features) + 
+    geom_aggrcoverage(aes(col = features)) + 
     labs(x = 'Distance to center of genomic features', y = 'Forward RNA-seq coverage') + 
     theme_bw() + 
     theme(legend.position = 'top')
@@ -538,9 +523,8 @@ CoverageExperiment(
     scale = TRUE, center = TRUE
 ) |> 
     aggregate() |> 
-    ggplot(aes(x = coord, y = mean)) + 
-    geom_ribbon(aes(ymin = ci_low, ymax = ci_high, fill = track), alpha = 0.2) + 
-    geom_line(aes(col = track)) + 
+    ggplot() + 
+    geom_aggrcoverage(aes(col = track)) + 
     facet_grid(track ~ .) + 
     labs(x = 'Distance from TSSs', y = 'Mean coverage') + 
     theme_bw() + 
@@ -586,9 +570,8 @@ AC |>
             stringr::str_detect(track, 'H3') ~ "H3"
         )
     ) |> 
-    ggplot(aes(x = coord, y = mean)) + 
-    geom_ribbon(aes(ymin = ci_low, ymax = ci_high, fill = track), alpha = 0.2) + 
-    geom_line(aes(col = track)) + 
+    ggplot() + 
+    geom_aggrcoverage(aes(col = track)) + 
     facet_grid(~histone) + 
     labs(x = 'Distance from TSSs', y = 'Mean histone PTM coverage') + 
     theme_bw() +