bio-causal-genomics-mediation-analysis

name	bio-causal-genomics-mediation-analysis
description	Decompose genetic effects into direct and indirect paths through mediating variables using the mediation R package. Tests whether gene expression, methylation, or other molecular phenotypes mediate the effect of genetic variants on disease. Use when testing whether a molecular...
license	MIT
author	AIPOCH

Source: https://github.com/aipoch/medical-research-skills

Version Compatibility

Reference examples tested with: R stats (base), ggplot2 3.5+

Before using code patterns, verify installed versions match. If versions differ:

• R: packageVersion('<pkg>') then ?function_name to verify parameters

If code throws ImportError, AttributeError, or TypeError, introspect the installed package and adapt the example to match the actual API rather than retrying.

Mediation Analysis

"Test whether gene expression mediates the effect of this variant on disease" → Decompose the total genetic effect into direct and indirect (mediated) paths through a molecular phenotype, estimating ACME, ADE, and proportion mediated with bootstrap confidence intervals.

• R: mediation::mediate() for causal mediation analysis

Framework

Causal mediation decomposes the total effect of a treatment (genotype) on an outcome (phenotype) into:

• ACME (Average Causal Mediation Effect) - Indirect effect through the mediator
• ADE (Average Direct Effect) - Direct effect not through the mediator
• Total effect = ACME + ADE
• Proportion mediated = ACME / Total effect

Typical genomic applications:

• SNP -> gene expression (mediator) -> disease
• SNP -> DNA methylation (mediator) -> gene expression
• SNP -> protein levels (mediator) -> clinical outcome

Basic Mediation with the mediation Package

Goal: Decompose a genetic effect into direct and indirect (mediated) paths through a molecular phenotype.

Approach: Fit separate models for mediator and outcome, then run mediate() with bootstrap to estimate ACME (indirect), ADE (direct), and proportion mediated.

library(mediation)

# --- Step 1: Fit mediator model ---
# How does the treatment (genotype) affect the mediator (expression)?
mediator_model <- lm(expression ~ genotype + age + sex + pc1 + pc2, data = dat)

# --- Step 2: Fit outcome model ---
# How do treatment and mediator jointly affect the outcome?
# For binary outcome, use glm with family = binomial
outcome_model <- glm(
  disease ~ genotype + expression + age + sex + pc1 + pc2,
  data = dat, family = binomial
)

