Bioinformatics, step by step · Workshop
Module 2 · Normalisation by hand
Videos: T3 Counts, lengths and TPM and T5 Size factors. Script: scripts/02_normalisation_by_hand.R
Goal. Compute TPM, CPM and DESeq2 size factors yourself, and check that you get the same numbers as the tools.
A. TPM from kallisto's numbers
TPM = (count ÷ length in kb) ÷ sum over all transcripts × 10⁶. kallisto uses the effective length, which is the number of positions where a fragment could start. For short transcripts, this is much smaller than the length.
Q2.1 Compute TPM for SRR6068402 from
est_countsandeff_length. Does it match kallisto'stpmcolumn?Q2.2 Repeat with
lengthinstead ofeff_length. Which transcripts change the most?
B. Why don't the gene-level TPM columns sum to one million?
Q2.3 Compute
colSums()ofdata/tpm_gene.txt. Why are the sums below 10⁶? What share of TPM did we drop?
C. CPM
Q2.4 Compute CPM for every gene and sample. CPM corrects for sequencing depth. What does it not correct for, compared with TPM? When does that not matter? (Hint: are you comparing genes, or samples?)
D. Size factors with the median of ratios
The my_size_factors() function in the script is the whole algorithm, in four lines:
- keep only genes with no zero counts;
- compute the log geometric mean of each gene: its pseudo-reference;
- for each sample, compute log(count) − log(pseudo-reference) for every gene;
- the size factor is exp(median of those values).
Q2.5 Do your size factors match
sizeFactors(estimateSizeFactors(dds))?Q2.6 Reproduce the worked example from video T5 (Actb, Gapdh, Myh6, Nppa and Postn; samples 02 and 04).
Q2.7 Which sample has the biggest gap between its size factor and its library-size ratio? Look at the share of reads that goes to Myh6 in each sample. Explain the gap in one sentence.
Q2.8 (theory) Why the median, and why the geometric mean? What would go wrong with the ordinary mean of the ratios?