收藏！看完联川这1万多字的问题解答，你就能从转录组小白变成大神

你是否困扰过，拿到转录组的结题报告和大量的表格无从下手？又会疑惑，拿到转录组的数据该怎么往下做验证？今年我们曾经出过一个专题，5千字带你了解2021年联川全转录组生信内容大升级！广受好评，今天我们会再次升级！

转录调控领域深耕16年的联川生物，每个月平均发表相关领域文章超过30篇，影响因子超过10分每个月有3-4篇以上。

今天为大家奉上最简单也是最容易被人忽略的转录组干货——转录组Q&A一百问。我们从生信、数据挖掘、后期验证、多组学联合分析等多个维度为大家解答各种可能做转录组测序会遇到的高频问题，希望能够为大家解决转录组的各种问题。

1.快速筛出候选基因

       转录组测序可以帮助研究者全局且快速地对某一生命过程中所有发生变化的基因进行筛查，然而只有从全局的基因变化分析到定位具体基因的讨论，研究者才能有的放矢开展后续的基因表型验证和机制探索。我们对快速筛出候选基因提供以下通用性的建议和方法，您也可以结合其他筛选指标一起筛选出关键基因。

1.1 筛选思路

       候选基因应该与您所研究的生物学问题密切相关，并在某一生命过程中发挥核心的生物学功能、参与重要的通路。基因功能的背景信息可以来源于参考基因组注释、基因功能数据库（如GO、KEGG）注释或者相关研究的文献和综述等，对于大且复杂的研究方向，比如“纤维化”、“植物盐胁迫”等，建议通过文献综述（可以在PubMed检索，https://pubmed.ncbi.nlm.nih.gov/）获取可能相关的基因、GO和KEGG通路信息，然后在自己的测序结果中进行筛选。无论通过怎样的方式进行数据筛选，人为的解读和挑选是必不可少的。有时候差异基因过多，还是不得不用肉眼去挖掘每个基因所代表的意义。这个工作不可省略，我们不得不去做这个繁琐但必要的工作！很多大课题组筛选数据也没有其他技巧，都是挨个儿看挨个敲除每个基因来验证。如果暂时没有思路，可以考虑一些调控上游的因素，如转录因子、组蛋白修饰、DNA结合蛋白、RNA结合蛋白等。

1.2 筛选方法

1）基于GO富集分析结果

       打开GO富集结果*_GO_enrichment.xlsx文件，在GO_Term列中，用关键词逐一筛选与您研究相关的生物学功能，通常设置pvalue小于0.05。获得相关显著富集的GO Term后建议在基因差异分析表格*_Gene_differential_expression.xlsx中的GO列进行筛选（使用GO ID “GO:0030335”或GO Term “positive regulation of cell migration”进行筛选，由于很多GO Term名称高度相似，建议优先使用GO ID进行筛选），可以得到注释到相关GO Term的所有基因，再筛选significant列显示为“yes”的即可获得注释到相关GO Term且差异显著的基因。由于*_Gene_differential_expression.xlsx包含所有基因的表达量信息和差异分析信息（差异倍数、p值、q值），因此在*_Gene_differential_expression.xlsx中进行筛选可以获得某GO Term中的基因表达全貌，有利于捕获那些不满足默认差异阈值条件（|log2FC|>=1 & q<0.05）但同样具有统计学差异（q<0.05或p<0.05）的基因。示例如下：

       （使用Excel或WPS的筛选功能，从GO富集结果文件中筛选“positive regulation of cell migration”功能，找到11个差异表达的基因，随后可以在*_Gene_differential_expression.xlsx中的GO列筛选获得这11个差异表达基因的信息。

详情请点击：

2）基于KEGG富集分析结果

打开KEGG通路富集结果*_KEGG_enrichment.xlsx文件，在pathway_name列中，用关键词逐一筛选与您研究相关的通路，通常设置pvalue小于0.05。获取关注通路的富集情况后可以参照GO筛选方法在*_Gene_differential_expression.xlsx中的KEGG列检索相关通路，即可获得注释到关注通路的全部基因和差异显著基因的信息。以查找“PI3K-Akt signaling pathway”为例，示例如下：

       （使用Excel或WPS的筛选功能，从KEGG通路富集结果文件中筛选“PI3K-Akt signaling pathway”通路，有17个差异表达基因与“PI3K-Akt signaling pathway”相关，随后可以在*_Gene_differential_expression.xlsx中的KEGG列筛选获得这17个差异表达基因的信息。）

       请注意：在某些实验处理条件下，差异表达基因可能较少，因而GO和KEGG通路富集结果较少，如果已经有关注的GO或KEGG通路，可直接从*_Gene_differential_expression.xlsx中的GO列或KEGG列进行筛选，获取候选基因信息。

详情请点击：

3）基于GSEA分析结果

       GSEA从基因集的富集角度出发，理论上更容易囊括细微但协调性的变化对生物通路的影响。我们可以基于GSEA分析结果筛选在实验组或对照组中富集的基因集（GO条目、KEGG通路等），然后可以关注CORE ENRICHMENT（对ES值有主要贡献的基因），即Leading edge subset。

详情请点击：

4）缩小范围

       综合上述两种方法，初步确定候选的基因范围。接下来通过以下条件，进一步缩小候选基因范围（仅供参考）。

       表达丰度：受系统噪音影响，低表达丰度基因（FPKM值<10）的差异倍数并不可靠，建议选择表达丰度中等水平的基因进行后续研究。

       基因新旧：在PubMed上查询候选基因的相关文章数量，建议选择较新的基因（文献数量<100）开展深入研究，请注意同一个基因常存在多个别名。

       分子大小：考虑到后续功能验证实验是否能开展，建议选择0.5~2.5K大小的基因。

       除了上述建议，您还可以关注测序结果里面的特殊基因，比如相关修饰酶、转录因子信息，其中动物转录因子可以前往AnimalTFDB（http://bioinfo.life.hust.edu.cn/AnimalTFDB#!/）查询，植物转录因子可以前往PlantTFDB（http://planttfdb.gao-lab.org/）查询。如果您研究的物种在AnimalTFDB和PlantTFDB中没有收录，可以通过序列比对近源物种或模式物种转录因子的方式获得可能的转录因子信息。如果有类似需求，可以将需求（附上项目编号）和需要参考的物种信息反馈至售后技术支持邮箱（[e-mail]）。

       通过蛋白质-蛋白质相互作用（PPI）网络分析可以在STRING数据库中查询关注基因或基因集的相关性蛋白，以及哪些蛋白处于网络的核心。

