Pathway Analysis Tutorial

MSI-RI

November 11, 2025

​

What is Pathway Analysis?

Analysis technique to summarize biological pathways/processes/cellular compartments/molecular functions, etc that are changing in an experiment

→ Exploratory!

What is Pathway Analysis?

Typically performed after RNA-seq differential gene expression (DGE) analysis

bmm Overexpression VS Wild Type (baseline)

What is Pathway Analysis?

Typically performed after RNA-seq differential gene expression (DGE) analysis

bmm Overexpression VS Wild Type (baseline)

Nice! *

Significance threshold of FDR < 0.05 yields 112 DEGs!

What is Pathway Analysis?

Typically performed after RNA-seq differential gene expression (DGE) analysis

bmm Overexpression VS Wild Type (baseline)

…

What biological processes are these significant genes related to?

Significance threshold of FDR < 0.05 yields 112 DEGs!

What is Pathway Analysis?

Typically performed after RNA-seq differential gene expression (DGE) analysis

bmm Overexpression VS Wild Type (baseline)

…

What biological processes are these significant genes related to?

I definitely don’t want to look up each gene individually……

Why do we use it?

1. Helps us derive meaning from results of high-throughput data, which can often be too difficult and overwhelming to interpret gene-by-gene

​

Why do we use it?

1. Helps us derive meaning from results of high-throughput data, which can often be too difficult and overwhelming to interpret gene-by-gene

​

1. Helps us validate/confirm expected results

If I knock out or overexpress a gene involved in regulating metabolism, I would expect that gene to be differentially expressed, but I may also expect that all the genes that are regulated by it in a metabolic pathway to also be perturbed by its change in expression.

​

​

​

Why do we use it?

1. Helps us derive meaning from results of high-throughput data, which can often be too difficult and overwhelming to interpret gene-by-gene

​

1. Helps us validate/confirm expected results

If I knock out or overexpress a gene involved in regulating metabolism, I would expect that gene to be differentially expressed, but I may also expect that all the genes that are regulated by it in a metabolic pathway to also be perturbed by its change in expression.

​

1. Helps generate new hypotheses and design subsequent experiments

Where do we find these biological pathways?

• Commonly-used Databases:
• Gene Ontology (GO)� http://geneontology.org/
• Kyoto Encyclopedia of Genes and Genomes (KEGG)� https://www.genome.jp/kegg/
• REACTOME� https://reactome.org/
• Molecular Signatures Database (MSigDB)� https://www.gsea-msigdb.org/gsea/msigdb
• There are many more!

Examples of pathway terms from databases

• Cell cycle
• Striated muscle contraction
• Androgen receptor signaling
• Huntington’s disease
• Kinase activity
• etc….

​

Terms may overlap across databases or may be unique to a database. If there is a particular pathway you are interested in, search the database websites to see if they have it and what genes are contained in it!

Examples of pathway terms from databases

​

• Cell cycle
• Striated muscle contraction
• Androgen receptor signaling
• Huntington’s disease
• Kinase activity
• etc….

​

Terms may overlap across databases or may be unique to a database. If there is a particular pathway you are interested in, search the database websites to see if they have it and what genes are contained in it!

Each pathway is defined by a set of genes ("gene set")

Databases: Gene Ontology (GO)

• Hierarchical, relational database of terms organized into 3 ontologies

​

Ontology: a set of concepts and categories in a subject area or domain that shows their properties and the relations between them

Databases: Gene Ontology (GO)

• Hierarchical, relational database of terms organized into 3 ontologies
• Each term has a GO term ID (e.g. GO:0019827 = “Stem cell maintenance”)

​

Databases: Gene Ontology (GO)

• Hierarchical, relational database of terms organized into 3 ontologies
• Each term has a GO term ID (e.g. GO:0019827 = “Stem cell maintenance”)
• Genes are annotated with GO terms based on different sources of evidence:
• experimental evidence
• phylogenetic evidence
• computational evidence
• author statements
• curatorial statements
• automatically generated annotations

See: http://geneontology.org/docs/guide-go-evidence-codes/

