WO2012009952A1

WO2012009952A1 - 关于基因表达的rna测序质控方法及装置

Info

Publication number: WO2012009952A1
Application number: PCT/CN2011/001158
Authority: WO
Inventors: 彭智宇; 韩祖晶; 李英睿; 王俊
Original assignee: BGI Shenzhen Co Ltd
Current assignee: BGI Shenzhen Co Ltd
Priority date: 2010-07-22
Filing date: 2011-07-13
Publication date: 2012-01-26
Anticipated expiration: 2013-01-22
Also published as: CN101914619A

Description

关于基因表达的 RNA测序质控方法及装置技术领域

本发明涉及生物技术领域，尤其涉及一种关于基因表达的核糖核酸 ( RNA, RiboNucleic Acid )测序质控方法及装置。背景技术

基因表达是指基因片段脱氧核糖核酸（ DNA , Deoxyribonucleic acid ) 转录成信使核糖核酸（mRNA, Messenger RNA )及 mRNA翻译成蛋白质的过程。随着人类基因组计划（ HGP , Human Genome Project ) 全部核苷酸测序的完成，人类基因组研究的重心逐渐进入后基因组时代 ( Postgenome Era ) 向基因的功能及基因的多样性倾斜。通过对个体在不同生长发育阶段或不同生理状态下大量基因表达的平行分析，研究相应基因在生物体内的功能，阐明不同层次多基因协同作用的机理，进而在人类重大疾病如癌症、心血管疾病的发病机理、诊断治疗、药物开发等方面的研究发挥巨大的作用。它将大大推动人类结构基因组及功能基因组的各项基因组研究计划。

一直以来，基于分子杂交的方法对基因表达进行分析，从经典的核酸分子杂交方法（southern , northern blotting ) 到目前的基因芯片技术，都是应用已知核酸序列作为探针与互补的靶核苷酸序列杂交，通过随后的信号检测进行定性与定量分析。

新一代高通量测序技术的发明对生物学特别是基因组学的研究来说具有划时代的意义，它的高通量性使得对一个物种的转录本和基因组进行细致全貌的分析成为可能。随着 Solexa 测序技术的出现，使得高通量，低成本测序成为可能，并且与芯片技术的模拟信号相比，基于 Solexa 测序技术的表达分析避免了芯片技术中的交叉杂交、分析模型复杂以及灵敏度低等缺点；但是，由于高通量测序读长的限制，使其在对未知基因组进行从头测序（de novo sequencing ) 的应用受到限制，这部分工作仍然需要传统测序手段（读取长度可达到 850 碱基）的协助。而这并不影响高通量测序技术在 mRNA表达谱、 microRNA表达谱、转录组测序、染色体免疫共沉淀（ ChlP-chip ， Chromatin 数字基因表达谱（DGE， Digital Gene Expression Profiling ) 和转录组分析（RNA-Seq ) 是利用新一代高通量测序技术和高性能计算分析技术对某一物种特定组织和状态下的基因表达情况进行序列捕捉和精确解析的新方法。随着新一代测序技术的不断发展，对基因表达的研究也会更加深入，因此，需要对基因表达的分析手段进行相关性评估，从而排除由于分析手段自身的不准确或不稳定性所造成的分析误差，从而选取具有较高可靠性的基因表达分析手段，以便真实反映基因测序的准确性，确保评估可靠，从而保证产业可行性和生产的稳定性。发明内容

本发明要解决的一个技术问题是提供一种关于基因表达的 RNA测序质控方法及装置，通过对基因表达的分析为基因测序提供质量控制方案。

本发明的一个方面提供了一种关于基因表达的 RNA 测序的质控方法，该方法包括：对测序技术得到的测序片段分别进行数字基因表达语分析（DGE ) 和转录组分析（RNA-Seq ); 数字基因表达谱分析的结果和转录组分析的结果分别与实时定量基因扩增荧光检测（qPCR， Real-time Quantitative PCR Detecting System；其中 PCR, Polymerase Chain Reaction, 聚合酶链式反应）的结果进行相关性分析；根据相关性分析结果，判断数字基因表达谱分析和转录组分析在基因表达定量上的差异，并从数字基因表达谱分析和转录组分析中选取一种测序分析方式；从所选取的测序分析方式获取的分析结果中选取一百万标签数据（ 1M reads ), 进行基因表达的测序稳定性分析。

本发明提供的关于基因表达的 RNA 测序的质控方法的一个实施例中，该方法还包括：釆用高通量测序技术进行关于基因表达的 RNA 测序；对数字基因表达谱分析的结果和转录组分析的结果分别进行去接头序列和去低质量序列的处理。

本发明提供的关于基因表达的 RNA 测序的质控方法的一个实施例中，通过高通量测序技术对样品片段的基因表达进行多次测序，并对多次测序的数据取平均值以获得实时定量基因扩增荧光检测的结果。

本发明提供的关于基因表达的 RNA 测序的质控方法的一个实施例中，数字基因表达谱分析的结果和转录组分析的结果分别与实时定量基因扩增荧光检测（qPCR ) 的结果进行相关性分析进一步包括：当参考基因不全时，将数字基因表达醤分析的结果和转录组分析的结果分别与实时定量基因扩增荧光检测的结果进行相关性分析；和 /或在相同测序量的情况下，比较数字基因表达谱分析的结果和转录组分析的结果检测到的基因数。

