diff --git a/README.md b/README.md index 86fd5697..7cbe8d46 100644 --- a/README.md +++ b/README.md @@ -34,12 +34,14 @@ data to biological insights. ### Filtering and Selection Functions -| Function | Description | -|-----------------------|-----------------------------------------------------| -| `identify_abundant()` | Identify abundant transcripts without removing them | -| `keep_abundant()` | Keep abundant transcripts | -| `keep_variable()` | Keep variable transcripts | -| `filterByExpr()` | Filter by expression | +| Function | Description | +|-----------------------------------------|------------------------------------------------------------------------| +| `identify_abundant()` | Identify abundant transcripts without removing them | +| `keep_abundant()` | Keep abundant transcripts | +| `identify_abundant_per_category()` | Identify abundant transcripts per experimental category | +| `keep_abundant_per_category()` | Keep transcripts abundant in at least one experimental category | +| `keep_variable()` | Keep variable transcripts | +| `filterByExpr()` | Filter by expression | ### Dimensionality Reduction Functions diff --git a/inst/NEWS.rd b/inst/NEWS.rd index b72a4f6c..c6009c82 100644 --- a/inst/NEWS.rd +++ b/inst/NEWS.rd @@ -4,6 +4,7 @@ \section{Changes in development version}{ \itemize{ \item \strong{BREAKING CHANGE:} \code{reduce_dimensions()} no longer accepts \code{.abundance} (a tidy-eval assay symbol). The assay must be given explicitly as a character string via \code{assay = "..."}, matching \code{assayNames(object)}. Omitting \code{assay} or passing an unknown assay name raises an error, so users must deliberately choose which abundance matrix to use (commonly a library-size scaled assay from \code{scale_abundance()}). Update all calls, tests, and vignettes that used \code{.abundance = counts} to \code{assay = "counts"} (or the appropriate assay name). + \item Added \code{identify_abundant_per_category()} and \code{keep_abundant_per_category()} as category-aware alternatives to \code{identify_abundant()} and \code{keep_abundant()}. These functions evaluate abundance within each experimental group independently and retain features that are sufficiently expressed in at least one category, making it possible to preserve group-specific or condition-specific signals that would be lost by global abundance filtering. The grouping is specified via \code{formula_design} (e.g. \code{~dex}), and the minimum number of categories required is controlled by \code{minimum_category} (default 1). }} \section{Changes in version 1.2.0, Bioconductor 3.12 Release}{ diff --git a/vignettes/introduction.Rmd b/vignettes/introduction.Rmd index 92bcdadb..fde0b711 100644 --- a/vignettes/introduction.Rmd +++ b/vignettes/introduction.Rmd @@ -47,6 +47,8 @@ Function | Description ------------ | ------------- `identify_abundant()` | Identify abundant transcripts without removing them `keep_abundant()` | Keep abundant transcripts +`identify_abundant_per_category()` | Identify abundant transcripts per experimental category +`keep_abundant_per_category()` | Keep transcripts abundant in at least one experimental category `keep_variable()` | Keep variable transcripts `filterByExpr()` | Filter by expression @@ -328,6 +330,95 @@ airway = airway_abundant_formula > **Tip:** Use `formula_design` for complex designs, and use the CPM threshold for library-size-aware filtering. +### Category-Aware Abundance Filtering + +Standard abundance filtering (`keep_abundant()`) asks whether a feature is sufficiently expressed across *all* samples. This works well when the signal of interest is globally abundant. However, it can silently discard features that are strongly expressed in only one experimental group — for example, a gene that is highly induced by treatment but nearly absent in the control condition. + +`keep_abundant_per_category()` addresses this by evaluating abundance *within each experimental category independently*. A feature is retained if it satisfies the abundance threshold in **at least one** category, regardless of its expression in the other categories. + +This approach is particularly important for: + +- **Case/control studies** where disease biomarkers are expressed only in the case group +- **Treatment-response experiments** where genes are induced by a specific stimulus +- **Multi-tissue or multi-cell-type datasets** with tissue-specific or cell-type-specific transcripts +- **Subgroup discovery** analyses where group-specific signals would be lost by global filtering + +#### Why this matters + +Consider a gene expressed at high levels in the `trt` group but near-zero in `untrt`. Global filtering may remove it because its average across all samples is low. But that gene is exactly what you want to detect in a differential expression test. Category-aware filtering preserves it. + +#### Run per-category filtering + +```{r filtering-per-category} +# Keep genes abundant in at least one treatment group +airway_abundant_per_category = airway_abundant_formula |> + keep_abundant_per_category( + formula_design = ~dex, # experimental grouping variable + minimum_counts = 10, # minimum raw count threshold + minimum_proportion = 0.7 # fraction of samples in a category that must pass + ) +``` + +#### Using CPM threshold instead of raw counts + +```{r filtering-per-category-cpm} +# Use counts-per-million threshold for library-size-aware filtering +airway_abundant_per_category_cpm = airway_abundant_formula |> + keep_abundant_per_category( + formula_design = ~dex, + minimum_count_per_million = 1, + minimum_proportion = 0.7 + ) +``` + +#### Require abundance in multiple categories + +By default (`minimum_category = 1`) a feature is kept if it passes the threshold in at least one category. Set `minimum_category = 2` to require abundance in at least two categories: + +```{r filtering-per-category-min2} +airway_abundant_per_category_min2 = airway_abundant_formula |> + keep_abundant_per_category( + formula_design = ~dex, + minimum_counts = 10, + minimum_proportion = 0.7, + minimum_category = 2 # must be abundant in both categories + ) +``` + +#### Mark features without filtering + +If you only want to annotate features with a logical `.abundant` column without immediately filtering, use the companion function `identify_abundant_per_category()`: + +```{r filtering-per-category-identify} +airway_marked = airway_abundant_formula |> + identify_abundant_per_category(formula_design = ~dex) +``` + +#### Compare global vs. category-aware filtering + +```{r filtering-compare-global-vs-category} +# Summary statistics +airway_abundant_formula |> summarise( + filtering = "global", + n_features = n_distinct(.feature) +) + +airway_abundant_per_category |> summarise( + filtering = "per category (>=1)", + n_features = n_distinct(.feature) +) + +airway_abundant_per_category_min2 |> summarise( + filtering = "per category (>=2)", + n_features = n_distinct(.feature) +) +``` + +> **When to use category-aware filtering:** +> Use `keep_abundant_per_category()` when your experimental design contains distinct biological groups and you want to retain features that are informative for at least one of them. Use global `keep_abundant()` when you are interested only in features that are broadly expressed across all samples. + + + ### Remove Redundant Transcripts Redundancy removal is a standard approach for reducing highly correlated features.