详情请点击：

5）组合筛选

有时候我们需要在多个细胞系中探索敲除（或敲减）某基因后共同差异变化的基因，或者是过表达某基因或敲减某基因的体系中变化相反的基因，或者是动物模型或细胞模型中用药后表达有回复的基因，亦或者是关注响应时间序列或浓度梯度变化的基因，上述目的往往需要结合两个比较组或多个比较组的差异分析信息，以下筛选方法供参考：

（1）韦恩图

       韦恩图（Venn），也叫温氏图、维恩图，是用于显示集合重叠区域的关系图。我们可以用韦恩图筛选不同比较组（不同基因集）共同的元素和特有元素，比如共同差异表达基因、特异差异表达基因等，并将这种筛选过程和结果可视化。灵活利用韦恩图可以通过一次筛选或多次筛选获取我们关注的核心基因集，然后可以对核心基因集进行进一步分析。比如利用韦恩图筛选某基因敲减组差异上调基因和过表达组差异下调基因的交集，韦恩图详细方法介绍参见：

2. OmicStudio开挂绘图 | UpSet图原理讲解+实操演示（扫描下方二维码）

（2）表达模式分析

       当样本为时间节点样本、浓度梯度样本、治疗或用药前后设置时，若进行两两比较筛选，当样本组较多时可能需要多次的相交筛选才能获得目标变化趋势的基因集。基于STEM（http://www.cs.cmu.edu/~jernst/stem/）的基因表达模式聚类分析（趋势分析）可以快速获得我们关注趋势的基因集或者哪种变化趋势是最显著的变化趋势。OmicStudio已经上线表达模式分析工具（https://www.omicstudio.cn/tool/37），输入基因表达量数据即可进行表达模式分析。表达模式分析需要将所有组（或自己关注的组）按照时间节点排序、浓度梯度排序或对照、模型、治疗的顺序排序，输入的是基因表达量（比如FPKM值）信息，每组需要使用一个表达量表征该组的基因表达量，通常使用表达量平均值（对该组所有样本取表达量平均值）。关于OmicStudio上的表达模式分析可以参考：

（3）加权共表达网络分析

       加权共表达网络分析，即WGCNA，旨在寻找协同表达的基因模块（module），并探索基因网络（模块）与关注的表型之间的关联关系以及网络中的核心基因（Hub基因）。WGCNA适用于复杂的数据模式，如果样本数大于15（建议5组以上），可以考虑开展WGCNA分析。关于WGCNA介绍可以参考：

（4）富集分析

       我们默认是对差异显著基因进行富集分析，在我们更改差异显著阈值、韦恩图筛选、趋势分析和WGCNA分析后可以获得相关关注基因集信息，随后可以对关注基因集进行富集分析，以获悉与这些筛选基因相关的GO或通路。OmicStudio已上线通用版富集分析工具，可以适用于几乎所有物种。

详情请点击：

2. qPCR验证候选基因

       开展后续实验前，通常会使用qPCR对4.2中的候选基因进行定量检测，以验证是否与测序结果（变化趋势）一致。由于不同检测技术具有不同的检测误差，建议选择同时满足以下条件的差异表达基因进行qPCR验证。

       1）兴趣是第一位的。选择自己感兴趣的基因进行验证，即使其表达倍数可能没有到达2倍（比如1.5倍）。如果没有明确感兴趣的基因，可以参考以下几条；

       2）优先选择表达量变化幅度较大的基因，比如差异倍数 ≥ 2，即fc ≥ 2 or ≤ 0.5，或是log2(fc) ≥ 1 or ≤ -1；

       3）优先选择具有统计学差异的基因，比如qvalue < 0.05；如果满足此q值条件的基因较少，可以使用pvalue < 0.05进行筛选，对于基因群体而言，使用q值可以降低假阳性；对于某基因个体而言，如果p值小于0.05，也可以考虑进行后续验证与功能机制探索。如果自己的样本类型异质性较高且在进行测序时生物学重复数偏低，即使p值大于0.05，在其表达量尚可且变化幅度尚可的情况下也是可以考虑的；

       4）优先选择表达量比较高的基因，在比较的样本组中，至少有一组样本的基因表达丰度（FPKM值）平均值大于10（不是绝对标准）；

       5）建议优先选择生物学重复样本间一致性较高的基因（对于异质性较高的样本，基因在生物学重复样本间可能也有明显差异，但差异倍数是实验组和对照组表达量平均值相除，一些异常表达量可能对差异倍数产生影响）。