本发明提供的关于基因表达的 RNA 测序的质控方法的一个实施例中，当参考基因不全时，将数字基因表达谱分析的结果和转录組分析的结果分别与实时定量基因扩增荧光检测的结果进行相关性分析的步骤进一步包括：将参考基因从 3，端到 5，端平均切成三份；对三份参考基因分别进行数字基因表达谱分析和转录组分析；将所获得的分析结果分别与实时定量基因扩增荧光检测的结果进行相关性分析。

本发明提供的关于基因表达的 RNA 测序的质控方法的一个实施例中，在相同测序量的情况下，比较数字基因表达谱分析的结果和转录组分析的结果检测到的基因数的步骤进一步包括：从高通量测序得到的测序片段中取出三百万标签数据（3M reads )分别进行数字基因表达谱分析和转录组分析；从高通量测序得到的测序片段中取出两百万标签数据 ( 2M reads )分别进行数字基因表达谱分析和转录组分析；和 /或从高通量测序得到的测序片段中取出一百万标签数据（1M reads )分别进行数字基因表达谱分析和转录组分析；在相同测序量的情况下，分别比较数字基因表达谱分析和转录组分析方法能够检测到的基因数。

本发明提供的关于基因表达的 RNA 测序的质控方法的一个实施例中，从所选取的测序分析方式获取的分析结果中选取一百万标签数据 ( 1M reads ), 进行基因表达的测序稳定性分析的步骤进一步包括：从数字基因表达谱分析结果中取出一百万标签数据（1M reads ), 并将其与全部的数字基因表达谱分析结果进行相关性分析；和 /或从转录组分析结果中取出一百万标签数据（1M reads ), 并将其与全部的转录组分析结果进行相关性分析。

本发明的另一个方面提供了一种关于基因表达的 RNA测序的质控装置，该装置包括：基因表达测算模块，用于对测序技术得到的测序片段分别进行数字基因表达谱分析（DGE )和转录组分析（RNA-Seq ); 相关性分析模块，用于将数字基因表达谱分析的结果和转录组分析的结果分别与实时定量基因扩增荧光检测（qPCR ) 的结果进行相关性分析；测序分析方式选取模块，用于根据相关性分析结果，判断数字基因表达谱分析和转录组分析在基因表达定量上的差异，并从数字基因表达谱分析和转录组分析中选取一种测序分析方式；测序稳定性分析模块，用于从所选取的测序分析方式获取的分析结果中选取一百万标签数据（ 1M reads ), 进行基因表达的测序稳定性分析。

本发明提供的关于基因表达的 RNA 测序的质控装置的一个实施例中，相关性分析模块进一步包括：第一相关性分析子模块，用于当参考基因不全时，将参考基因从 3，端到 5，端平均切成三份；对三份参考基因分别进行数字基因表达谱分析和转录组分析；将所获得的分析结果分别与实时定量基因扩增荧光检测的结果进行相关性分析；第二相关性分析子模块，用于在相同测序量的情况下，从高通量测序得到的测序片段中取出三百万标签数据（3M reads ) 分别进行数字基因表达谱分析和转录组分析，从高通量测序得到的测序片段中取出两百万标签数据（ 2M reads ) 分别进行数字基因表达谱分析和转录组分析；和 /或从高通量测序得到的测序片段中取出一百万标签数据（1M reads ) 分别进行数字基因表达谱分析和转录组分析；以及在相同测序量的情况下，分别比较数字基因表达谱分析和转录组分析方法能够检测到的基因数。

本发明提供的关于基因表达的 RNA 测序的质控装置的一个实施例中，测序稳定性分析模块进一步包括：第一测序稳定性分析子模块，用于从数字基因表达谱分析结果中取出一百万标签数据（1M reads ), 并将其与全部的数字基因表达谱分析结果进行相关性分析；第二测序稳定性分析子模块，用于从转录组分析结果中取出一百万标签数据（ 1M reads ), 并将其与全部的转录组分析结果进行相关性分析。

本发明提供了一种关于基因表达的 R A 测序质控方法及装置，通过对基因表达的分析手段进行相关性分析和综合性评估，从而选取具有较高可靠性的基因表达分析手段，真实反映基因测序的准确性，保证产业可行性，为生产的稳定性提供质量控制方案。附图说明

图 1示出本发明实施例提供的一种关于基因表达的 R A测序的质控方法的流程图；

图 2示出了本发明两样品的 DGE分析的结果与 qPCR结果的相关性分析结果的示意图，其中图 2 ( a )示出了样品 UHRR的 DGE分析的结果与 qPCR结果的分析结果示意图，图 2 ( b )示出了样品 HBRR的 DGE分析的结果与 qPCR结果的分析结果示意图；图 3示出了本发明两样品的 RNA-Seq分析的结果与 qPCR结果的相关性分析结果的示意图，其中图 3 ( a )示出了样品 UHRR的 RNA-Seq分析的结果与 qPCR 结果的分析结果示意图，图 3 ( b ) 示出了样品 HBRR 的 RNA-Seq分析的结果与 qPCR结果的分析结果示意图；

图 4示出本发明提供的关于基因表达的 RNA测序的质控方法的另一个实施例的流程图；

