Different results of spearman correlation between TPM and FPKM

Question

TPM and FPKM of RNA-Seq data form GDC TCGA calculated based STAR were retrieved, respectively. The correlation between a specific gene, e.g. HIF1A, and other genes were calculated based on TPM and FPKM, respectively. And the significant genes were filtered based on P value < 0.05. Here is the intersection of these two significant gene set shown in venn diagram. IMO, the two gene sets should not be obiviously different. Ang advice for this discrepance would be appreciated!

Answer 1

This shouldn't be surprising that you see different correlations between gene expression data when expressed in different units. To see why, let's look at how these units are defined.

Let's denote the absolute number of fragments we sequence from an RNA-seq measurement as $X_{gi}$ for gene $g$ in sample $i$, and $l_g$ as the length of gene $g$ in kilobases. Then our different expression values will be:

$$ f_{gi} = \frac{ X_{gi} \cdot 10^6 }{ l_g \cdot \sum_h X_{hi} } \text{ in FPKM} $$

or

$$ t_{gi} = \frac{ X_{gi} \cdot 10^6 }{ l_g \cdot \sum_h \frac{ X_{hi} }{ l_h } } \text{ in TPM} $$

You convert from FPKM to TPM units like so:

$$ t_{gi} = 10^6 \cdot \frac{ f_{gi} }{ \sum_h f_{hi} } $$

As you can see, within any single sample, $i$, this last denominator is fixed. But it is not fixed across samples (i.e. $\sum_h f_{hi} \ne \sum_h f_{hj}$ if $i \ne j$). This means that the ordering for a single gene's expression across samples is not preserved when you switch between units (i.e. if $f_{gi} > f_{gj}$, it is not guaranteed that $t_{gi} > t_{gj}$).

Because Spearman correlation is based on the order of each gene across samples in your dataset, and the order of genes in FPKM will be somewhat different from the order of genes in TPM, you will find different correlations between HIF1A and other genes when using the different units.

Here's a simplified, concrete example. Imagine that your transcriptome is made of up 3 genes (A, B, and C) with lengths 1 kb, 2 kb, and 3 kb, respectively. And imagine you're doing RNA-seq in 4 samples. If your absolute count matrix is this

$$ X = 10^6 \cdot \begin{bmatrix} 1 & 2 & 9 & 7 \\ 2 & 2 & 6 & 5 \\ 3 & 2 & 2 & 5 \end{bmatrix} $$

Then the expression matrix in FPKM is

$$ F = \begin{bmatrix} 166\,666.7 & 333\,333.3 & 529\,411.76 & 411\,764 \\ 166\,666.7 & 166\,666.7 & 176\,470.59 & 147\,058.82 \\ 166\,666.7 & 111\,111.1 & 39\,215.69 & 98\,039.22 \end{bmatrix} $$

and in TPM its

$$ T = \begin{bmatrix} 333\,333.3 & 545\,454.5 & 710\,526.32 & 626\,865.7 \\ 333\,333.3 & 272\,727.3 & 236\,842.11 & 223\,880.6 \\ 333\,333.3 & 181\,818.2 & 52\,631.58 & 149\,253.7 \end{bmatrix} $$

If you calculate the Spearman correlation for genes B and C against A, then you get these values.

Gene	Correlation with A (FPKM)	Correlation with A (TPM)
B	0.316	-0.8
C	-1	-1

You can see that gene C is anti-correlated with A, but B actually switches signs and goes from moderately correlated in FPKM to strongly anti-correlated in TPM.

So your observation that the same comparison on the different units gives pretty different results is expected. If you made a scatter plot with Spearman's $\rho$ from TPM along the x-axis and Spearman's $\rho$ from FPKM along the y-axis for all the genes that had a correlation p-value < 0.05, you'll probably have a lot of genes that also flip signs like I showed with this toy example.

This paper by Wagner et al. is probably a good paper to review. It can give some context about RPKM and how to think about count-based sequencing data. Correlation-based analyses on data like this are typically used as a first-pass approach to get an idea about how things are behaving. But it's not very precise, as you've experienced. If you're doing exploratory work, I'd be inclined to pick TPM units, look at the results you have, then try to come up with a differential analysis question that can be handled in a more statistically robust way.