Skip to contents

Chemical Enrichment Analysis Using CTD

Source: R/enrichment_analysis.R
enrichment_CTD.Rd

Identifies chemicals whose known gene targets are significantly enriched in a user-supplied gene list or expression matrix, using data from the Comparative Toxicogenomics Database (CTD).

Four methods are available:

ORA

Over-Representation Analysis (default). Tests whether the overlap between your gene list and each chemical's target genes is larger than expected by chance. Uses enricher. Input: data frame with column EntrezID (character or numeric Entrez gene IDs) and an optional numeric value column.

GSEA

Gene Set Enrichment Analysis. Uses a ranked gene list (ranked by the numeric column in the input, e.g. p-values or fold changes) to detect chemicals whose targets cluster toward the top or bottom of the ranking. Uses fgsea.

CAMERA

Competitive gene-set test accounting for inter-gene correlation. Uses camera. Input: a numeric expression matrix (genes x samples) plus a design matrix and a contrast.

GSVA

Gene Set Variation Analysis (sample-level scoring). Returns a matrix of per-sample enrichment scores for each chemical, suitable for downstream clustering or association testing. Uses gsva. Input: a numeric expression matrix (genes x samples).

Usage

enrichment_CTD(
  x,
  method = c("ORA", "GSEA", "CAMERA", "GSVA"),
  design = NULL,
  contrast = NULL,
  id_type = NULL,
  pAdjustMethod = "BH",
  interaction_types = NULL,
  gene_id_type = c("symbol", "entrez"),
  ...
)

Arguments

x

The input. Its expected type depends on method:

  • For "ORA" and "GSEA": a data frame with at least two columns, EntrezID (character or numeric Entrez gene IDs) and a numeric value column (e.g. p-value). For GSEA, an optional column named stat can be added with a signed ranking statistic (e.g. the moderated t-statistic from limma::eBayes()); when present it is used directly for ranking, preserving directionality and avoiding ties. When absent, the second column is transformed via -log10() with a warning. The second column is ignored by ORA.

  • For "CAMERA" and "GSVA": a numeric expression matrix with genes in rows and samples in columns. rownames(x) must be either Entrez IDs or HGNC SYMBOLs.

method

Character. Enrichment method: "ORA" (default), "GSEA", "CAMERA", or "GSVA".

design

Design matrix (required when method = "CAMERA").

contrast

Contrast specification for camera (column number, column name, or numeric vector). Required when method = "CAMERA".

id_type

Either "entrez", "symbol", or NULL (default) for auto-detection from rownames(x). Only used when method is "CAMERA" or "GSVA".

pAdjustMethod

Character. Method for multiple testing correction: one of "BH" (Benjamini-Hochberg, default), "bonferroni", "fdr" (alias for BH), or "none". Not used for method = "GSVA" (which returns scores rather than p-values).

interaction_types

Character vector of CTD InteractionActions values to retain when building gene sets, or NULL (default) to use all cached interactions. Values follow the verb\^{}noun convention used by CTD, e.g. "increases\^{}expression", "decreases\^{}expression", "affects\^{}binding". A gene is included in a chemical's set if any of its recorded interaction actions matches one of the specified types. Requires that import_CTD() has been run (the filter is applied to the cached ctd_interactions.rda file). Restricting to expression interactions is recommended for RNA-seq analyses to improve biological specificity.

gene_id_type

Character. Identifier type used in the EnrichedGenes output column: "symbol" (default) returns HGNC gene symbols with Entrez ID as fallback for unmapped genes; "entrez" skips the symbol lookup and returns Entrez IDs directly. Only applies to "ORA" and "GSEA"; ignored by "CAMERA" and "GSVA".

...

Additional arguments forwarded to the underlying engine: enricher for ORA (e.g. universe, minGSSize, maxGSSize), fgseaMultilevel for GSEA (e.g. minSize, maxSize, nproc), camera for CAMERA, gsva for GSVA (e.g. minSize, maxSize).

Value

  • For "ORA", "GSEA", and "CAMERA": a data frame of enrichment results sorted by PValueAdjusted ascending. All three methods share the leading columns ChemicalID, ChemicalName, Method, PValue, PValueAdjusted; method-specific extras follow (see the package vignette for the full per-method schema).

  • For "GSVA": a numeric matrix of enrichment scores with chemicals (CTD chemical IDs) in rows and samples in columns.

Details

Before calling this function you must import the CTD data once with import_CTD. If the cached data is not found, the function stops with an informative error message.

Data Licensing Disclaimer

This package does not bundle or redistribute any CTD data. The Comparative Toxicogenomics Database is maintained by NC State University and its data are subject to specific licensing terms. Users must download the data directly from https://ctdbase.org and comply with the CTD Terms of Service (https://ctdbase.org/about/legal.jsp).

See also

import_CTD to import and cache the CTD data; plot_CTD to visualize results.

Examples

# Import the bundled sample data first:
sample_file <- system.file(
    "extdata", "CTD_chem_gene_ixns_sample.csv",
    package = "ctdR"
)
import_CTD(sample_file)
#> Reading CTD chemical-gene interactions from: /home/runner/work/_temp/Library/ctdR/extdata/CTD_chem_gene_ixns_sample.csv
#> Filtered to 86 human interactions
#> Mapping genes for 10 chemicals...
#> Warning: 10 ChemicalID(s) appear with more than one ChemicalName in the CTD file; only the first name per ID is retained. Affected IDs: D000082, D001564, D002104, D003907, D004958 ... (and 5 more)
#> CTD data cached successfully in: ~/.cache/ctdR
#>   10 chemicals | 17 unique genes | 0 s

# ORA / GSEA: prepare a gene list with Entrez IDs and a numeric value
genes <- data.frame(
    EntrezID = c("7124", "3569", "7157", "672", "1956"),
    pvalue = c(0.001, 0.003, 0.01, 0.02, 0.05)
)
ora_results <- enrichment_CTD(genes, method = "ORA")
#> Warning: column name ‘FoldEnrichment’ is duplicated in the result
gsea_results <- enrichment_CTD(genes, method = "GSEA")
#> Warning: GSEA: no 'stat' column found in input. Falling back to -log10(second column) for ranking, which loses directionality and may produce ties at non-significant genes. Add a signed ranking statistic (e.g. the moderated t-statistic from limma::eBayes()) as a column named 'stat' to suppress this warning.
#> Warning: All values in the stats vector are greater than zero and scoreType is "std", maybe you should switch to scoreType = "pos".

# CAMERA / GSVA: expression matrix + design + contrast.
# Uses the bundled GSE311566 subset
# (Dex vs DMSO, female PBMCs; see inst/extdata/README.md).
gse <- readRDS(system.file(
    "extdata", "GSE311566_subset.rds", package = "ctdR"
))
expr <- gse$expr
grp  <- gse$coldata$group
d    <- model.matrix(~ grp)
camera_results <- enrichment_CTD(expr, method = "CAMERA",
    design = d, contrast = 2)

# GSVA: per-sample enrichment scores
gsva_scores <- enrichment_CTD(expr, method = "GSVA")
#> ℹ GSVA version 2.6.2
#> ℹ Searching for rows with constant values
#> ℹ Calculating GSVA ranks
#> ℹ kcdf='auto' (default)
#> ℹ GSVA dense (classical) algorithm
#> ℹ Row-wise ECDF estimation with Gaussian kernels
#> ℹ Calculating row ECDFs
#> ℹ Calculating column ranks
#> ℹ GSVA dense (classical) algorithm
#> ℹ Calculating GSVA scores for 10 gene sets
#> ✔ Calculations finished