图 5示出了本发明样品 UHRR三等分参考基因序列的 DGE分析的结果与 qPCR结果的相关性分析结果的示意图，其中图 5 ( a ) 示出了样品 UHRR第一段的 DGE分析的结果与 qPCR结果的分析结果示意图，图 5 ( b ) 示出了样品 UHRR第二段的 DGE分析的结果与 qPCR结果的分析结果示意图，图 5 ( c ) 示出了样品 UHRR第三段的 DGE分析的结果与 qPCR结果的分析结果示意图；

图 6示出了本发明样品 UHRR三等分参考基因序列的 RNA-Seq分析的结果与 qPCR结果的相关性分析结果的示意图，其中图 6 ( a ) 示出了样品 UHRR第一段的 RNA-Seq分析的结果与 qPCR结果的分析结果示意图，图 6 ( b ) 示出了样品 UHRR 第二段的 RNA-Seq 分析的结果与 qPCR 结果的分析结果示意图，图 6 ( c ) 示出了样品 UHRR 第三段的 RNA-Seq分析的结果与 qPCR结果的分析结果示意图；

图，是本发明样品 UHRR在相同测序量下， DGE和 RNA-Seq检测到的基因数的示意图；

图 8示出本发明提供的关于基因表达的 R A测序的质控方法的另一个实施例的流程图；

图 9示出本发明实施例提供的一种关于基因表达的 RNA测序的质控装置的结构示意图；

图 10示出本发明提供的关于基因表达的 RNA测序的质控装置的另一个实施例的结构示意图；

图 11示出本发明提供的关于基因表达的 RNA测序的质控装置的另一个实施例的结构示意图。具体实施方式

下面参照附图用本发明的示例性实施例对本发明进行更全面的描述及说明。

图 1示出本发明实施例提供的一种关于基因表达的 RNA测序的质控方法的流程图。

如图 1 所示，关于基因表达的 RNA测序的盾控方法 100 包括步骤 102, 对测序技术得到的测序片段分别进行数字基因表达傳分析（DGE )和转录组分析（ R A-Seq )。本发明实施例中，测序方法可以采用高通量测序技术，例如采用 Illumina GA Solexa测序技术； Solexa是一种基于边合成边测序技术（SBS， Sequencing-By-Synthesis ) 的新型测序方法，通过利用单分子阵列实现在小型芯片（Flow Cell )上进行桥式 PCR反应。新的可逆阻断技术可实现每次只合成一个碱基，不需要标记荧光基团，再利用相应的激光激发荧光基团捕获激发光，从而读取碱基信息。实验可以釆用 36 Single End测序平台，对 RNA标准品 /实验样品分别进行双酶切测序和随机打断测序。

步骤 104, 数字基因表达讲分析的结果和转录组分析的结果分别与实时定量基因扩增荧光检测（qPCR )的结果进行相关性分析。稍后对关于 DGE 和 RNA-Seq分析的结果与 qPCR的结果的相关性分析方法做进一步的详细介绍。

步骤 106, 根据相关性分析结果，判断数字基因表达傳分析和转录组分析在基因表达定量上的差异，并从数字基因表达谱分析和转录组分析中选取一种测序分析方式。例如，综合分析数字基因表达谱和 RNA-Seq在基因表达定量（涉及基因数和基因表达量）上的差异，具体来说，可以包括分析正常测序量时比较数字基因表达傳和 RNA-Seq分析的结果与 qPCR结果的相关性，分析参考基因不全时比较数字基因表达谱和 RNA-Seq分析的结果与 qPCR结果的相关性，以及在相同测序量下比较数字基因表达傳和 RNA-Seq 能检测到的基因数中的至少任意一种方式。根据前述综合分析结果，得出 DGE和 RNA-Seq在基因表达定量上的差异，从而选取合适的测序分析方式。

步骤 108, 从所选取的测序分析方式获取的分析结果中选取一百万标签数据 ( 1M reads ), 进行基因表达的测序稳定性分析。例如，根据前述综合分析，如果 RNA-Seq分析方式所获取的基因表达定量更为准确（即 RNA- Seq得到的基因表达量更接近于 qPCR得到的基因表达量），则从 RNA-Seq 分析方式所获取的分析结果中随机选取 1M reads, 并将其与全部的转录组分析结果进行相关性分析；所述随机的选取方式可以是将所有测序得到的 read 完全打乱，再从中任意取出 1M的 reads; 如果 DGE和 RNA-Seq分析方式所获取的基因表达定量相当，则可以从中任选一种，以所选取的方式所获得的分析结果中选取 1M reads, 并将其与全部的转录组分析结果进行相关性分析；从而根据分析结果对生产测序的稳定性进行检测和评估以确保测序工作的准确性（其中关于 "检测和评估" 主要是通过分析测序结果的重复性，由于 1M reads 的基因数目和表达量是确定的，如果某次测序与确定结果重复性不好就说明该次测序不稳定不正确）。