# --- Step 3: Run mediation analysis ---
# treat: name of treatment variable (genotype)
# mediator: name of mediator variable (expression)
# boot = TRUE: Use nonparametric bootstrap for CIs
# sims: Number of bootstrap simulations (1000 minimum for publication)
med_result <- mediate(
  mediator_model, outcome_model,
  treat = 'genotype', mediator = 'expression',
  boot   sims  


summarymed_result

Interpreting Results

# Extract key quantities
acme <- med_result$d0           # Indirect (mediated) effect
acme_ci <- med_result$d0.ci     # 95% CI for ACME
ade <- med_result$z0            # Direct effect
total <- med_result$tau.coef    # Total effect
prop_med <- med_result$n0       # Proportion mediated

cat('ACME (indirect):', round(acme, 4), '\n')
cat('ACME 95% CI:', round(acme_ci[1], 4), 'to', round(acme_ci[2], 4), '\n')
cat('ADE (direct):', roundade  
cat total  
cat prop_med  

eQTL Mediation

Goal: Test whether gene expression mediates the effect of an eQTL on a disease outcome across multiple genes.

Approach: Wrap the mediation workflow in a function, loop over candidate genes, and adjust p-values for multiple testing.

library(mediation)

run_eqtl_mediation <- function(dat, snp_col, expr_col, outcome_col, covariates) {
  covar_formula <- paste(covariates, collapse = ' + ')

  med_formula <- as.formula(paste(expr_col, '~', snp_col, '+', covar_formula))
  out_formula <- as.formula(paste(outcome_col, '~', snp_col, '+', expr_col, '+', covar_formula))

  med_model <- lm(med_formula, data = dat)

  if (length(unique(dat[[outcome_col]]   
    out_model  glmout_formula data  dat family  binomial
    
    out_model  lmout_formula data  dat
  

  result  mediate
    med_model out_model
    treat  snp_col mediator  expr_col
    boot   sims  
  

  data.frame
    snp  snp_col gene  expr_col
    acme  resultd0 acme_p  resultd0.p
    ade  resultz0 ade_p  resultz0.p
    total  resulttau.coef total_p  resulttau.p
    prop_mediated  resultn0
  



genes    
covars      

mediation_results  do.callrbind lapplygenes g 
  run_eqtl_mediationdat  g  covars



mediation_resultsacme_fdr  p.adjustmediation_resultsacme_p method  

Multi-Omics Mediation

Goal: Test cascading mediation chains across multiple molecular layers (e.g., SNP -> methylation -> expression -> disease).

Approach: Fit sequential models for each link in the chain and run separate mediation analyses for each mediator-outcome pair.

# Test mediation chains: SNP -> methylation -> expression -> disease
library(mediation)

# Step 1: SNP -> methylation
mod_meth <- lm(methylation ~ genotype + age + sex, data = dat)

# Step 2: methylation -> expression (controlling for genotype)
mod_expr <- lm(expression ~ methylation + genotype + age + sex, data = dat)

# Step 3: expression -> disease (controlling for methylation and genotype)
mod_disease <- glm(
  disease ~ expression + methylation + genotype + age + sex,
  data = dat, family = binomial
)

# Test methylation as mediator of SNP -> expression
med_meth_expr <- mediate(mod_meth, mod_expr, treat = 'genotype', mediator = 'methylation',
                         boot = TRUE sims  


med_expr_disease  mediatemod_expr mod_disease treat   mediator  
                            boot   sims  

High-Dimensional Mediation (HDMA)

Goal: Test thousands of potential mediators simultaneously (e.g., all CpG sites) to identify which mediate a genetic effect.

Approach: Use HIMA's penalized regression to jointly select significant mediators from a high-dimensional mediator matrix and estimate their indirect effects.

# For testing many potential mediators simultaneously (e.g., all CpG sites)
# install.packages('HIMA')
library(HIMA)

# X: treatment (genotype), M: high-dimensional mediators, Y: outcome
# HIMA uses penalized regression to select significant mediators

result <- hima(
  X = dat$genotype,
  Y = dat$disease,
  M = as.matrix(dat[, mediator_cols]),
  COV.XM = as.matrix(dat[, covariate_cols]),
  Y.family = 'binomial',
  M.family = 'gaussian',
  penalty = 'MCP'    # Minimax concave penalty (default)
)

# Results: significant mediators with estimated indirect effects
significant_mediators <- result[result$BH.FDR < 0.05, ]

Assumptions and Diagnostics

# --- Sequential ignorability assumption ---
# 1. No unmeasured confounders between treatment and mediator
# 2. No unmeasured confounders between mediator and outcome
# 3. No unmeasured confounders between treatment and outcome
# This assumption is UNTESTABLE but can be probed with sensitivity analysis

# --- Sensitivity analysis ---
# Tests how robust results are to unmeasured confounding
sens <- medsens(med_result, rho.by = 0.1, effect.type = 'indirect', sims = 1000)
summary(sens)

# rho: Correlation between residuals of mediator and outcome models
# At what rho does ACME cross zero? (larger |rho| = more robust)
# rho at which ACME = 0 is called the sensitivity parameter
# |rho| > 0.3: Reasonably robust to unmeasured confounding

plot(sens)

Visualization

library(ggplot2)

plot_mediation_diagram <- function(acme, ade, total, prop_med) {
  cat('Mediation Path Diagram:\n\n')
  cat('  Genotype ---[a]---> Mediator ---[b]---> Outcome\n')
  cat('      |                                     ^\n')
  cat('      +----------[c\' (ADE)]----------------+\n')
  cat('\n')
  cat('  Indirect (a*b = ACME):', round(acme, 4), '\n')
  cat('  Direct (c\' = ADE):', round(ade, 4), '\n')
  cat('  Total (c):', round(total, 4) 
  cat prop_med  


plot_mediation_results  results_df 
  results_dfgene  factorresults_dfgene levels  results_dfgeneorderresults_dfprop_mediated

  ggplotresults_df aesx  gene y  prop_mediated 
    geom_colfill   alpha   
    geom_hlineyintercept   linetype   color   alpha   
    coord_flip 
    labsx   y   title   
    theme_minimal

Related Skills

• mendelian-randomization - Causal inference using genetic instruments
• colocalization-analysis - Test if signals share a causal variant
• population-genetics/association-testing - GWAS for treatment-outcome associations
• multi-omics-integration/mofa-integration - Multi-omics data for mediation chains

Input Validation

This skill accepts requests that match the documented purpose of bio-causal-genomics-mediation-analysis and include enough context to complete the workflow safely.

Do not continue the workflow when the request is out of scope, missing a critical input, or would require unsupported assumptions. Instead respond:

bio-causal-genomics-mediation-analysis only handles its documented workflow. Please provide the missing required inputs or switch to a more suitable skill.