Functional annotation and pathway enrichment

Practical 10: Genomics, systems biology and multi-omics

Learning objectives

Preserve feature-to-gene mappings; define an eligible background; calculate a hypergeometric tail and verify its equivalent Fisher test; assess background sensitivity.

Biological question and data

Is a selected gene list enriched for a supplied pathway under an explicitly defined selection null?

The gene/pathway table is fictional. Every background gene is eligible, and the null treats a fixed-size selected set as uniformly drawn from that background. Selected genes are supplied for teaching, not discovered by a calibrated genome-wide analysis.

No biological download is required. Numerical inputs are hypothetical, not observed animal data. The runnable examples use base R. Count and abundance units are specified where used.

1. Define the eligible universe

Use the 24 assayed gene IDs and a prespecified six-gene list.

Show R code 
helper_dir <- if (dir.exists("R")) "R" else "../R"
stopifnot(file.exists(file.path(helper_dir,"teaching_data.R")))
source(file.path(helper_dir,"teaching_data.R"))
d <- teaching_data()
background <- rownames(d$counts)
selected <- sprintf("gene%02d",1:6)
stopifnot(!anyDuplicated(background),all(selected%in%background))
annotation <- d$annotation
stopifnot(!anyDuplicated(annotation$gene),setequal(annotation$gene,background))
table(annotation$pathway)

pathwayA pathwayB pathwayC 
       8        8        8

The measured universe is known. Counting annotations outside it would change the null population.

2. Calculate enrichment

For each set, compute overlap, null expectation and an upper-tail hypergeometric p-value. Confirm the first result with a 2-by-2 Fisher test.

Show R code 
one_set <- function(pathway) {
 member <- annotation$gene[annotation$pathway==pathway]
 N <- length(background); K <- length(member); n <- length(selected)
 x <- sum(selected%in%member)
 c(N=N,K=K,n=n,x=x,expected=n*K/N,p=phyper(x-1,K,N-K,n,lower.tail=FALSE))
}
result <- as.data.frame(t(vapply(unique(annotation$pathway),one_set,numeric(6))))
result$BH <- p.adjust(result$p,"BH")
result
          N K n x expected          p           BH
pathwayA 24 8 6 6        2 0.00020803 0.0006240899
pathwayB 24 8 6 0        2 1.00000000 1.0000000000
pathwayC 24 8 6 0        2 1.00000000 1.0000000000
Show R code 
z <- result[1,]
tab <- matrix(c(z$x,z$n-z$x,z$K-z$x,z$N-z$K-z$n+z$x),2,byrow=TRUE)
stopifnot(abs(fisher.test(tab,alternative="greater")$p.value-z$p)<1e-12)

PathwayA has six selected members in a background of eight pathway members. The enrichment model concerns gene-set representation, not measured pathway activity.

3. Change the background and explain the effect

Keep the pathway count fixed while adding unmeasured nonmembers to the universe. This is deliberately an inappropriate comparison for the same assay.

Show R code 
data.frame(background=c("eligible panel","inflated universe"),
 N=c(24,100),p=c(phyper(5,8,16,6,lower.tail=FALSE),
 phyper(5,8,92,6,lower.tail=FALSE)))
         background   N           p
1    eligible panel  24 2.08030e-04
2 inflated universe 100 2.34889e-08

Adding many nonmembers makes the same overlap appear more surprising. That is a change in assumptions, not more biological evidence.

4. Distinguish pathway size from evidence strength

The overlap count should be compared with its null expectation. A pathway containing many eligible genes can have a large overlap without unusual representation. Compare two fictional pathways using the same selected-list size.

Show R code 
pathway_check <- data.frame(pathway=c('small','large'),
 eligible_total=100,pathway_size=c(10,60),selected_total=10,overlap=c(4,7))
pathway_check$expected <- with(pathway_check,selected_total*pathway_size/eligible_total)
pathway_check$tail <- with(pathway_check,
 phyper(overlap-1,pathway_size,eligible_total-pathway_size,selected_total,lower.tail=FALSE))
pathway_check$BH <- p.adjust(pathway_check$tail,'BH')
pathway_check
  pathway eligible_total pathway_size selected_total overlap expected
1   small            100           10             10       4        1
2   large            100           60             10       7        6
         tail         BH
1 0.008224876 0.01644975
2 0.374237582 0.37423758

5. Report the tested universe

The eligible universe, mapping and selection rule define the test. Overlapping pathways can share most of their selected genes, so several discoveries need not mean several independent processes. Inspect the actual genes driving each overlap.

A binary selected list also discards effect direction and ranking. Ranked-set methods address different questions and assumptions; they should not be described as the same hypergeometric calculation. Neither type directly measures biochemical pathway activity.

Interpretation challenge

Can the smallest enrichment p-value identify the pathway with the largest physiological change?

No. It describes representation under the supplied selection null. Pathway size, background, mapping, overlap and biological measurements must be considered separately.

Interpretation and limits

Hypergeometric enrichment assumes an appropriate equal-selection null. Expression level, feature length, detectability, correlation and annotation overlap can require a more suitable method. Document mapping, set version, selection rule, background and multiplicity.

Exercises

  1. Why exclude unmeasured genes from this universe?
  2. Does a small p-value establish pathway activation?
  3. What if one measured gene maps to three proteins?
  1. They could not have entered the selected list under this assay.
  2. No; it is evidence about representation under the supplied null.
  3. Define the gene/protein analysis unit and avoid counting repeated gene evidence as independent entries.

Continue learning

B03 Databases and B14 Pathways and networks.

Continue with networks. The Multi-Omics pathway reuses this practical rather than creating a duplicate enrichment analysis. Return to the practical index or the unified map.