数字基因表达语分析 ( DGE ) 实验部分主要包括：样本制备实验和测序实验。主要试剂耗材为 Illumina Gene Expression Sample Prep Kit和 Solexa 测序芯片（ flowcell )，主^ 器是 Illumina Cluster Station ( Illumina公司）和 Illumina Genome Analyzer ( Illumina公司 ) 系统。具体实验流程：提取 6 g总 R A，利用 Oligo(dT)磁珠吸附纯化 mRNA, 并以 Oligo(dT)引导反转录合成双链 cDNA。标签 5，末端的产生可用两种内切酶实现： Nlalll或者 DpnII, 通常我们使用 Nlalll, 它识别并切断 cDNA上的 CATG位点，利用磁珠沉淀纯化带有 cDNA3，端的片段，将其 5，末端连接 Illumina接头 1 (即序列： ACAGGTTCAGAGTTCTACAGTCCGACATG )。 Illumina接头 1与 CATG位点的结合处是 Mmel的识别位点， Mmel是一种识别位点与酶切位点分离的内切酶，酶切 CATG位点下游 17bp处，这样就产生了带有接头 1 的 Tag。通过磁珠沉淀去除 3，片段后，在 Tag3，末端连接 Illumina接头 2

(即序列： CAAGCAGAAGACGGCATACGANN ), 从而获得两端连有不同接头序列的 21bp标签 library. 经过 15个循环的 PCR线性扩增后，通过 6%TBE PAGE胶电泳纯化 85碱基条带，解链后，单链分子被加到 Solexa 测序芯片（flowcell )上并固定，每条分子经过原位扩增成为一个单分子簇

( cluster )测序模板，加入 4色荧光标记的 4种核苷酸，采用边合成边测序法（sequencing by synthesis, SBS ) 测序。每个通道将产生数百万条原始 Read, Read的测序读长为 35bp。利用 OligodT的 beads富集总 RNA 中 mRNA, 并逆转录为双链 cDNA, 釆用 4 碱基识别酶 Nlalll, 酶切双链 cDNA, 链接 Illumina接头 1，利用 Mmel酶切 3，端 CATG下游 17bp碱基，并在 3，端链接 Illumina接头 2。再加入 Primer GX1和 Primer GX2进行 PCR扩增。扩增后样本通过 6% TBE PAGE胶回收 85碱基条带，纯化后通过 Illumina基因表达测序法测序。转录组分析 ( RNA-Seq ) 实验部分测序基本过程包括：提取样品总 RNA后，用带有 Oligo(dT)的磁珠富集真核生物 mRNA (若为原核生物，则用试剂盒去除 rRNA 后进入下一步）。加入 fragmentation buffer将 mRNA打断成短片段，以 mRNA为模板，用六随机引物（random hexamers )合成第一条 cDNA链，然后加入緩冲液、 dNTPs、 RNase H 和 DNA polymerase I合成第二条 cDNA 链，再经过 QiaQuick PCR试剂盒（ Qiagen公司生产）纯化并加 EB緩冲液洗脱之后做末端修复并连接测序接头，然后用琼脂糖凝胶电泳进行片段大小选择，最后进行 PCR扩增，使用建好的测序文库进行测序。

接下来对 DGE和 R A-Seq分析的结果与 qPCR的结果的相关性分析方法做伴细介绍：

数字基因表达语分析 ( DGE )的结果与 qPCR的结果的相关性分析，主要涉及 DGE标准分析中表达量 TPM ( Transcripts Per Million clean reads ) 的计算方式，具体来说： TPM =每个基因包含的原始 Clean Tags t/该样本中总 clean Tags 数 ^x 1,000,000 (参见 Deep sequencing-based expression analysis shows major advances in robustness, resolution and inter-lab portability over five microarray platforms, Peter A.C. 't Hoen, Yavuz Ariyurek, et al" Nucleic Acids Research, 15 October 2008， Vol.36, No.21 )。

图 2示出了本发明两样品的 DGE分析的结果与 qPCR结果的相关性分析结果的示意图。通常来说，如 DGE数据产量为 3M reads, 可以随机从样品测序数据中取 3M reads来进行 DGE结果准确性的分析；所述随机的选取方式可以是将所有测序得到的 reads 完全打乱，再从中任意取出 3M 的 reads. 由于 UHRR和 HBRR是 RNA标准样品，能下载获取的是该 RNA 标准样品的 qPCR结果，而样品测序数据是不能下载的，需要自行进行测序。图 2 )示出了样品 UHRR的 DGE分析结果与 qPCR结果的分析结果示意图，图 2 ( b )示出了样品 HBRR的 DGE分析结果与 qPCR结果的分析结果示意图；其中本发明使用的 UHRR是 Stratagene公司的 Universal Human Reference RNA ( UHRR ) 标准品， HBRR 是 Ambion 公司的 Human Brain Reference RNA ( HBRR )标准品。如图 2 ( a )所示，样品 UHRR的 DGE分析的结果与 qPCR结果的相关系数约为 0.3，如图 2 ( b ) 所示，样品 HBRR的 DGE分析的结果与 qPCR结果的相关系数约为 0.53 (其中 UHRR和 HBRR样品在 DGE分析中能检测到的基因数都是 716， UHRR和 HBRR样品在 qPCR 中能检测到的基因数都是 687 )。

转录组分析 ( RNA-Seq ) 的结果与 qPCR结果的相关性分析，主要涉及 RNA-Seq标准分析中表达量 RPKM ( Reads Per Kb per Million reads ) 的计算方式，具体来说： RNA-Seq 标准分析中表达量 RPKM 的算法（参见 Mapping and quantifying mammalian transcriptomes by RNA-Seq, Ali Mortazavi et al. , 30 May 2008， Nature Methods! Advance Online Publication 如下所示:

其中， RPKM(A)是关于基因 A的表达量， C为唯一比对到基因 A的 reads数， N为唯一比对到基因组的总 reads数， L为基因 A编码区的 ½ 数。