注意：

       1）由于测序结果只对测序样本负责，因此qPCR验证测序结果，请使用测序同批样本；

       2）您只需对候选基因进行qPCR验证，我们不建议您对非候选的差异表达基因或非差异表达基因进行qPCR验证，因为测序技术已较为成熟，不必要对技术可靠性进行验证；

       3）考虑到qPCR验证测序结果的吻合率通常为60~70%，如可能，建议选5~10个或更多的候选基因进行验证；

       4）如果样本类型异质性比较高，比如临床样本、动物模型、组织液等，那么注意验证样本的选择或者使用更多的样本（远大于3）进行qPCR验证，以判断某基因是否在两组中存在统计学差异。异质性较高的样本的组内表达平均值容易受偏离值或异常值的影响，从而影响差异倍数的计算，因此对于此种情况，不要只看差异倍数和p值，也应该重点关注组内表达量的一致性或两组间的整体差异；

       5）低表达基因（比如所有样本中FPKM值小于1）无论是测序定量和qPCR定量，结果可靠性都会降低，因此不建议选择低表达基因；如果您对其感兴趣、并且使用多对引物（针对不同外显子设计引物）在多于三次重复实验的基础上得到一致的结果，以您的qPCR结果为准。

       6）由于可变剪切的存在，一个基因可能存在多个转录本，针对于不同外显子设计引物的qPCR结果可能不同（甚至相反）。基因定量是基因座所有转录产物的综合定量结果，在设计引物时可以考虑引物设计在大部分转录本的共有序列区间。当然如果对某条mRNA（比如经典转录形式）或lncRNA感兴趣，引物尽可能设计在其特有序列区间。

3.常见问题

Q：3.1 原始数据FASTQ是什么？

A：fastq是一个文本格式，用于储存生物学序列及其相应质量值（通常是核酸序列的）。为了方便储存及可读这些信息，这些序列以及质量信息使用ASCII字符标示。该格式最初由Sanger开发，目的是将FASTA序列与质量数据放到一起，目前已经成为高通量测序结果的事实标准。通常fastq文件中每一个序列含有4行信息，第一行：以‘@’开头，是这一条read的名字，这个字符串是根据测序时的状态信息转换过来的，中间不会有空格，它是每一条read的唯一标识符，同一份FASTQ文件中不会重复出现，甚至不同的FASTQ文件里也不会有重复；第二行：表示序列信息，制表符或者空格不允许出现。一般是明确的DNA或者RNA字符，由A，C，G，T和N这五种字母构成，N代表的是测序时那些无法被识别出来的碱基；第三行：用于将测序序列和质量值内容分离开来。以‘+’开头，后面是描述信息等，或者什么也不加。在旧版的FASTQ文件中会直接重复第一行的信息，但现在一般什么也不加（节省存储空间）；第四行：测序read的质量值，每个字符与第二行的碱基一一对应，按照一定规则转换为碱基质量得分，进而反映该碱基的错误率，因此字符数必须和第二行保持一致，它描述的是每个测序碱基的可靠程度，用ASCII码表示。第四行中每个字符对应的ASCII值减去64，即为对应第二行碱基的测序质量值。如果测序错误率用E表示，Illumina的碱基质量值用Qphred表示，则有下列关系：Qphred=-10log10E。

 为了节省空间，原始数据fastq一般以压缩形式fastq.gz（或fq.gz）储存、分发和上传公共数据库，如无必要，无需解压。

Q：3.2 一般来说，FPKM值多大算表达量高的？

A：FPKM值和物种、不同的发育生长时期都是有关系的，没有一个绝对和统一的标准。根据文献中的报道，一般认为0.5以下的受测序误差影响比较大，可以认为是没有表达，0.5-1属于低丰度表达基因，1-10属于中低丰度表达基因，10-50以上属于中高丰度表达基因，50以上属于高丰度表达基因，该标准仅供参考。

Q：3.3 MD5值是什么？

A：MD5全称是Message-Digest Algorithm 5，是一种被广泛使用的密码散列函数，可以产生出一个128位（16字节）的散列值（hash value），用于确保信息传输完整一致。MD5值等同于文件的ID，它的值是唯一的。如果文件内容（不是文件名）发生改变，那么MD5值会发生变化。由于二代测序原始数据一般较大，为了保证数据在传输过程没有发生损坏，可以使用MD5值进行文件校验。Windows平台可以使用notepad++生成文件的MD5值（工具>>>MD5>>>从文件生成...），Mac平台可以打开终端，定位到文件的位置（比如如果文件在桌面，可以使用指令cd Desktop定位到桌面），然后使用指令MD5 Control1_Data1.fq.gz生成Control1_Data1.fq.gz文件的MD5值，而使用MD5 *可以生成此目录下所有文本文件的MD5值。

Q：3.4 为什么分析结果中没有类似于人、小鼠中的基因名？

A：分析结果中的基因ID、基因名（symbol）和转录本ID都来源于参考基因组，一般只有部分物种（如小鼠、大鼠等）有类似人中的基因名（如GAPDH、Gapdh），当然即使是大鼠目前也有部分基因没有基因名而以基因ID代替。其他研究较多的物种可能有自己专门的网站和基因命名方式，比如拟南芥以AT1G开头的基因名，番茄以Solyc01g开头的基因名等。另外部分物种可能有多个数据库，比如水稻、茶叶等，不同数据库间的基因ID、基因名和转录本ID形式可能都不同。此外对于研究较少物种或新组装的基因组可能基因ID、转录本ID和基因名三者相同。

Q：3.5 怎么查找基因序列进行后续实验，比如qPCR？

A：所有转录本的序列信息已放入3_transcripts_sequence文件夹中，您可以使用基因对应的转录本ID在merged.fa中查询其序列。OmicStudio已上线序列批量提取工具，提交merged.fa文件和转录本ID列表，即可批量提取关注转录本的序列。另外对于Windows平台，可以使用notepad++软件（https://notepad-plus.en.softonic.com/）打开此文件（merged.fa），使用快捷键（Ctrl+F）检索转录本所在位置，将序列复制出来即可；对于Mac平台，可以使用CotEditor打开此文件，使用快捷键（command+F）检索转录本所在位置，将序列复制出来即可。

Q：3.6 为什么同一个基因（名）有多个表达量并且有的上调、有的下调？

