Dr Cheng ZhangDigital Liver Group · translational systems medicine

Bioinformatics, step by step · Workshop

Download handouts + scripts Download the data

Transcriptomics workshop: from reads to pathways

A hands-on workshop that goes with the video series Bioinformatics, step by step (Foundations F1–F8 and Transcriptomics T1–T10). Each module has a short handout (this folder), a runnable R script (scripts/), and questions; the answers are in the answers (available from your instructor). Every module starts by pointing to the video that explains the theory.

Data. Public RNA-seq from mouse hearts after myocardial infarction (MI) or sham surgery, collected at day 1 and day 3, two replicates per group: PRJNA412001. Paper: HIF1 mediates a switch in pyruvate kinase isoforms after myocardial infarction (PubMed 29652636).

sample condition sample condition
SRR6068402 SHAM_1D SRR6068406 SHAM_3D
SRR6068403 SHAM_1D SRR6068407 SHAM_3D
SRR6068404 MI_1D SRR6068408 MI_3D
SRR6068405 MI_1D SRR6068409 MI_3D

Learning outcomes

By the end, participants can:

  1. read FASTQ, FASTA, GTF and SAM files by hand, and explain each field (quality scores, flags, CIGAR);
  2. explain what an RNA-seq count measures, and how kallisto turns reads into estimated counts;
  3. compute CPM, TPM and DESeq2 size factors by hand, and say when to use each;
  4. explore samples with PCA, and explain why a variance stabilising transformation is needed first;
  5. explain why counts are modelled with a negative binomial distribution, and what dispersion means;
  6. run DESeq2, and reproduce a Wald test and a Benjamini–Hochberg adjustment by hand;
  7. run over-representation and gene set enrichment analyses, and sanity check the results;
  8. write a design formula with more than one factor, and interpret its coefficients.

Modules

# module video(s) script time
0 Setup – 00_install.R before the day
F The files of an RNA-seq experiment F1–F8 00b_files_and_formats.R 45 min
1 Reads to counts T1, T2 01_consolidate_kallisto.R 45 min
2 Normalisation by hand T3, T5 02_normalisation_by_hand.R 45 min
3 Exploring samples with PCA T4 03_exploration_pca.R 40 min
4 What kind of noise? Count models T6 04_count_models.R 35 min
5 Differential expression T7 05_differential_expression.R 60 min
6 Multiple testing T8 06_multiple_testing.R 35 min
7 Functional analysis T9 07_functional_analysis.R 60 min
8 More than one factor T10 08_multifactor.R 40 min

Suggested one-day schedule (about 7 hours of contact time). Participants watch videos F1–F8 and T1–T3 before the day.

time
09:00 welcome; module F (files)
09:45 modules 1–2
11:15 break
11:30 modules 3–4 (watch videos T4 and T6 together first)
12:45 lunch
13:45 modules 5–6
15:20 break
15:35 modules 7–8
17:15 wrap-up: what would you do differently with your own data?

For a two-session version, split after module 4. For a short NGS-basics session, use videos F1–F8 and module F on their own. For a flipped version, assign the videos as homework and spend all contact time on the scripts and questions.

How to run

Open RStudio, set the working directory to this workshop folder (Session → Set Working Directory), then open the scripts in order and run them line by line. Results go to results/. Modules 6–8 read the DESeq2 output that module 5 writes, so run module 5 first.

setwd("path/to/transcriptomics tutorial/workshop")
source("scripts/00_install.R")    # once

Folder contents

workshop/
  README.md               this file
  00_setup.md, 00b_files_and_formats.md … 08_multifactor.md   handouts, one per module
  scripts/                one R script per module (tested with R 4.4.1, DESeq2 1.44)
  solutions/answers.md    answers, with the numbers you should get
  data/
    kallisto/             kallisto 0.44.0 output for all 8 samples (abundance.tsv, run_info.json)
    KallistoResults.zip   the same, zipped (as in the 2020 workshop)
    mart_export.txt       Ensembl BioMart: protein coding transcript -> gene
    count_gene.txt        gene-level estimated counts (21,814 genes x 8 samples)
    tpm_gene.txt          gene-level TPM
    metadata.txt          sample -> condition
    metadata_multi.txt    sample -> MI (0/1), day (1/3)
    KEGG.gmt, GOBP.gmt    gene set collections for module 7
    foundations/          first 10,000 read pairs of SRR6068402 (FASTQ), Pkm FASTA + GTF (Ensembl 99), for module F
  img/                    screenshots for module 1

What changed from the 2020 version

This workshop started from the PhD course material at sysmedicine/phd2020/transcriptomics. Main changes:

  • New modules: F (opening every file format by hand) and the theory modules 2, 4 and 6 that reproduce each calculation by hand: TPM, size factors, Poisson vs negative binomial, the Wald test and Benjamini–Hochberg.
  • PCA now compares raw TPM, log TPM and DESeq2's vst: on raw TPM, PC1 (87.6%) is dominated by a handful of very highly expressed genes.
  • Multifactor bug fixed. The 2020 script built model.matrix(~MI+day) but then still fitted design = ~conds with the old coldata. Module 8 now fits ~ day + MI, and adds an interaction model.
  • Functional analysis adds over-representation by hand (phyper), fgsea, and a sanity check on a suspicious result (olfactory transduction called "up" in heart tissue).
  • Question typo fixed: the 2020 question 4 asked about day 1 twice; it now asks about day 3.
  • Transcript-to-gene aggregation uses rowsum(), which is faster than merge() + aggregate() and gives the same table.