测序公司给的分析报告中为什么没有Ensembl GTF注释版本？

导读：

在做RNA-seq等测序分析时，通常会用到Ensembl的GTF基因注释文件，它有不同的版本。我们调研了十几家测序公司的报告，令人惊讶的是，在报告中明确说明GTF版本的公司几乎没有，甚至个别公司连基因组版本都没有明确说明。本文以human为例，比较了不同版本GTF所包含的基因情况，以帮助大家更好地理解基因注释的变化。

在进行RNA-seq定量分析时，常用的定量软件，例如featureCounts，HTSeq等都会要求输入一个GTF基因注释文件。不同的数据库和研究机构（如NCBI、Ensembl、GENCODE等）使用不同的算法和标准进行基因注释，因此即使是同一基因组版本，不同机构提供的注释文件也可能存在差异。

随着科学技术的进步和研究的深入，对基因组的认识不断更新，新的基因、转录本和其他基因组特征的发现需要反映在注释文件中，因此会定期更新GTF注释文件以包含最新的生物学信息。Ensembl GTF约3-6个月更新一次，每次更新都会删除部分基因或者添加部分基因。因此，不同的GTF版本中包含的基因数不一样。

以人类hg38基因组的GTF注释为例，从2014年7月的v76到2024年8月的v113版，前前后后总共有38个版本，基因数（包括protein_coding，lncRNA等）从58051个（v85）到78932个（v113），极差为20881个基因。

在对数据进行重新分析时，使用的GTF版本不匹配，会造成部分基因的缺失。例如GEO数据库中很多RNA-seq上传的都是Ensembl基因的原始count矩阵，在将原始count转成FPKM或者TPM时，必需知道基因的长度，如果GTF版本不匹配，就会出现部分基因长度未知，无法计算该基因的FPKM或者TPM了。

例如，在使用rnanorm将count矩阵转成TPM时，会弹出如下警告信息：

UserWarning: X contains 120 genes that are not in GTF <_io.TextIOWrapper name='Homo_sapiens.GRCh38.104.chr.gtf' mode='r' encoding='UTF-8'>. This will result in NaN values for missing genes in the output

最终生成的结果中，有120个基因的tpm是空值。

120个基因影响不大，约占全部基因的2‰左右。然而，参考上述统计结果，最大的差别可以高达2w个基因，是不是瞬间觉得“非常可怕”？

此问题由来已久，我们这里给出了一个预测Ensembl GTF版本的解决方案：

1，下载所有版本的GTF文件

2，提取其包含的基因

3，将原始表达矩阵里边的基因与各个版本的GTF基因取交集

4，若提供的表达矩阵里边的基因是某个版本GTF基因的子集，则表达矩阵使用的GTF就有可能是该GTF版本

Python代码如下（相关文件加好友索取）：

因此，在出报告时必需明确GTF注释版本，最好在报告中直接给出分析时所使用的GTF文件。这样可以确保结果的可重复性，同时也为进一步信息挖掘提供方便。

小伙伴们赶紧看下你手上的报告，找一找基因组版本和注释版本吧！

微生信助力高分文章，用户220000+，谷歌学术4400+