Skip to contents

Tests all other groups against a single reference element using pairwise Wilcoxon tests with Benjamini-Hochberg (BH) correction. This is the statistically appropriate approach for targeted biological hypotheses (e.g., "Is each region type enriched compared to intergenic?").

Usage

av_compare_against_reference(
  obj,
  group_by,
  compare_by,
  reference_level,
  abs_score = TRUE,
  min_n = 5
)

Arguments

obj

An AlphaVarSet object.

group_by

Column name defining the groups (e.g., "region_type").

compare_by

Column name defining the items (e.g., "gene_name"). Must be a categorical column (character or factor).

reference_level

Character. The reference group label to test against (e.g., "intergenic"). Must exist in the group_by column.

abs_score

Logical. If TRUE (default), uses absolute scores.

min_n

Minimum sample size per group per item (default 5).

Value

The updated AlphaVarSet object with results stored in obj$stats$comparisons[[comparison_name]] (named "<group_by>_vs_<reference_level>_by_<compare_by>"). Structure matches av_compare_groups() output (method, global_p, results_tidy, params, dropped_groups).

Details

Why Many-to-One instead of All-vs-All? When the biological question targets enrichment relative to a neutral baseline (e.g., intergenic), testing all \(k(k-1)/2\) pairs inflates the multiple testing burden unnecessarily. Many-to-One tests only \(k-1\) comparisons per item, maximizing statistical power under BH correction.

Effect Size: Rank-Biserial Correlation (r). Positive r indicates the non-reference group is stochastically larger than the reference.

Examples

example_file <- system.file("extdata", "bcl11a_rs1427407_example.csv", package = "AlphaVaR")
obj <- av_import_csv(example_file)
#> Reading bcl11a_rs1427407_example.csv ...

# Cell lines as the baseline: k - 1 tests per modality instead of all pairs.
obj <- av_compare_against_reference(
  obj,
  group_by = "biosample_type",
  compare_by = "output_type",
  reference_level = "cell_line"
)
#> Warning: Dropped 4 of 2676 rows with NA in 'biosample_type', 'output_type', or 'quantile_score'.
#> Warning: 1 items skipped (reference group n < 5 ): SPLICE_SITES 
res <- obj$stats$comparisons[["biosample_type_vs_cell_line_by_output_type"]]
res$method
#> [1] "Many-to-One Wilcoxon vs cell_line (BH)"
head(res$results_tidy)
#> # A tibble: 6 × 16
#>   item         group_1 group_2   n_1   n_2 median_1 median_2 median_diff   p_val
#>   <chr>        <chr>   <chr>   <int> <int>    <dbl>    <dbl>       <dbl>   <dbl>
#> 1 ATAC         tissue  cell_l…    97   198    0.668    0.852     -0.184  7.54e-5
#> 2 DNASE        primar… cell_l…    87   108    0.692    0.529      0.163  1.03e-2
#> 3 CHIP_HISTONE in_vit… cell_l…    34   115    0.533    0.655     -0.123  6.53e-2
#> 4 DNASE        tissue  cell_l…   108   108    0.459    0.529     -0.0698 6.33e-2
#> 5 PROCAP       primar… cell_l…     8    40    0.397    0.596     -0.199  1.28e-1
#> 6 CHIP_TF      tissue  cell_l…    21   286    0.569    0.468      0.101  1.13e-1
#> # ℹ 7 more variables: effect_size <dbl>, effect_type <chr>, method <chr>,
#> #   epsilon_sq <dbl>, kw_p <dbl>, p_adj <dbl>, significance_stars <chr>