kegg-pathway-analysis

KEGG Pathway Enrichment Analysis Guide

Overview

KEGG (Kyoto Encyclopedia of Genes and Genomes) pathway enrichment analysis identifies biological pathways that are statistically over-represented among differentially expressed genes. This guide covers the two main enrichment approaches (ORA and GSEA), critical workflow decisions such as splitting genes by directionality, tool selection between R clusterProfiler and Python gseapy, and strategies for handling the notoriously unreliable KEGG REST API. It addresses recurring failure modes that produce incorrect pathway counts or stalled analyses.

The three most common errors in KEGG pathway analysis are: (1) combining up-regulated and down-regulated genes into a single enrichment run, which masks true pathway signals; (2) analysis failures caused by KEGG REST API timeouts with no fallback strategy; and (3) delaying result reporting while attempting cosmetic pathway name lookups that may never complete. This guide provides concrete solutions for each.

Key Concepts

ORA vs GSEA

Over-Representation Analysis (ORA) and Gene Set Enrichment Analysis (GSEA) are the two primary methods for pathway enrichment, and they differ in both input and statistical approach.

ORA takes a pre-filtered gene list (e.g., genes with padj < 0.05 and |log2FC| > 1.5) and tests whether KEGG pathway members are over-represented in that list relative to a background universe. ORA uses a hypergeometric test (Fisher's exact test). It is straightforward but discards magnitude information and depends heavily on the significance cutoff chosen.

GSEA takes a ranked list of all genes (typically ranked by log2 fold change or a signed significance statistic) without any cutoff. It computes a running enrichment score by walking down the ranked list and identifies pathways whose members cluster toward the top or bottom of the ranking. GSEA captures subtle coordinated changes that ORA may miss.

In practice, ORA via enrichKEGG() (clusterProfiler) or gp.enrichr() (gseapy) is the more common starting point. GSEA via gseKEGG() or gp.prerank() is preferred when you want to avoid arbitrary cutoffs or when effect sizes are small.

Directionality in Enrichment

When performing ORA, gene directionality -- whether a gene is up-regulated or down-regulated -- is critical. A single pathway can contain genes regulated in opposite directions. If up-regulated and down-regulated genes are combined into one list, their opposing signals cancel out, diluting the enrichment signal and masking genuinely enriched pathways. Running enrichment separately for up-regulated and down-regulated gene sets produces more accurate and interpretable results. This splitting is mandatory for ORA. GSEA inherently handles directionality through the signed ranking, though interpreting leading-edge genes by direction is still important.

KEGG Organism Codes

KEGG uses three-letter (or four-letter) organism codes to identify species-specific pathway databases. Using the wrong code silently returns empty results. Common codes:

| Organism | Code | |---|---| | Human | hsa | | Mouse | mmu | | Rat | rno | | Zebrafish | dre | | Drosophila | dme | | C. elegans | cel | | E. coli K-12 | eco | | P. aeruginosa PA14 | pau | | P. aeruginosa PAO1 | pae | | S. cerevisiae | sce | | A. thaliana | ath |

Gene ID format also varies by organism: eukaryotic species typically require Entrez gene IDs, while bacterial species use locus tags. Mismatched ID types are a silent failure mode.

KEGG API Reliability

The KEGG REST API (rest.kegg.jp) is rate-limited, frequently slow, and prone to timeouts. Both clusterProfiler::enrichKEGG() and direct HTTP requests to KEGG can fail unpredictably. Planning for API failures is not optional -- it is a necessary part of any KEGG-based workflow. Strategies include pre-fetching and caching pathway data, using offline gene set databases bundled with gseapy, and implementing retry logic with timeouts.

Decision Framework

Question: What enrichment analysis do you need?
|
+-- Have a pre-filtered DEG list (with cutoffs applied)?
|   +-- Yes --> ORA
|   |   +-- Using R? --> clusterProfiler::enrichKEGG()
|   |   +-- Using Python? --> gseapy.enrichr()
|   |   +-- KEGG API failing? --> gseapy with offline gene sets
|   +-- No, want cutoff-free analysis --> GSEA
|       +-- Using R? --> clusterProfiler::gseKEGG()
|       +-- Using Python? --> gseapy.prerank()
|
+-- Need to split by direction?
|   +-- ORA --> YES, always split up/down (mandatory)
|   +-- GSEA --> No split needed (direction encoded in ranking)
|
+-- KEGG API unreliable?
    +-- Try cached/pre-fetched data first
    +-- Fall back to gseapy offline databases
    +-- Use retry logic with short timeouts

| Scenario | Recommended Approach | Rationale | |----------|---------------------|-----------| | Standard ORA with R | clusterProfiler::enrichKEGG(), split by direction | Most widely used, integrates with Bioconductor ecosystem | | Standard ORA with Python | gseapy.enrichr() with KEGG_2021_Human | Offline gene sets avoid API dependency | | Cutoff-free enrichment | GSEA via gseKEGG() or gp.prerank() | Captures subtle coordinated changes, no arbitrary threshold | | KEGG API is down | Switch to gseapy offline databases | gseapy bundles KEGG gene sets locally | | Comparing conditions | Run separate up/down enrichment per condition | Enables direction-aware set operations across conditions | | Non-model organism | Verify organism code, use KEGGREST to check availability | Wrong code silently returns empty results |

Best Practices

1. Always split ORA by gene direction. Run enrichKEGG() or gp.enrichr() separately for up-regulated and down-regulated genes. Combining them inflates the gene list, dilutes enrichment signal, and produces incorrect pathway counts. Report the union of significant pathways from both directions.

1. Specify the background universe explicitly. Set the universe to all tested genes (the full set from your differential expression analysis), not just the significant ones. Omitting the universe defaults to all genes in the KEGG database, which inflates significance for well-studied pathways.

1. Pre-fetch and cache KEGG data before running enrichment. Download pathway-gene mappings at the start of the analysis and save them locally. This avoids mid-analysis failures when the KEGG API becomes unresponsive and makes the analysis reproducible.

1. Report the numeric answer before resolving pathway names. Once you have computed the count or list of significant pathway IDs, emit that result immediately. Resolving IDs to human-readable names via additional KEGG API calls is cosmetic and can timeout, losing the primary result.

1. Apply multiple testing correction consistently. Use adjusted p-values (p.adjust < 0.05, typically BH method) rather than raw p-values. Both clusterProfiler and gseapy apply correction by default, but always verify the cutoff is on the adjusted value.

1. Verify gene ID format matches the organism. Eukaryotic KEGG pathways expect Entrez gene IDs; bacterial species expect locus tags. A mismatch silently returns zero enriched pathways. Use bitr() in clusterProfiler or equivalent ID conversion if your input uses gene symbols.

1. Use gseapy as a fallback when clusterProfiler fails. When enrichKEGG() fails due to KEGG API issues, gseapy's enrichr() function with bundled offline gene sets (e.g., KEGG_2021_Human) provides equivalent ORA results without any network dependency.

Common Pitfalls

1. Combining up-regulated and down-regulated genes into a single enrichment run. Pathways with genes regulated in opposite directions cancel out, producing fewer significant pathways than the true count. Results from combined lists are unreliable. Example of the anti-pattern:

```r # WRONG: combining up and down genes into one list all_sig_genes <- rownames(subset(res, padj < 0.05