图 3示出了本发明两样品的 RNA-Seq分析的结果与 qPCR结果的相关性分析结果的示意图。通常来说，如 RNA-Seq数据产量为 3M reads, 可以从样品测序数据中随机取 3M reads来进行 DGE结果准确性的分析，所述随机的选取方式可以是将所有测序得到的 reads 完全打乱，再从中任意取出 3M 的 reads。图 3 ( a ) 示出了样品 UHRR 的 RNA-Seq 分析的结果与 qPCR结果的分析结果示意图，图 3 ( b )示出了样品 HBRR的 RNA-Seq分析的结果与 qPCR结果的分析结果示意图；其中样品 UHRR的 RNA-Seq分析的结果与 qPCR结果的相关系数约为 0.91，样品 HBRR的 RNA-Seq分析的结果与 qPCR结果的相关系数约为 0.86 (其中 UHRR和 HBRR样品在 RNA-Seq 分析中能检测到的基因数都是 872， UHRR 和 HBRR 样品在 qPCR 中能检测到的基因数都是 851 )。此外，需要说明的是：从样品 UHRR和 HBRR中抽取 3M reads进行 RNA-Seq与 qPCR相关性分析，与用全部数据计算得到的 RNA-Seq和 qPCR的相关系数相同，都分别是 0.91 和 0.86。由此可见，对于基因测序的数据量来说，其对 RNA-Seq 的定量分析几无影响或者说影响甚微。

本发明提供的关于基因表达的 RNA 测序的质控方法的一个实施例中，在对测序技术得到的测序片段分别进行 DGE和 RNA-Seq分析步骤之前，对数字基因表达谱分析的结果和转录组分析的结果分别进行去接头序列；进一步地，也可以对去接头序列的结果再进行去低质量序列的处理，从而获取能够用于标签的数据（clean tag ) 以进行后续分析。

本发明提供的关于基因表达的 RNA 测序的质控方法的一个实施例中，通过高通量测序技术对样品片段的基因表达进行多次测序，并对多次测序的数据取平均值以获得实时定量基因扩增荧光检测的结果。例如， UHRR和 HBRR样品的 qPCR数据是从 GEO (高通量基因表达， Gene Expression Omnibus ) 上下栽的，具体来说其下载路径： http://www.ncbi.nlm. nih.gov/geo/quer /acc.cgi?acc=GSE5350 , 其中 UHRR的 accession number是 GSM129638, 公开日是 2006年 9 月 8 日； HBRR的 accession number是 GSM129645, 公开日是 2006年 9月 8 日。对 UHRR和 HBRR样品分别进行多次测序（如 4次）的平行实验，并对该 4 次平行实验的关于基因数和基因表达量的结果取平均值以作为 qPCR定量结果。

本发明提供的关于基因表达的 RNA测序的质控方法，基于对测序片段进行 DGE和 RNA-Seq分析，并对 DGE和 RNA-Seq分析的结果与 qPCR 结果进行相关性的综合分折，从而选取适宜的测序分析方式进行基因表达的测序稳定性分析。通过本发明提供的关于基因表达的 RNA测序的质控方法的一个实施例，其能够真实反映基因测序的准确性，保证产业可行性，为生产的稳定性提供质量控制方案。

图 4示出本发明提供的关于基因表达的 RNA测序的质控方法的另一个实施例的流程图。

如图 4 所示，关于基因表达的 RNA 测序的质控方法 400 包括步骤 402、 404、 405、 406、 408，其中步楝 402、 406和 408可以分别执行与可以分别执行与图 1所示的步驟 102、 106和 108相同或相似的技术内容，为筒洁起见，这里不再赘述其技术内容。

如图 4所示，在步骤 402后，执行步骤 404, 当参考基因不全时，将数字基因表达傳分析的结果和转录组分析的结果分别与实时定量基因扩增荧光检测的结果进行相关性分析。参考基因都是现有数据库中已经拼接好的核酸序列（http：〃 www.ncbi.nlm.nih.gov/ ), 这些核^^列有艮多版本（由不同的研究机构，数据中心等等单位发布的），每个机构由于其技术水平的限制，所以发布的结果与基因的真实情况是有不同的，因此可能存在参考基因不全或不完整的情形。例如，在用 DGE分析的结果和 RNA-Seq分析的结果分别与实时定量基因扩增荧光检测（qPCR ) 的结果进行相关性分析时，当参考基因不全 /不完整时，可以采用如下方式进行相关性的分析。