A：联川生物针对定量和差异分析提供了基因层面和转录本层面两个维度的信息。summary/4_transcript_expression/5_isoforms_fpkm_expression.xlsx为转录本定量结果，summary/6_differential_expression/*/*_Transcript_differential_expression.xlsx（*为差异比较组）为转录本差异分析结果。由于可变剪切，一个基因可能转录多条转录本，一般情况下，来自于同一基因的转录本具有相同的基因名和基因ID，如果我们在转录本分析结果中查询基因名，查到的是基因对应转录本的表达量信息和差异信息，不同转录本的表达量和差异信息往往是不同的。

当然在基因定量结果（summary/4_transcript_expression/4_genes_fpkm_expression.xlsx、summary/6_differential_expression/*/*_Gene_differential_expression.xlsx）中可能检索基因名（极少情况）也能检索到多个结果，原因是不同的基因座具有相同的名称。

Q：3.7 为什么研究的是A组织（如肝脏），富集分析却富集到B组织（如神经）相关通路或疾病的？

A：一般基因是多功能性的，其在不同发育时期、不同组织和不同状态下可能具有不同的功能。基因富集分析是将所有显著性差异基因向Gene Ontology数据库（或其他数据库）的各Term映射，也就是会统计显著性差异基因在数据库中收录的所有功能（或涉及条目、通路）。如果一个基因差异显著，那么其涉及的所有通路都会被富集，只是通路的富集程度不同。因此如果出现显著富集通路与自己研究不相关或相关性不大的情况，可以直接忽略；或者，如果有兴趣的话，也可以探索其他组织中有明确功能报道的基因是否在我们研究的模型中发挥作用以及发挥怎样的作用。

除了上述情况，极少情况我们会遇到富集分析结果中出现不属于本物种的通路，比如动物样本富集分析结果中出现光合作用（Photosynthesis），原因可能是在进行基因的注释时没有区分物种。一般如果参考基因组中没有基因的GO、KEGG注释，同时GO、KEGG官网又没有收录此物种的情况下，通常策略是将基因序列与GO、KEGG收录物种的蛋白序列进行比对，以获得相关注释，然后基于此进行富集分析。如果在比对时没有区分物种或大类，可能会出现动物基因获得了植物的通路注释，从而在富集分析结果中体现出来。如果遇到此类问题，可以重新进行基因注释和富集分析。

Q：3.8 在填写立项信息和分析要求时将实验组和对照组填反了，会不会影响分析结果？

A：比较组格式为：<实验组>VS<对照组>，即分析流程判断VS前为实验组，VS后为对照组。实验组和对照组的顺序只是影响差异基因上下调的判定（上调和下调刚好对调），对p值和q值无影响，对差异基因的判定和富集分析结果无影响，只需要将差异倍数（fc）变为倒数(1/fc)，log2fc变为相反数（-log2fc），up和down对调即可。上下调是相对而言的（比较FPKM的大小可以判断在哪组或哪个样本中的表达量高），A相对于B下调也可以描述为B相对于A上调。

Q：3.9 怎样筛选关注通路（如NF-kappa B signaling pathway）的差异显著基因并绘制热图？

