<p>If I read you correctly, you want to find the gene with the highest expression?</p>
<p>With RNA-seq data from multiple samples, I'd find this by ranking each gene by TPM/FPKM and taking the mean rank across samples. This method will tend to mitigate any large differences in a gene's expression between samples.</p>
<p>As you are comparing (RNA-seq) expression across genes you must take into account transcript length.</p>

CC BY-SA 4.0
2021-05-14T12:20:03.213Z
2021-05-14T12:20:03.213Z
Greg
<p>Transcript abundance quantification is a tricky topic since a read often could belong to several transcripts, so any &quot;count&quot; is a best guess as to which transcript it actually originates from. That being said, there are tools that can help you here:</p>
<ol>
<li><code>salmon</code> (as you mentioned) to quanitfy. Run it with <code>--numGibbsSamples 50</code> (or higher if your computer has the hardware. Also run with <code>-d</code> for the next tool</li>
<li><code>terminus</code> is a tool that also originates from the COMBINE-lab, It's relatively experimental, but it applies well to your use case. It will take the salmon output and assign counts to transcript groups. If a set of transcripts can be clearly quantified above a certain confidence, that transcript will stay unique. If it is confounded with other transcripts, those will be placed in a transcript group and their counts will be summarised. More info <a href="https://github.com/COMBINE-lab/terminus" rel="noreferrer">here</a></li>
</ol>
<p>You are right in that TPM, FPKM etc. are not appropriate measures to compare across samples. TMM, or DESeq2's median ratio method (what you get with <code>counts(dds, normalized=TRUE)</code> are a step in the right direction, but are still inappropriate across samples. I'd suggest this <a href="https://haroldpimentel.wordpress.com/2014/05/08/what-the-fpkm-a-review-rna-seq-expression-units/" rel="noreferrer">often posted blog post</a> in the topic.</p>
<p>There are no perfect solutions here since RNA-Seq is inherently a relative quanitifcation method so this will always be an issue. IMO, use <code>DESeq2::varianceStabilizingTransformation()</code>, <code>DESeq2::rlog()</code> or <code>limma::voom()</code> after terminus and rank your transcripts that way.</p>
<p>Just keep the inherent transcript incertainties in mind.</p>

CC BY-SA 4.0
2021-05-14T16:11:20.177Z
I did not know about `terminus`, it looks really interesting. Could you explain why TMM/median ratio method are still inappropriate across samples? I thought they were designed for between sample normalization.
Tomas Bencomo
CC BY-SA 4.0
2021-05-15T19:22:30.347Z
Don't get me wrong, I think both TMM and median-ratio are good enough for your approach, I just don't like to throw around the term normalization lightly. IMO normalization implies aligning dsitributions to allow for statistical testing, which those two don't _fully_ do. If there is a huge disparity in the number and magnitude of expressed genes between conditions, or lots of zero inflation, those methods do not produce a clean normalization either.
Bastian Schiffthaler
CC BY-SA 4.0
2021-05-14T07:37:31.873Z
2021-05-14T07:37:31.873Z
Bastian Schiffthaler
<p><code>vcf2maf</code> uses <code>VEP</code> to annotate variants, and I believe selects the default Ensembl transcript to use for annotation. Sometimes the transcript that <code>VEP</code> selects is not the transcript I'm interested in, usually because the selected transcript is not the most highly expressed transcript in my tissue of interest (skin). <code>vcf2maf</code> allows you to provide a transcript override list so that VEP annotates the variant using the specified transcripts instead.</p>
<p>I have several skin samples sequenced with bulk RNA-Seq. I want to estimate the average abundance for each transcript across all samples and then use these abundances to rank transcripts from most to least abundant. Then I will use the most abundant transcript as the default VEP transcript. I plan to use <code>salmon</code> or <code>kallisto</code> to quantify transcript abundance. Should I use TPM or normalized counts to calculate average expression?</p>
<p>My initial thought is to use normalized counts (generated by DESeq2 from raw counts). Are there any problems with this approach? GTEx displays transcript abundance with average TPM, but I thought TPM was inappropriate to use across samples because it doesn't account for between sample differences.</p>
<p>Update: I forgot to mention I also tried using TPM ranks like @ATpoint describes. I haven't fully compared how this compares to transcripts identified by normalized counts, but the initial genes I checked showed good concordance between methods</p>

<rna-seq><deseq2><transcriptome><normalization><salmon>
Calculating most abundant transcript from RNA-Seq data
# Calculating most abundant transcript from RNA-Seq data

(Asked by: username_0 on 2021-05-14)

`vcf2maf` uses `VEP` to annotate variants, and I believe selects the default Ensembl transcript to use for annotation. Sometimes the transcript that `VEP` selects is not the transcript I'm interested in, usually because the selected transcript is not the most highly expressed transcript in my tissue of interest (skin). `vcf2maf` allows you to provide a transcript override list so that VEP annotates the variant using the specified transcripts instead.

I have several skin samples sequenced with bulk RNA-Seq. I want to estimate the average abundance for each transcript across all samples and then use these abundances to rank transcripts from most to least abundant. Then I will use the most abundant transcript as the default VEP transcript. I plan to use `salmon` or `kallisto` to quantify transcript abundance. Should I use TPM or normalized counts to calculate average expression?

My initial thought is to use normalized counts (generated by DESeq2 from raw counts). Are there any problems with this approach? GTEx displays transcript abundance with average TPM, but I thought TPM was inappropriate to use across samples because it doesn't account for between sample differences.

Update: I forgot to mention I also tried using TPM ranks like @ATpoint describes. I haven't fully compared how this compares to transcripts identified by normalized counts, but the initial genes I checked showed good concordance between methods

**Tags:** rna-seq, deseq2, transcriptome, normalization, salmon

---

## 2 answers

### Answer by username_2 on 2021-05-14. Score: 5 | Accepted answer

Transcript abundance quantification is a tricky topic since a read often could belong to several transcripts, so any "count" is a best guess as to which transcript it actually originates from. That being said, there are tools that can help you here:

1. `salmon` (as you mentioned) to quanitfy. Run it with `--numGibbsSamples 50` (or higher if your computer has the hardware. Also run with `-d` for the next tool
2. `terminus` is a tool that also originates from the COMBINE-lab, It's relatively experimental, but it applies well to your use case. It will take the salmon output and assign counts to transcript groups. If a set of transcripts can be clearly quantified above a certain confidence, that transcript will stay unique. If it is confounded with other transcripts, those will be placed in a transcript group and their counts will be summarised. More info [here](https://github.com/COMBINE-lab/terminus)

You are right in that TPM, FPKM etc. are not appropriate measures to compare across samples. TMM, or DESeq2's median ratio method (what you get with `counts(dds, normalized=TRUE)` are a step in the right direction, but are still inappropriate across samples. I'd suggest this [often posted blog post](https://haroldpimentel.wordpress.com/2014/05/08/what-the-fpkm-a-review-rna-seq-expression-units/) in the topic.

There are no perfect solutions here since RNA-Seq is inherently a relative quanitifcation method so this will always be an issue. IMO, use `DESeq2::varianceStabilizingTransformation()`, `DESeq2::rlog()` or `limma::voom()` after terminus and rank your transcripts that way.

Just keep the inherent transcript incertainties in mind.

> **Comments:**
> **username_0**: I did not know about `terminus`, it looks really interesting. Could you explain why TMM/median ratio method are still inappropriate across samples? I thought they were designed for between sample normalization.  
> **username_2**: Don't get me wrong, I think both TMM and median-ratio are good enough for your approach, I just don't like to throw around the term normalization lightly. IMO normalization implies aligning dsitributions to allow for statistical testing, which those two don't _fully_ do. If there is a huge disparity in the number and magnitude of expressed genes between conditions, or lots of zero inflation, those methods do not produce a clean normalization either.  

---

### Answer by username_1 on 2021-05-14. Score: 0

If I read you correctly, you want to find the gene with the highest expression?

With RNA-seq data from multiple samples, I'd find this by ranking each gene by TPM/FPKM and taking the mean rank across samples. This method will tend to mitigate any large differences in a gene's expression between samples.

As you are comparing (RNA-seq) expression across genes you must take into account transcript length.


---