具体来说，在非模式生物中，有理由怀疑参考基因序列的不全会造成 DGE 定量不准确；首先，将完整的参考基因序列（如 NCBI 中人的 refseq 基因）从 3，端开始进行三等分，然后将三等分的基因序列当作完整的参考基因序列，分别进行 DGE分析的结果与 qPCR结果的相关性分析。图 5示出了本发明样品 UHRR三等分参考基因序列的 DGE分析的结果与 qPCR结果的相关性分析结果的示意图，其中图 5 ( a ) 示出了样品 UHRR第一段的 DGE分析的结果与 qPCR结果的分析结果示意图，图 5 ( b ) 示出了样品 UHRR第二段的 DGE分析的结果与 qPCR结果的分析结果示意图，图 5 ( c ) 示出了样品 UHRR第三段的 DGE分析的结果与 qPCR结果的分析结果示意图；分析发现样品 UHRR的这三部分序列，其 DGE分析结果与 qPCR结果的相关系数分别约为 0.71， 0.39和 0.33 (用完整基因序列作分析时 DGE分析结果与 qPCR结果的相关系数为 0.76 ), 其在 DGE分析中能检测到的基因数分别是 774， 596和 435。同样地，对于采用 RNA-Seq分析方式进行基因序列表达的，也是先将完整的参考基因序列从 3，端开始进行三等分，然后将三等分的基因序列当作完整的参考基因序列，分别进行 RNA-Seq与 qPCR相关性分析。图 6 示出了本发明样品 UHRR三等分参考基因序列的 RNA-Seq分析的结果与 qPCR 结果的相关性分析结果的示意图，其中图 6 ( a ) 示出了样品 UHRR第一段的 RNA-Seq分析的结果与 qPCR结果的分析结果示意图，图 6 ( b ) 示出了样品 UHRR第二段的 RNA-Seq分析的结果与 qPCR结果的分析结果示意图，图 6 ( c ) 示出了样品 UHRR 第三段的 RNA-Seq 分析的结果与 qPCR结果的分析结果示意图；分析发现样品 UHRR的这三部分序列，其 RNA-Seq 分析结果与 qPCR 结果的相关系数分别约为 0.85 0.91和 0.84 (用完整基因序列作分析时 RNA-Seq与 qPCR的相关系数为 0.91 ), 其在 RNA-Seq分析中能检测到的基因数分别是 917、 911 和 896。

DGE 由于其自身的缺点，它无法检测出不含 CATG (或 GATC )位点的基因，倾向于得到每条 mRNA 最靠近 3，端的 Tag作为该 mRNA的标签，因此它对参考基因的要求比较严格，参考基因序列不完整对 DGE 结果的影响很大；而 RNA-Seq 是对 mRNA 进行随机打断，所以每条 mRNA 能够得到很多标签，对参考基因的依赖性不是很强，在参考基因不完整的情况下也能够得到比较准确的表达量信息。由此可见，基因序列不完整对 DGE分析结果影响较大，而对于 RNA-Seq 的分析结果影响不大；也就是说，对于参考基因不全时，如果使用 DGE进行分析，则优选釆用基因表达 3，端开始的笫一段；而进一步地，优选采用 RNA-Seq分析方式对基因片段进行表达分析。

步骤 405，在相同测序量的情况下，比较数字基因表达谱分析的结果和转录组分析的结果检测到的基因数。例如，在相同测序量下，比较 DGE和 RNA-Seq两种分析方式检测到的基因数具体可以包括：从高通量测序得到的测序片段中随机取出 3M的 reads数分别进行数字基因表达语和 RNA-Seq 分析，从高通量测序得到的测序片段中随机取出 2M 的 reads 数分别进行数字基因表达谱和 RNA-Seq 分析，从高通量测序得到的测序片段中随机取出 1M的 reads数分别进行 DGE和 RNA-Seq分析；所述随机的选取方式可以是将所有测序得到的 reads 完全打乱，再从中任意取出对应数量的 reads. 在相同测序量的情况下，可以从前述三种方式中任选至少一种方式来分别比较数字基因表达谱和 RNA-Seq 能检测到的基因数。图 7是本发明样品 UHRR在相同测序量下， DGE和 RNA-Seq检测到的基因数的示意图。如图 7所示，在相同测序量时 RNA-Seq能检测到的基因比 DGE多约 1000个基因。

图 8示出本发明提供的关于基因表达的 RNA测序的质控方法的另一个实施例的流程图。

如图 8 所示，关于基因表达的 RNA 测序的质控方法 800 包括步驟 802、 804、 806、 808、 809, 其中步骤 802、 804、 806 可以分别执行与可以分别执行与图 1所示的步骤 102、 104、 106相同或相似的技术内容，为筒洁起见，这里不再赘述其技术内容。

如图 8所示，在步骤 806后，执行步骤 808，从数字基因表达谱分析结果中随机取出一百万标签数据（1M reads ), 并将其与全部的数字基因表达谱分析结果进行相关性分析。所述随机的选取方式可以是将所有测序得到的 reads完全打乱，再从中任意取出 1M的 reads。

步骤 809，从转录组分析结果中随机取出一百万标签数据（ 1M reads ), 并将其与全部的转录组分析结果进行相关性分析。所述随机的选取方式可以是将所有测序得到的 reads 完全打乱，再从中任意取出 1M 的 reads。

关于步骤 808 和步骤 809 中的相关性分析方法，可以采用如下方式：每个 lane加入一个已知库（只测 1M左右 reads ), 通过比较已知库的测序结果检测测序稳定性。例如，本发明实施例中采用 1M测序量标准品相关性分析方法，两次 1M测序量标准品 UHRR与重复实验样品 UHRR 相关性对照表，如表 1所示。

表 1 1M测序量标准品与重复实验的数据相关性对照表表 1 表明，对 1M 测序量标准品与重复实验的数据进行分析比对发现， Gene相关性无论是 Spearman相关系数还是 Pearson相关系数都很高；由此可以说明测序结果是正常、可信任的，使用标准品检测测序稳定性的方法具有可行性的，能够通过对基因表达的分析为基因测序提供质量控制方案。用 1M reads 进行基因表达质控的方法在 RNA测序中的应用能够评估生产稳定性。

此外，需要说明的是，对于本发明的质控方法和装置而言，不论是从定量准确性的角度来说，还是从检测到的基因数目来说，又或者对参考基因的依赖性等角度来比较，在质控方案中采用 RNA-Seq分析方法具有比 DGE更准确地反映基因表达的优势。