GO graphs show relationships of terms

Kogelman et al., 2011; https://doi.org/10.1186/1471-2164-12-607

GO graphs show relationships of terms

3 ontologies:

MF, BP, and CC

Kogelman et al., 2011; https://doi.org/10.1186/1471-2164-12-607

GO graphs show relationships of terms

Parent term (General)

Child term (Specific/

Specialized)

Thousands of genes annotated with term

Tens to hundreds of genes annotated with term

Kogelman et al., 2011; https://doi.org/10.1186/1471-2164-12-607

GO graphs show relationships of terms

No over-representation

Kogelman et al., 2011; https://doi.org/10.1186/1471-2164-12-607

Parent term (General)

Child term (Specific/

Specialized)

Thousands of genes annotated with term

Tens to hundreds of genes annotated with term

GO graphs show relationships of terms

GO graphs show relationships of terms

Significantly over-represented in our test set

Kogelman et al., 2011; https://doi.org/10.1186/1471-2164-12-607

Parent term (General)

Child term (Specific/

Specialized)

Tens to hundreds of genes annotated with term

No over-representation

Thousands of genes annotated with term

GO graphs show relationships of terms

GO graphs show relationships of terms

Significantly over-represented in our test set

Kogelman et al., 2011; https://doi.org/10.1186/1471-2164-12-607

Parent term (General)

Child term (Specific/

Specialized)

Thousands of genes annotated with term

Tens to hundreds of genes annotated with term

No over-representation

GO graphs show relationships of terms

No over-representation

Kogelman et al., 2011; https://doi.org/10.1186/1471-2164-12-607

Parent term (General)

Child term (Specific/

Specialized)

Thousands of genes annotated with term

Tens to hundreds of genes annotated with term

GO graphs show relationships of terms

No over-representation

Significantly over-represented in our test set

Side note: Plots such as this can be created using the goplot() function in the clusterProfiler R package.

Kogelman et al., 2011; https://doi.org/10.1186/1471-2164-12-607

Databases: KEGG

https://www.genome.jp/kegg/

KEGG pathview

This is an example of the “Cell Cycle” pathway in KEGG. This was produced using KEGG pathview (pathview package in R)

KEGG pathview

This is an example of the “Cell Cycle” pathway in KEGG. This was produced using KEGG pathview (pathview package in R)

​

What should we look at?

• The number of genes perturbed in the pathway

​

KEGG pathview

This is an example of the “Cell Cycle” pathway in KEGG. This was produced using KEGG pathview (pathview package in R)

​

What should we look at?

• The number of genes perturbed in the pathway
• What directions the genes are changing

red = upregulation

green = downregulation

​

KEGG pathview

This is an example of the “Cell Cycle” pathway in KEGG. This was produced using KEGG pathview (pathview package in R)

​

What should we look at?

• The number of genes perturbed in the pathway
• What directions the genes are changing

red = upregulation

green = downregulation

• What are those genes roles in the pathway? Are they activators or repressors of the pathway?
• For example, downregulation of a repressor could mean activation of the pathway

Databases: REACTOME

Covers “classical intermediary metabolism, signaling, transcriptional regulation, apoptosis and disease”

​

REACTOME may* have more disease-related terms and pathways than KEGG and GO.

https://reactome.org/

*: There is a lot of overlap between these databases! Be sure to search through the databases to find the pathways you care about analyzing!

Databases: Molecular Signatures

H - hallmark

C1 - positional gene sets

C2 - curated gene sets

• includes KEGG and REACTOME!

C3 - regulatory target gene sets

C4 - computational gene sets

C5 - ontology gene sets

• includes GO!

C6 - oncogenic signature sets

C7 - immunologic signature sets

C8 - cell type signature gene sets

https://www.gsea-msigdb.org/gsea/msigdb

Databases: Molecular Signatures

H - hallmark

C1 - positional gene sets

C2 - curated gene sets

• includes KEGG and REACTOME!

C3 - regulatory target gene sets

C4 - computational gene sets

C5 - ontology gene sets