A：summary/6_differential_expression/*/*_Gene_differential_expression.xlsx是某比较组的差异分析结果，在表格中的KEGG列使用Excel或WPS筛选”NF-kappa B signaling pathway“即可筛选与NF-kappa B signaling pathway相关的基因，结合“significant”为“yes”的筛选可以获得与NF-kappa B signaling pathway相关的差异显著基因，然后复制基因名（gene_name列）和样本表达量列（FPKM值列）到新的Excel，将此Excel上传至OmicStudio的热图绘制模块，即可进行热图绘制（对于有生物学重复的情况，在“绘图前数据处理”部分选择“按行中心化”和“按行标准化”，数据处理可以选择“不处理”或其他两个选项）。

其他常见问题FAQ

1. 为什么差异分析表格中只有 100 行，99 个基因的信息？

答：我们的分析结果分为三个部分，即结题报告（网页版）、完整分析结果（单独的summary）和原始数据。结题报告是对结果的导读以及文件结构的说明，其支持文件夹 src 中包含完整分析结果的完整结构，只是针对大的文件表格取其前 100 行进行展示；如果您查看的是 src 中summary_part 中的结果，表格只有 100 行信息。数据挖掘可以基于网页版报告的说明在完整分析结果中进行。

2. 样本间的相关性有何意义，是如何计算的？

答：样本的相关性反映样本间的相似程度，样本相关系数越接近于 1，样本间的相似度越高。相关性系数的计算方法有三种， 即 pearson、spearman 和Kendall。联川生物使用 R 语言进行样本间 pearson 相关系数的计算，基于所有表达基因，不建议基于筛选的基因进行 pearson 相关系数的计算去反映生物学重复和样本或组间差异。

通常我们认为生物学重复间的相关性系数要大于生物学重复外的样本的相关系数，但由于异质性的存在，部分样本类型（如临床样本、动物模型、组织液、细胞pool 等）可能存在部分样本与其他组更为接近的情况，即整体相关性系数都很高，没有理想的按组分群，部分病人相关性可能与健康人更接近。上述情况是可能且正常的，如果相关性系数都很高的情况下，不是必须要剔除组内偏离的样本（一般剔除样本需谨慎），这也是对于这类样本我们一般建议增加生物学重复数的原因（比如 n 大于 6 或 10），即通过设置更多的生物学重复筛选更具代表性的差异表达基因并降低异质性的影响；如果某个体的确严重偏离，也可以将其剔除。

3. 差异基因数目是m 个，为什么 GO、KEGG 富集分析中的 TS 数目小于 m?

答：GO、KEGG 富集分析是针对具有GO 或KEGG 注释基因而言的，没有 GO 或KEGG 的基因不会作为GO 富集分析或KEGG

富集分析的前景和背景的，TS 指的是所有差异显著基因（m）中具有GO 或 KEGG 注释的，因此 TS 不大于 m。

4. 热图使用什么数据绘制的，从哪个文件中可以找到绘图数据?

答：热图使用表达量绘制的，如果我们需要绘制所有样本的基因热图可以前往 summary\4_transcript_expression，

4_genes_fpkm_expression.xlsx 是所有样本基因表达量文件，其中基因名列和所有 FPKM 列（表达量，其他组学也是上传表达量矩阵），将其复制到新 Excel 即可得到基因表达量矩阵。从中筛选我们关注的基因即可针对我们关注的基因进行热图绘制。

summary\6_differential_expression\*\*_Gene_differential_expression.xlsx 是*（差异比较组）的差异分析文件，取关注基因的基因名列和所有 FPKM 列（表达量）至新的Excel 中即可获得这些基因的热图绘制文件。热图绘制可以在 OmicStudio（https://www.omicstudio.cn/tool/4）进行绘制。

需要强调的是对于生物学重复的样本，我们一般对表达量取 Z 值，然后对 Z 值进行赋色，以直观比较基因在不同样本中的表达量高低。Z 值计算公式为：

Zsample-i =[(FPKMsample-i)-Mean(FPKM of all samples)] / [Standard deviation(FPKM of all samples)]

当对数据进行如下处理时等同于 Z 值，建议在绘图时参照如下勾选。其他的颜色、图片比例可以基于个性化选择进行调整。

备注：

1） 如果绘制热图时基因数目太多（比如超过 100 个），不建议显示基因名称，因为如果需要显示清楚每个基因的名称，图形可能很长。如果是此情况，可以不显示行名，并通过图片高度和宽度调整比例以体现选择基因的表达量变化趋势；

2） 一般聚类可以将表达模式相似的基因聚在一起，从而更好的通过热图直观看出表达趋势。如果需要输出的热图和输入的基因顺序或样本顺序一致，可以不对行（一般是基因）或列（一般是样本）进行举例。

5. 差异基因GO/KEGG 富集分析散点图怎么看？数据筛选的依据是什么？

答：以 GO 富集分析为例，GO 富集分析散点图是根据富集的显著性（pvalue）取 top20 的 GO Term 进行绘图，如果老师您要根据 GO 富集分析结果进行数据筛选，我们这边的建议是先看是否显著富集（pvalue < 0.05），再看富集程度（Rich factor=S gene number/B gene number），再看 S gene number。显著富集（pvalue < 0.05），富集程度高（Rich factor 值越大），并且 S gene number 数目多是最理想的目标（没有明确目标时的建议），S gene number 比较多的也可以列入筛选的范围，Rich factor 和S gene number两者需要兼顾，没有绝对的优先级。另外差异基因 GO 和 KEGG 富集分析，除了看富集的显著性、富集程度、富集的差异基因数目之外，重点还是要关注具体的 GO 和 KEGG 注释是否和老师您的研究方向一致（自己的研究兴趣是第一位的）。

6.差异基因GO 富集分析柱状图为什么始终都是 biological_process（25 个），cellular_componen（15个）molecular_function（10 个）？数据筛选的依据是什么？

答：GO 富集性分析结果柱状图反映在生物学过程（biological process）、细胞组分（cellular component）和分子功能（molecular

function）富集的GO_Term 上差异基因的个数分布情况。由于在生物学过程、细胞组分、分子功能这三种 GO_function 上富集的 GO_Term 数目比较多，无法把所有的注释结果都展示在一张图中，因此三种 GO_function 我们分别挑选 Top25、Top15、

Top10 进行展示，数据筛选的依据是：先筛选 GO_function（以 biological process 为例），再根据注释到 GO_Term 的差异基因数目（S gene number）从大到小降序排列，筛选出 Top25。当筛选到 Top25 时，如果出现几个GO_Term 对应的 S gene number 数目相等时，代码会随机选一个进行绘图，另外 S gene number 数目相等的 GO_Term 在柱状图上是随机排列的，没有先后顺序，所以 GO 富集分析柱状图横坐标跟 GO 富集分析结果表中 S gene number 降序排列的结果不是完全一一对应的。

7. 请问为什么差异表达分析结果中，最后significant 都是no?

答：差异表达谱中包含差异显著基因和非差异基因，老师您先选中首行，点击筛选，然后在 significant 那一列选“yes”，筛选之后的基因即为差异表达基因。在分析结果中同时满足 log2(fc) 绝对值大于等于 1 且 q 值小于 0.05 的基因标为 yes，否者标为 no。

8. 请问为什么差异表达分析结果中，有的基因q 值小于 0.05，定义为差异显著基因，有的 q 值小于 0.05，定义为非差异显著

基因?

答：分析结果中同时满足 log2(fc) 绝对值大于等于 1 且q 值小于 0.05 的基因定义为差异显著基因，如果某基因（您关注的基因）q 值小于 0.05，但差异倍数没到 2 倍，也是可以考虑下游分析验证的。在生信分析时需要设置阈值，那么如果某基因满足 q 小于 0.05，但是差异倍数为 1.99，那么分析会认为其差异不显著，如果差异倍数为 2.01，那么分析会认为其差异显著。事实上在生物学实验中，2.01 和 1.99 没有明显差别，如果对此类基因感兴趣，也是可以考虑将此基因纳入候选。

9. Coverage 怎么计算的？0-1,2-5 数值表示什么意思?

答：Coverage 是测序序列覆盖整个转录本的绝对深度，即在考虑测序深度的同时考虑了转录本各个位置被测序reads 覆盖的广度，0-1、2-5 等表示不同的覆盖度区间。

我们的测序结果中采用 StringTie 软件进行转录本组装、合并、定量，最后得到 Coverage 和 FPKM 值，Coverage 的计算方法描述如下：

an exon coverage value is calculated like this: all the read alignments intersecting the exon are considered and the coverage values for each base (i.e. the number of the read alignments covering that base) are added up and then divided by the length of that exon.

另外 Coverage 是没有进行标准化处理的测序原始数据，FPKM 值是标准化处理后的表达量，可直接用于后续分析。

10. 你们云平台上的通用版富集分析结果中的 GO 富集分析柱状图的纵坐标百分比是怎么计算的?

答：柱状图是根据百分比结果列出的，纵坐标 Percent of genes 计算公式如下：

percent=s_gene_number*fold(BP/CC/MF)/(max(s_gene_number)*1.1)

s_gene_number：表示该功能条目下显著基因的个数

max(s_gene_number)：显著基因数目最多的条目的显著基因数目

max(BP/CC/MF)：表示三大功能中各自的 s_gene_number 最大值

fold(BP/CC/MF)：表示 max(s_gene_number)的 1.1 倍与 max(BP/CC/MF)的比值再取整fold(BP/CC/MF)=int(max(s_gene_number)*1.1/max(BP/CC/MF)

11. KEGG 富集的表格中那个 p 值有什么用？是越小越好？大于一定值这个结果就不能用了吗？

答：关于 KEGG 富集的 p 值是利用超几何检验计算的，在散点图中 p 值和 Rich facrtor(S_gene_number/B_gene_number)以及

S_gene_number 都是可以作为参数筛选使用，选择依据推荐是 pvalue——Rich factor——S_gene_number。不能说 p 值越大越不能用，建议您根据您所关心的功能及三个筛选依据进行综合考量。

12. RNA-Seq 测序中，增加生物学重复个数和单个样本数据量对差异分析的影响是什么？

答：增加生物学重复个数和单个样品数据量，都可以改善定量的结果。随着生物学重复数（n）的增加，差异分析的假阳性率（FPR）变化不大，但真阳性率（TPR）在不断提高。即提高生物学重复数，差异表达基因的检测更加敏感；随着生物学重复数的增加，差异分析的真阳性率（TPR）在不断提高；而测序深度的提高对真阳性率（TPR）的提高没有生物学重复增加明显。建议在实验设计时，如果允许多设几个重复，特别是对于异质性较高的样本类型。

13. 链特异性建库是什么，有什么优势？

答：联川生物在转录组建库时采用链特异性建库（fr-firstrand），链特异性转录组测序（strand-specific RNA-seq/ssRNA-seq） 可以保留转录组测序时转录本的方向信息，即可以确定转录本是来源于基因组上面的正义链还是反义链。其构建文库的方法有多种，其中用的最普遍的即是 dUTP 方法。相对于传统转录组测序而言，链特异性文库在基因结构的确定，non-coding 转录本（例如lncRNA 和 antisense transcript）的鉴定，原核生物的操纵子（operon）鉴定以及转录本的基本定量方面，都具有绝佳的优势。

链特异性建库的关键就在于合成 cDNA 的第二链时，由 dUTP 代替dTTP，然后用 UDG 处理，第二链就会降解，而第一链保留下来, 继而测序。因此，测序得到的转录本序列信息，只是来源于第一链的。

链特异性建库有如下优势：

1） 定量更准确

由于链特异性测序方法可以区分转录本的来源，因此在计算某些转录本的表达量时，可以排除来自其互补链的转录本。

2） 可变剪切事件的检测更准确

因为链特异性文库可以排除反义链上 antisense 转录本的影响，可变剪接事件的检测假阳性更低。

3） Non-coding transcript 的检测

链特异性文库可以显著提高 non-coding transcript 的检出效率。对于 antisense 的 non-coding 转录本，如果用普通文库，是无法区分的；如果是基因间的 non-coding 转录本，普通文库无法确定转录本的方向。

4） 原核生物操纵子（operon）的预测

原核生物的基因是多顺反子的结构，反义转录本上的基因，如果不加区分，那么对应位置的基因表达量会计算不准确，并且预测 operon 以及基因结构也更不准确。

5） 组装结果更真实

一般的转录组组装出来的 unigene 既包括编码转录本，也包括一些非编码转录本（比如lncRNA），但是如果不区分正反链， 那么有互补配对关系的编码与非编码转录本会被组装成一条转录本。

链特异性文库详细情况，请点击：

参考文献：Parkhomchuk D, Borodina T, Amstislavskiy V, Banaru M, Hallen L, Krobitsch S, Lehrach H, Soldatov A. Transcriptome

analysis by strand-specific sequencing of complementary DNA. Nucleic Acids Res. 2009 Oct;37(18):e123. doi: 10.1093/nar/gkp596. Epub 2009 Jul 20. PMID: 19620212; PMCID: PMC2764448.

 14.GO、Pathway 富集分析中是否一定需要选择显著富集的通路?

答：建议不要仅仅基于 Pathway 富集分析的结果解读数据，人为的解读和挑选是必不可少的。因为生物数据的解读，在现阶段更多是生物学问题，而不是数学问题。原因大体如下：

1） 基因调控是个系统，不要仅仅看成一个一个孤立的 Pathway

基因调控是个系统，可以从两个层面进行解读：

a）1 个基因的改变可以造成整个系统的改变；举几个例子：

把 1 个生命活动必须的蛋白敲除后，整个细胞会发生紊乱。而植物抗病应激，也往往是 1 个受体蛋白识别了病原的外源蛋白，然后导致整个细胞系统的变化。

b）1 个基因往往有多个功能，但执行具体的功能往往是不同蛋白复合物共同作用的结果。

例如。基因 X 理论上在不同情况下，有可能参与 A、B、C 通路。在某个生物处理下，或许基因 X 只在A 通路里起作用。但如果进行基因注释的话，X 同样也会被注释到 B、C。所以，富集分析的结果总是会涉及特别多的通路。例如，研究人的项目，无论什么研究背景，常常会富集到帕金森综合症通路。不是你的材料真的得了帕金森综合症，只是那些与你实验处理相关的基因，在一定条件下也可以参与到帕金森综合症的过程，所以被注释到了这个通路里。

2） Pathway 富集分析的统计假设，并非在任何情况下都适用

Pathway 富集分析，在生物学上的假设是：1 个 Pathway 上游基因的改变，会导致下游相关基因改变，从而改变通路中大量基因的表达，达到统计学上富集的效果。但很多 Pathway 中，基因 A、B、C 并不是相互调控的关系，而是共同参与某个过程的不同部分。

例如，代谢物 X 的合成修饰。基因 A、B、C 分步骤参与合成的 3 个步骤。基因 A 给X 前体加了羟基，然后传递到下游；基因 B 又给 X 前体加了苯环，再传递到下游；基因 C 又给X 的前体加了个乙酰基 ，完成 X 的合成。那么，基因 A、B、C 是参与了的相同的通路。如果基因 A 发生表达量变化，会直接调控影响 B、C 的表达量变化吗？看来很有可能不会，所以从RNA-seq 差异分析的富集分析结果中，这个通路是不显著的。那么基因A 的表达变化是否有生物学意义？当然有，因为代谢物 X 的合成的确受影响了。

类似的例子，理论上 DNA 差异甲基化的结果，就不能看Pathway 富集分析的结果。1 个Pathway 中的 1 个基因的 DNA 甲基化变化，就足以改变这个通路的基因表达，而不需要整个通路的甲基化都发生变化。DNA 甲基化、组蛋白 CHIP-seq 的结果，其实只看功能注释或通路注释就足够了，不需要考虑富集。

所以，我们还是要观察、理解某个核心 Pathway 中基因的相互作用，才能判断其中的基因变化是否有生物学意义，而不

仅仅看富集分析的 p 值或 q 值。

3） 目前的 Pathway 是不完整的

目前 KEGG 等数据库收录的是已有的研究结果，但这些Pathway 的信息，远没有到达完善的水准。大部分通路只是了解1 个大概的调控途径，而中间有什么转录因子参与、是否还有其他代谢物的生成，都是不知道的。这些通路的完整性，也会影响 Pathway 富集分析结果。例如，基因 A 发生变化了，看起来下游基因没有变化。也许是还有其他的调控在起作用，只是这些调控作用现在还不知道而已。

总结：Pathway 和 GO 富集分析结果的解读，应该从生物学意义的角度出发，p 值和 q 值只是个参考而已，那些不显著的通路也值得解读（从功能注释的角度解读，而不是从富集分析的角度解读）。只要结果可以解释，有意义，不用太迷信 p 值。

 15.1. TPM 是什么，标准化方法比 FPKM 更好？

答：随着测序的发展，为了解决比较不同样本基因表达量的问题，发展出不同的基因表达量标准化方法，比如 RPKM、FPKM、TPM、SRPBM(用于 circRNA 标准化)等，其中最常见的三种方法是 RPKM、FPKM、TPM。无论是 RPKM、FPKM 还是 TPM都是基于基因的 reads count 标准化而来，都是基因的相对定量结果。

一．RPKM、FPKM、TPM 的影响因素

RPKM、FPKM 和 TPM 都是基于基因的 reads count 标准化得到的，而影响 reads count 除了基因的客观表达水平外（我们关注的），还包括基因的长度和样本的测序深度（简单理解为测序数据量）。

在同一个样本中，基因越长，随机打断得到的片段就越多，该基因被测到的概率就越大，比对到该基因的reads 就越多，从而

reads count 越大；不同样本里，样本的测序深度越高，同一基因被测到的次数越多，比对到该基因的 reads 数就越多，从而 reads

count 越大。

如果我们想比较同一样本不同基因的表达量，只需要标准化基因长度影响即可；如果我们想比较同一基因不同样本中的表达量，只需要标准化测序深度的影响即可；如果我们既想比较同一基因在不同样本中的表达量，又想比较同一样本不同基因的表达量，那么需要同时对基因长度和测序深度进行标准化。

二．RPKM、FPKM、TPM 的区别

早期的转录组测序为单端测序，以 RPKM（Reads Per Kilobase Million）值作为基因的标准化方法，其先将测序深度标准化， 然后将基因长度标准化。而随着双端测序技术（现在的转录组测序一般均为双端测序）的出现与普及，逐渐使用FPKM 替代

RPKM 进行基因的标准化。FPKM 使用 Fragments 概念代替 Reads，是目前 RNA 测序最广泛使用的标准化方法之一。

对于双端测序，paired-reads 表示一个 Fragments 从两端测序产生的 reads，通常以 read1 和 read2 表示。如果一对 paired-reads 都比对上了，那么这一对 pair-reads 称为一个 fragment；如果一个比对上了，另一个没比对上，那么这个比对上的 reads 就称为一个 fragment。

FPKM 与 RPKM 类似，同样是先将测序深度标准化，然后将基因长度标准化。当为单端测序时，FPKM 和 RPKM 是相同的。但是按照 FPKM 值进行标准化的方式有个特点，即每个样本的FPKM 值总和并不一致。基于此有研究者提出了 TPM 的标准化方法，以解决 FPKM 总和不一致的情况。

TPM 表示 Transcripts Per Million，其是先对基因长度标准化，再对测序深度标准化，与FPKM 正好相反。基因的 TPM 值等于基因 FPKM 值在样本中的占比（百分数）乘以 10⁶，每个样本的基因的 TPM 值总和均为10⁶。

三．一定需要使用 TPM 替换FPKM？

那么我们是否一定需要以 TPM 代替FPKM 进行标准化呢？答案是否定的，原因如下：

1） 当两个样品中的绝大多数基因都发生相同方向变化时（概率比较小），TPM 和 FPKM 均不可靠。特别是当超高表达量的某基因 A 表达量发生变化时，在测序结果中会测得基因 A 更多的 reads，其他基因的 reads 就会相对较少，使用 FPKM 值进行标准化时，会得到其他基因的表达量都降低的假象；此时如果使用 TPM 的标准化方法，其标准化表达量会进一步被影响。

2） TPM 的标准化方式影响的是基因的表达量数值，但是对于差异分析的 p 值和 q 值无影响，因为 p 值和 q 值的计算是基于的基因的 reads count，不是基于 FPKM 或TPM。

3） TPM 影响的是差异分析中基因的差异倍数。TPM 是 FPKM 的百分数，虽然各样本的 FPKM 值总和不一致，但 TPM 和 FPKM值是线性相关的，基于 FPKM 和TPM 值对于大部分基因的变化趋势判断是一致的。

4） 研究者基于 FPKM 值的定量结果能够验证出差异并基于此进行下游机制研究（这也是一直沿用的重要原因），而且目前大量文章依然采用 FPKM 作为基因表达量的标准化方式，特别是顶级期刊依然不断有采用FPKM 值方式的文章发表，采用FPKM 或 TPM 不是制约文章发表的因素。

5） 目前尚没有生信研究文章从整体上表明 TPM 在定量上比 FPKM 有什么优势，TPM 只是解决了 FPKM 总和不一致的问题。如果科学界内有广泛共识认为 TPM 比FPKM 更有优势，结果更准确、更可靠，那么提供相关测序和分析服务的公司肯定会快速跟上的，但是目前并没有广泛共识。

6） TPM 和 FPKM 可能对那些变化微弱的基因影响较大，因为两者定量后可能得出的结论有差异，比如上下调（差异倍数和

1 的比较；0.99 和 1.01 的结论是不同的）以及显著性（在满足统计学阈值的情况下，差异倍数与 2 或 0.5 的比较；1.99 和 2.01 的结论是不同的）。但是基于某个基因的验证结果，我们不能得出 TPM 和FPKM 谁更好的结论，因为换个基因您可能会得出相反的结论。比如在p < 0.05 或 q < 0.05 的基础上TPM 计算出的差异倍数为 2.01，FPKM 计算出的差异倍数为 1.99，以默认的差异显著阈值来看，TPM 分析出是差异显著的，FPKM 分析出是差异不显著的。但是 1.99 和 2.01 没有本质差别。又比如 TPM 算出某基因变化倍数为 1.8 倍，FPKM 算出变化倍数为 1.5 倍，不能因为 TPM 算出的差异倍数可能与qPCR 算出的差异倍数更接近，说 TPM 更准确。可能产生差异的根本原因是二代测序和 qPCR 这两个技术本身的差异。

测序是基于比对到基因的 reads 去衡量基因的表达水平，也就是基于基因的所有外显子；而 qPCR 只是基于基因的一两个外显子去衡量基因的表达水平。

7） 不用过度纠结于 FPKM 还是 TPM，应该将精力放在数据的深度挖掘和实验验证上。如果您依然想采用 TPM 值，那么可以基于 FPKM 值进行TPM 值的计算。某样本 S（Sample S）中基因 A（geneA）的TPM 值为：TPMgeneA = （FPKMgeneA / Total FPKMSample S） X10⁶（使用Excel 或 WPS 求取样本S 的 FPKM 值总和并计算TPM）。

附上两篇使用 FPKM 进行标准化的参考文献：

Factor DC, Barbeau AM, et al. Cell Type-Specific Intralocus Interactions Reveal Oligodendrocyte Mechanisms in MS. Cell. 2020 Apr 16;181(2):382-395.e21. doi: 10.1016/j.cell.2020.03.002. Epub 2020 Apr 3. PMID: 32246942; PMCID: PMC7426147.

Guo CJ, Ma XK, Xing YH, et al. Distinct Processing of lncRNAs Contributes to Non-conserved Functions in Stem Cells. Cell. 2020 Apr 30;181(3):621-636.e22. doi: 10.1016/j.cell.2020.03.006. Epub 2020 Apr 6. PMID: 32259487.

16. 构建基因 A 敲除的细胞系或动物模型，为什么测序结果中还能够检测到基因 A 或基因 A 的表达水平没有明显变化？

答：我们构建基因敲除体系时一般是在染色体（DNA）上进行操作，主要针对基因（蛋白编码基因）的某个外显子，通过局部的 DNA 序列测序（PCR 加一代测序），确定基因外显子产生了移码突变，从而确定基因被敲除，并且通过 WB 检测蛋白水平显著的降低。一般通过 DNA（基因局部的 DNA 序列）的一代测序和 WB 结果即可确定目标基因被敲除。但是 RNA-Seq 反映的是转录水平，如果我们对于目标基因的操作并没有显著影响基因转录，那么在 RNA-Seq 中是可以检测到目标基因的， 同时敲除组的基因 FPKM 值可能相对于对照组下调、稍上调或者没有明显差异。

如果我们针对的是目标基因的启动子区域（一般针对 RNA 基因），那么理论上，无论是转录水平还是蛋白水平都是下调的。但是如果针对的是外显子，那么比较转录水平是没有意义的，有 WB 结果支持即可。

那么为什么在转录组中检测到敲除基因的表达，但是蛋白却又明显下调呢？原因是虽然移码突变不一定会影响转录，但是转录出来的 RNA 不能翻译（相当于 mRNA 变成了 lncRNA），从而影响了蛋白水平。

如果我们研究的基因敲除模型并不是自己构建的，可以先咨询敲除是针对基因的哪个区域，再去分析可能的原因。

17. GSEA 分析结果和 GO、KEGG 富集分析结果结论不一致，我该相信哪个？

答：富集分析和 GSEA 分析的原理不同，因而分析结果或结论是很可能不同的。富集分析（包括 GSEA）只是通过分析获得相关功能或机制研究的参考，基于此提出相关假设，之后需要通过实验验证去证明或证否。如果结论一致或不一致，您可以选择其中的一个分析结果来引出您的假设。

相关阅读

点击下方图片进入云平台资料汇总：

所见即所得，绘图高规格

联川云平台，让科研更自由