图 9示出本发明实施例提供的一种关于基因表达的 RNA测序的质控装置的结构示意图。

如图 9所示，一种关于基因表达的 RNA测序的质控装置 900包括：基因表达测算模块 902、相关性分析模块 904、测序分析方式选取模块 906和测序稳定性分析模块 908。

其中，基因表达测算模块 902, 用于对测序技术得到的测序片段分别进行数字基因表达谱分析（DGE )和转录组分析（RNA-Seq );

相关性分析模块 904，用于将数字基因表达谱分析的结果和转录组分析的结果分别与实时定量基因扩增荧光检测（qPCR ) 的结果进行相关性分析。

测序分析方式选取模块 906, 用于根据相关性分析结果，判断数字基因表达谱分析和转录组分析在基因表达定量上的差异，并从数字基因表达谱分析和转录组分析中选取一种测序分析方式。

测序稳定性分析模块 908，用于从所选取的测序分析方式获取的分析结果中选取一百万标签数据（1M reads ), 进行基因表达的测序稳定性分析。

图 10示出本发明提供的关于基因表达的 RNA测序的质控装置的另一个实施例的结构示意图。

如图 10所示，一种关于基因表达的 RNA测序的质控装置 1000 包括：基因表达测算模块 1002、相关性分析模块 1004、测序分析方式选取模块 1006和测序稳定性分析模块 1008，其中基因表达测算模块 1002、测序分析方式选取模块 1006和测序稳定性分析模块 1008 可以是与图 9 所示基因表达测算模块 902、测序分析方式选取模块 906和测序稳定性分析模块 908相同或相似的功能模块。为简洁起见，这里不再赘述。

如图 10所示，相关性分析模块 1004进一步包括：第一相关性分析子模块和第二相关性分析子模块；其中

第一相关性分析子模块 10041，用于当参考基因不全时，将参考基因从 3，端到 5，端平均切成三份；对三份参考基因分别进行数字基因表达谱分析和转录组分析；将所获得的分析结果分别与实时定量基因扩增荧光检测的结果进行相关性分析。

第二相关性分析子模块 10042, 用于在相同测序量的情况下，从高通量测序得到的测序片段中随机取出三百万标签数据（3M reads )分别进行数字基因表达谱分析和转录组分析，从高通量测序得到的测序片段中随机取出两百万标签数据（2M reads )分别进行数字基因表达谱分析和转录组分析；和 /或从高通量测序得到的测序片段中随机取出一百万标签数据（1M reads )分别进行数字基因表达谱分析和转录组分析；所述随机的选取方式可以是将所有测序得到的 reads完全打乱，再从中任意取出对应数量的 reads. 以及在相同测序量的情况下，分别比较数字基因表达傳分析和转录组分析方法能够检测到的基因数。

图 11示出本发明提供的关于基因表达的 RNA测序的质控装置的另一个实施例的结构示意图。

如图 11 所示，一种关于基因表达的 RNA测序的质控装置 1100 包括：基因表达测算模块 1102、相关性分析模块 1104、测序分析方式选取模块 1106和测序稳定性分析模块 1108。其中基因表达测算模块 1102、相关性分析模块 1104、测序分析方式选取模块 1106可以是与图 9所示基因表达测算模块 902、相关性分析模块 904、测序分析方式选取模块 906相同或相似的功能模块。为筒洁起见，这里不再赘述。

如图 11所示，测序稳定性分析模块 1108进一步包括：第一测序稳定性分析子模块 11081和第二测序稳定性分析子模块 11082, 其中

第一测序稳定性分析子模块 11081，用于从数字基因表达谱分析结果中随机取出一百万标签数据（1M reads ), 并将其与全部的数字基因表达谱分析结果进行相关性分析。所述随机的选取方式可以是将所有测序得到的 reads完全打乱，再从中任意取出 1M的 reads。

第二测序稳定性分析子模块 11082, 用于从转录组分析结果中随机取出一百万标签数据（1M reads ), 并将其与全部的转录组分析结果进行相关性分析。所述随机的选取方式可以是将所有测序得到的 reads完全打乱，再从中任意取出 1M的 reads。

本发明提供的关于基因表达的 RNA测序的质控装置，通过基因表达测算模块对基因片段进行分析，并通过相关性分析模块和测序分析方式选取模块进行相关性分析和综合性评估，从而选取具有较高可靠性的基因表达分析手段，真实反映基因测序的准确性，为生产的稳定性提供质量控制方案。

参考前述本发明示例性的描述，本领域技术人员可以清楚的知晓本发明提供的关于基因表达的 R A测序的质控方法及装置所具有的前述优点，本发明提供的质控方案适用于高通量测序技术，能够有效地评估 RNA 测序的稳定性，确保测序工作的准确性。

本发明的描述是为了示例和描述起见而给出的，而并不是无遗漏的或者将本发明限于所公开的形式。很多修改和变化对于本领域的普通技术人员而言是显然的。本发明中描述的功能模块以及功能模块的划分方式仅为说明本发明的思想，本领域技术人员根据本发明的教导以及实际应用的需要可以自由改变功能模块的划分方式及其模块构造以实现相同的功能；选择和描述实施例是为了更好说明本发明的原理和实际应用，并且使本领域的普通技术人员能够理解本发明从而设 ϋ

Claims

权利要求书

1. 一种关于基因表达的 R A测序的质控方法，其特征在于，所述方法包括：