• includes GO!

C6 - oncogenic signature sets

C7 - immunologic signature sets

C8 - cell type signature gene sets

https://www.gsea-msigdb.org/gsea/msigdb/human/genesets.jsp?collection=CGP

Choose a database that fits your needs!

Pathway Analysis Techniques

​

• Over-Representation Analysis (ORA)

“Are more genes associated with gene set X than would be expected (over-represented) in my significant genes?”

​

• Gene Set Enrichment Analysis (GSEA)

“Based on the expression responses of all of the genes in my dataset, is a particular set of genes overexpressed (or under) in one treatment group or the other?”

Over-Representation Analysis (ORA)

Over-representation Analysis (ORA)

“Are more genes associated with pathway term X than would be expected (over-represented) in my significant gene set?”

• Input required is a list of genes (i.e. “significant”, “test set”,“query set”):
• Significantly differentially expressed genes
• Genes in a network module
• Other “interesting” set of genes
• You gene list size matters:
• Ideally, ~150 to ~600 genes (for whole transcriptome)
• Too many or too few lead to unreliable results
• Not a hard rule!

Recall: GO Term Enrichment

No over-representation

Significantly over-represented in our test set

Kogelman et al., 2011; https://doi.org/10.1186/1471-2164-12-607

Side note: Plots such as this can be created using the goplot() function in the clusterProfiler R package.

Statistics of ORA

List of 94 DEGs

ORA

results

• GeneRatio: Of the 94 of the genes in our test/query set, 10 are annotated with “response to bacterium”

​

• BgRatio: Of the 12,175 genes in our reference set (not including the query gene set), 270 are annotated with “response to bacterium”

​

List of all 12445 genes tested

Your data

gene sets to test

Statistics of ORA

GO: 0009617 (5.7e-12)

response to bacterium

10/94 | 270/12175

• 10/94 of the genes in our test set are annotated with “response to bacterium”

​

• 19/12175 of the reference set (does not include query gene set) are annotated with “response to bacterium”

​

GO term ID (p-value significance)

GO term description

Fraction of test genes annotated with term

Fraction of reference genes annotated with term

Statistics of ORA: confusion matrix

significant genes:

genes examined:

Term

No Term

10

84

270

10 + 84 = 94 sig. genes

270 + 12175 = 12445 "universe"

12175

Fisher’s Exact Test or Hypergeometric test

GO: 0009617 (5.7e-12)

response to bacterium

10/94 | 270/12175

sig. genes with term

sig. genes without term

examined genes with term

examined genes without term

Statistics of ORA: confusion matrix

Term

No Term

10

84

270

10 / 94 = 0.106

270 / 12445 = 0.022

12175

Fisher’s Exact Test or Hypergeometric test

GO: 0009617 (5.7e-12)

response to bacterium

10/94 | 270/12175

sig. genes with term

sig. genes without term

examined genes with term

examined genes without term

is

this greater than this?

significant genes:

genes examined:

Statistics of ORA

• 8/69 of the genes we in our test set are annotated with “threonine-type endopeptidase activity”

​

• 19/11138 of the reference set (does not include query gene set) are annotated with “threonine-type endopeptidase activity”

​

GO term ID (p-value significance)

GO term description

Fraction of test genes annotated with term

Fraction of reference genes annotated with term

GO: 0004298 (5.7e-12)

threonine-type endopeptidase activity

8/69 | 19/11138

Statistics of ORA

test set

reference

Term

No Term

8

61

19

69

11,138

11,119

11,207 total genes in universe

GO: 0004298 (5.7e-12)

threonine-type endopeptidase activity

8/69 | 19/11138

Fisher’s Exact Test or Hypergeometric test

Statistics of ORA

Is this term over-represented in our test set?

test set

reference

Term

No Term

8

61

19

69

11,138

11,119

11,207 total genes in universe

Fisher’s Exact Test or Hypergeometric test

GO: 0004298 (5.7e-12)

threonine-type endopeptidase activity

8/69 | 19/11138

Reading ORA Results Table

Visualizating results: dotplot