对测序技术得到的测序片段分别进行数字基因表达语分析 ( DGE )和转录组分析 ( RNA-Seq ); 量基因扩增 ^检测（qPCR )的结果进行相关性分析； ' 、根据相关性分析结果，判断数字基因表达谱分析和转录组分析在基因表达定量上的差异，并从所述数字基因表达谱分析和转录组分析中选取一种测序分析方式；

从所选取的测序分析方式获取的分析结果中选取一百万标签数据 ( 1M reads ), 进行基因表达的测序稳定性分析。

2. 如权利要求 1所述的方法，其特征在于，所述方法还包括：采用高通量测序技术进行关于基因表达的 RNA测序；

对所述数字基因表达谱分析的结果和所述转录组分析的结果分别进行去接头序列和去低质量序列的处理。

3. 如权利要求 1 所述的方法，其特征在于，通过高通量测序技术对样品片段的基因表达进行多次测序，并对多次测序的数据取平均值以获得所述实时定量基因扩增荧光检测的结果。

4. 如权利要求 1 所述的方法，其特征在于，所述数字基因表达傳分析的结果和所述转录组分析的结果分别与实时定量基因扩增荧光检测

( qPCR )的结果进行相关性分析进一步包括：

当参考基因不全时，将数字基因表达语分析的结果和转录组分析的结果分别与实时定量基因扩增荧光检测的结果进行相关性分析；和 /或

在相同测序量的情况下，比较数字基因表达傳分析的结果和转录组分析的结果检测到的基因数。

5. 如权利要求 4所述的方法，其特征在于，所述当参考基因不全时，将数字基因表达谱分析的结果和转录组分析的结果分别与实时定量基因扩增荧光检测的结果进行相关性分析的步骤进一步包括：

将所述参考基因从 3，端到 5，端平均切成三份；

对所述三份参考基因分别进行数字基因表达谱分析和转录组分析；将所获得的分析结果分别与实时定量基因扩增荧光检测的结果进行相关性分析。

6. 如权利要求 4所述的方法，其特征在于，所述在相同测序量的情况步骤进一步包括：

从高通量测序得到的测序片段中取出三百万标签数据 ( 3M reads )分别进行数字基因表达语分析和转录组分析，从高通量测序得到的测序片段中取出两百万标签数据（2M reads )分别进行数字基因表达谱分析和转录组分析；和 /或从高通量测序得到的测序片段中取出一百万标签数据 ( 1M reads ) 分别进行数字基因表达谱分析和转录组分析；

在相同测序量的情况下，分别比较数字基因表达傅分析和转录组分析方法能够检测到的基因数。

7. 如权利要求 1 所述的方法，其特征在于，所述从所选取的测序分析方式获取的分析结果中选取一百万标签数据 ( 1M reads ), 进行基因表达的测序稳定性分析的步骤进一步包括：

从数字基因表达醤分析结果中取出一百万标签数据 ( 1M reads ), 并将其与全部的数字基因表达谱分析结果进行相关性分析；和 /或

从转录组分析结果中取出一百万标签数据 ( 1M reads ), 并将其与全部的转录组分析结果进行相关性分析。

8. —种关于基因表达的 RNA测序的质控装置，其特征在于，所述装置包括：

基因表达测算模块，用于对测序技术得到的测序片段分别进行数字基因表达谱分析 ( DGE )和转录组分析 ( R A-Seq )；

相关性分析模块，用于将所述数字基因表达傳分析的结果和所述转录组分析的结果分别与实时定量基因扩增荧光检测（qPCR )的结果进行相关性分析；

测序分析方式选取模块，用于根据相关性分析结果，判断数字基因表达傳分析和转录组分析在基因表达定量上的差异，并从所述数字基因表达傳分析和转录组分析中选取一种测序分析方式；

测序稳定性分析模块，用于从所选取的测序分析方式获取的分析结果中选取一百万标签数据 ( 1M reads ), 进行基因表达的测序稳定性分析。

9. 如权利要求 8所述的装置，其特征在于，所述相关性分析模块进一步包括：第一相关性分析子模块，用于当参考基因不全时，将所述参考基因从 3，端到 5，端平均切成三份；对所述三份参考基因分别进行数字基因表达谱分析和转录组分析；将所获得的分析结果分别与实时定量基因扩增荧光检测的结果进行相关性分析；

第二相关性分析子模块，用于在相同测序量的情况下，从高通量测序得到的测序片段中取出三百万标签数据 ( 3M reads )分别进行数字基因表达语分析和转录组分析；从高通量测序得到的测序片段中取出两百万标签数据 ( 2M reads )分别进行数字基因表达谱分析和转录组分析；和 /或从高通量测序得到的测序片段中取出一百万标签数据 ( 1M reads )分别进行数字基因表达傳分析和转录组分析；以及在相同测序量的情况下，分别比较数字基因表达谱分析和转录组分析方法能够检测到的基因数。

10. 如权利要求 8所述的装置，其特征在于，所述测序稳定性分析模块进一步包括：

第一测序稳定性分析子模块，用于从数字基因表达语分析结果中取出一百万标签数据 ( 1M reads ), 并将其与全部的数字基因表达语分析结果进行相关性分析；

第二测序稳定性分析子模块，用于从转录组分析结果中取出一百万标签数据 ( 1M reads ), 并将其与全部的转录组分析结果进行相关性分析。