• Shows the top X significant terms (configurable)
• Dot size indicates the number of genes of the test set that are annotated with the term
• GeneRatio is as defined in the table (see previous slides) the total number of genes in the test set with term, divided by test set size
• Dot color indicates significance based on adjusted p-value

bmm Overexpression VS Wild Type (baseline)

Visualizating results: cnetplot

cnetplots are network plots that show the relationships of terms with each other and the significant genes that belong to those terms.

​

​

Visualizating results: cnetplot

cnetplots are network plots that show the relationships of terms with each other and the significant genes that belong to those terms.

​

• “Hub” nodes (grey) are the pathway terms. Node size indicates the number of genes enriched in the term
• Other nodes (in color) are genes. Node color indicates the fold change of the genes.

​

​

​

“Hub” node

Other nodes

Visualizating results: emapplot

Comparing/Contrasting Types of ORA Plots

• dotplot shows:
• Only the significantly enriched terms
• How much of your test set overlaps with various terms
• cnetplot shows:
• Relationship of terms to genes in your test set
• Directionality of gene expression response for each term
• Genes that are shared across terms
• emapplot shows:
• Relationship between significantly enriched terms
• “Biothemes” that are represented in your test set (e.g., response to various stimuli)

These three plot types are from the clusterProfiler package, specifically. Other pathway analysis tools may have similar or complementary visualizations.

Other ways to use ORA

• On a gene module from Weighted Gene Coexpression Network Analysis (WGCNA)
• Use on top or bottom 200-300 loadings from a principal component axis to determine pathways highly variable genes are related to
• Really, any way you could make an “interesting list of genes” from your data

Gene Set Enrichment Analysis (GSEA)

Gene Set Enrichment Analysis (GSEA)

1. rank all genes by evidence of change in expression
2. ask whether a gene set's genes are nearer the top, or bottom, of the ranked list

Gene Set Enrichment Analysis (GSEA)

• rank all genes by evidence of change in expression
• ask whether a gene set's genes are nearer the top, or bottom, of the ranked list

• left is heatmap of expression, with rows sorted:
• greatest evidence of positive change (A-B) are near top
• greatest negative evidence near bottom
• scant evidence in middle
• ​

Gene Set Enrichment Analysis (GSEA)

• rank all genes by evidence of change in expression
• ask whether a gene set's genes are nearer the top, or bottom, of the ranked list

• left is heatmap of expression, with rows sorted:
• greatest evidence of positive change (A-B) are near top
• greatest negative evidence near bottom
• scant evidence in middle
• Where do we find genes in a gene set?

e.g. Gene set 1 (cell cycle): P53, CycA, CDK1, SCF, SKP2,...

Gene Set Enrichment Analysis (GSEA)

“Based on the expression responses of all of the genes in my dataset, is a particular set of genes overexpressed (or under) in one treatment group or the other?”

Gene set 1 (cell cycle): P53, CycA, CDK1, SCF, SKP2,...

Gene Set Enrichment Analysis (GSEA)

How to rank?

• -log(p-value) == statistical evidence only

​

• -log(p-value) x direction of fold change

== positive and negative evidence are separated (this graph uses this rank)

​

• logFC == measurement of change without statistical evidence

Subramanian et al., 2005; https://doi.org/10.1073/pnas.0506580102

rank

rank value

Gene Set Enrichment Analysis (GSEA)

Ranked on one of these:

• -log(p-value)
• -log(p-value) x direction of FC

Subramanian et al., 2005; https://doi.org/10.1073/pnas.0506580102

Gene Set Enrichment Analysis (GSEA)

Ranked on one of these:

• -log(p-value)
• -log(p-value) x direction of FC

Note: Genes can be tied of rank when they have equal values. Tied genes will be randomly ordered during the statistical test, impacting the GSEA “random walk”. Rerunning the ranking can change the result, so it is good to export the ranked gene list.

Subramanian et al., 2005; https://doi.org/10.1073/pnas.0506580102

Gene Set Enrichment Analysis (GSEA)

Response to hypoxia

Subramanian et al., 2005; https://doi.org/10.1073/pnas.0506580102

each horizontal line is the position of a "response to hypoxia" gene

Gene Set Enrichment Analysis (GSEA)

Response to hypoxia

Subramanian et al., 2005; https://doi.org/10.1073/pnas.0506580102

-log10(p) * sign(logFC)

Ranked Gene List

Gene Set positions

Gene Set Enrichment Analysis (GSEA)

Response to hypoxia

Subramanian et al., 2005; https://doi.org/10.1073/pnas.0506580102

-log10(p) * sign(logFC)

Ranked Gene List

Gene Set positions

Score (calculated left to right)

Gene Set Enrichment Analysis (GSEA)

Subramanian et al., 2005; https://doi.org/10.1073/pnas.0506580102

-log10(p) * sign(logFC)

Ranked Gene List

Gene Set positions

Score (calculated left to right)

Enrichment Score (ES)

• Equal to the score at the peak of the running sum
• If negative, peak is where the rank metric is negative (e.g. down-regulation)
• used in permutation test

​

​

Leading Edge

• The genes that were "in front" of the position where the ES occurred

​

leading edge genes

Gene Set Enrichment Analysis (GSEA)

Response to hypoxia

Subramanian et al., 2005; https://doi.org/10.1073/pnas.0506580102

Gene Set Enrichment Analysis (GSEA)

Response to hypoxia

Subramanian et al., 2005; https://doi.org/10.1073/pnas.0506580102

GSEA

Enrichment Score (ES)

• Equal to the score at the peak
• If negative, peak is where the rank metric is negative (e.g. down-regulation)
• not yet ready for statistical test

​

Subramanian et al., 2005; https://doi.org/10.1073/pnas.0506580102

Leading Edge

• The genes that were "in front" of the position where the ES occurred

​

Gene Set Enrichment Analysis (GSEA)

Response to hypoxia

Subramanian et al., 2005; https://doi.org/10.1073/pnas.0506580102

e.g. -log10(p) * sign(logFC)

GSEA

Enrichment Score (ES)

• Equal to the score at the peak

​

Subramanian et al., 2005; https://doi.org/10.1073/pnas.0506580102

GSEA

Subramanian et al., 2005; https://doi.org/10.1073/pnas.0506580102

Enrichment Score (ES)

• Equal to the score at the peak
• Calculated using running-sum statistic that increases or decreases with each hit or miss to a gene set. Magnitude of increase/decrease is correlated with phenotype (i.e. metric used for ranking).

GSEA

Enrichment Score (ES)

• Equal to the score at the peak
• Calculated using running-sum statistic that increases or decreases with each hit or miss to a gene set. Magnitude of increase/decrease is correlated with phenotype (i.e. metric used for ranking).
• Degree to which a gene set is overrepresented at the top or bottom of a ranked list

Subramanian et al., 2005; https://doi.org/10.1073/pnas.0506580102

GSEA statistical test: permutation test

https://www.gsea-msigdb.org/gsea/doc/GSEAUserGuideFrame.html

ACTUAL ranked gene list

Same data, but permuted gene names

compute ES score on permuted data,

re-permute many times

label swap

permutation

Calculating significance

Kolmogorov Smirnov permutation test

​

• significance is proportional to rank of actual ES / number of permutations

​

​

https://www.pathwaycommons.org/guide/primers/data_analysis/gsea/

actual ES

mean of ES permutations

Calculating significance

Kolmogorov Smirnov test

​

• “What is the probability that this collection of samples could have been drawn from a specific probability distribution?”

​

​

https://www.pathwaycommons.org/guide/primers/data_analysis/gsea/

GSEA

Normalized Enrichment Score (NES)

​

• Enrichment score is normalized by mean ES of all permutations (random ordering of genes)

​

https://www.gsea-msigdb.org/gsea/doc/GSEAUserGuideFrame.html

NES =

actual ES

mean (ESs against all permutations of the dataset)

Calculating significance

Kolmogorov Smirnov test

​

• “What is the probability that this collection of samples could have been drawn from a specific probability distribution?”

​

• Distribution is made by permuting genes in the ranked list (or by permuting samples)

https://www.pathwaycommons.org/guide/primers/data_analysis/gsea/

GSEA

Normalized Enrichment Score (NES)

​

• Enrichment score is normalized by mean ES of all permutations (random ordering of genes)
• Primary statistic for examining gene set enrichment results

​

NES =

actual ES

mean (ESs against all permutations of the dataset)

https://www.gsea-msigdb.org/gsea/doc/GSEAUserGuideFrame.html

GSEA

Normalized Enrichment Score (NES)

​

• Enrichment score is normalized by mean ES of all permutations (random ordering of genes)
• Accounts for gene set size, correlations between genes and gene set. Can be compared across gene sets
• Primary statistic for examining gene set enrichment results

​

NES =

actual ES

mean (ESs against all permutations of the dataset)

https://www.gsea-msigdb.org/gsea/doc/GSEAUserGuideFrame.html

Examples of GSEA plots

Examples of GSEA plots

Test vs control

High in test

Low in test

Examples of GSEA plots

Leading edge

Leading edge

Test vs control

High in test

Low in test

GSEA Results Table

Table of significant pathway terms

Example table from clusterProfiler wrapper

GSEA Results Table

# genes associated with pathway term

GSEA Results Table

Enrichment score (ES)

GSEA Results Table

Normalized enrichment score (NES)

​

(influenced by the permutation method, the number of permutations, and the size of the expression set)

GSEA Results Table

Nominal p-value

​

(estimates the statistical significance of the enrichment score for a single gene set)

GSEA Results Table

Family Wise Error Rate (FWER)

​

(a more conservatively estimated probability that the normalized enrichment score represents a false positive finding)

GSEA Results Table

False Discovery Rate (FDR)

​

(estimated probability that a gene set with a given NES represents a false positive finding)

GSEA Results Table

Rank at Max

​

(position in the ranked list at which the maximum enrichment score occurred)

GSEA Results Table

Leading edge subset statistics

​

(% tags = % gene set hits in leading edge (coverage of GS)

% list = % ranked genes in leading edge (coverage of ranked list)

% signal = )

GSEA Results Table

Leading edge subset statistics

​

% tags = % gene set hits in leading edge (coverage of GS)

% list = % ranked genes in leading edge (coverage of ranked list)

% signal: complicated, but

Directly proportional to % tags

Inversely proportional to % list

Inversely proportional to gene set size

Leading edge

“% Signal” Example

“Tag %” is the percentage of gene set hits in leading edge (coverage of GS), “Gene %” is the percentage of ranked genes in leading edge (called “list” in the GSEA results table), N is the size of your ranked list, and N_h is the size of the gene set being tested.

​

For “REACTOME_DNA_REPAIR” (tag=39%,list=20%):

When to use ORA vs. GSEA:

When to use GSEA:

• If you have no differentially genes, GSEA might be able to detect the subtle shifts in expression in many genes of a pathway
• ​

Factors to consider when running GSEA

• What databases (gene sets) to run against
• What metric to use to rank the gene list (i.e. log2FC, p-value, etc)
• Type of permutations performed (samples vs genes)
• Note: MSI script only permutes genes. Desktop GSEA can permute both
• # permutations to perform (> 10000)
• Size of the expression set
• Gene set size (min=25, max=500, for 22k-25k tested genes)

​

GSEA

• Use the full gene differential expression results tables from edgeR or DESeq2
• Choose to rank genes on:
• log2FC
• -log(p-value) * direction of the fold change (preferred!)

​

​

​

​

​

​

GSEA paper: https://www.pnas.org/doi/10.1073/pnas.0506580102

-6.41

6.16

-5.49

5.49

…

…

…

​

-log(FDR) * direction of FC

Next: A Hands-on Guide (approx 50 minutes)

​

The later section of this tutorial will show you how to use a MSI-developed script for pathway enrichment analysis (ORA and GSEA) and interpreting the results.

This script is run from the command line.

A command line refresher module will be included