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

Abstract

本发明公开了一种关于基因表达的RNA测序的质控方法及装置，该方法包括：对测序技术得到的测序片段分别进行DGE和RNA-Seq分析；DGE分析的结果和转录组分析的结果分别与qPCR的结果进行相关性分析；根据相关性分析结果，判断DGE分析和转录组分析在基因表达定量上的差异，并从DGE分析和转录组分析中选取一种测序分析方式；从所选取的测序分析方式获取的分析结果中选取一百万标签数据(1M reads)，进行基因表达的测序稳定性分析。本发明通过对基因片段进行相关性分析和综合性评估，从而选取具有较高可靠性的基因表达分析手段，确保测序工作的准确性，为生产的稳定性提供质量控制方案。

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表达谱、转录组测序、染色体免疫共沉淀(ChIP-chip，Chromatin Immunoprecipitation)以及DNA甲基化等方面的应用。

数字基因表达谱(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)，进行基因表达的测序稳定性分析的步骤进一步包括：从数字基因表达谱分析结果中取出一百万标签数据(1Mreads)，并将其与全部的数字基因表达谱分析结果进行相关性分析；和/或从转录组分析结果中取出一百万标签数据(1M reads)，并将其与全部的转录组分析结果进行相关性分析。

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

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

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

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

附图说明

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

图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结果的分析结果示意图；

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

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

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

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

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

具体实施方式

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

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

如图1所示，关于基因表达的RNA测序的质控方法100包括步骤102，对测序技术得到的测序片段分别进行数字基因表达谱分析(DGE)和转录组分析(RNA-Seq)。本发明实施例中，测序方法可以采用高通量测序技术，例如采用Illumina GA Solexa测序技术；Solexa是一种基于边合成边测序技术(SBS，Sequencing-By-Synthesis)的新型测序方法，通过利用单分子阵列实现在小型芯片(Flow Cell)上进行桥式PCR反应。新的可逆阻断技术可实现每次只合成一个碱基，不需要标记荧光基团，再利用相应的激光激发荧光基团捕获激发光，从而读取碱基信息。实验可以采用36Single 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分析方式所获取的分析结果中随机选取1Mreads，并将其与全部的转录组分析结果进行相关性分析；所述随机的选取方式可以是将所有测序得到的reads完全打乱，再从中任意取出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总RNA，利用Oligo(dT)磁珠吸附纯化mRNA，并以Oligo(dT)引导反转录合成双链cDNA。标签5’末端的产生可用两种内切酶实现：NlaIII或者DpnII，通常我们使用NlaIII，它识别并切断cDNA上的CATG位点，利用磁珠沉淀纯化带有cDNA3’端的片段，将其5’末端连接Illumina接头1(即序列：ACAGGTTCAGAGTTCTACAGTCCGACATG)。Illumina接头1与CATG位点的结合处是MmeI的识别位点，MmeI是一种识别位点与酶切位点分离的内切酶，酶切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碱基识别酶NlaIII、酶切双链cDNA，链接Illumina接头1，利用MmeI酶切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和RNA-Seq分析的结果与qPCR的结果的相关性分析方法做详细介绍：

数字基因表达谱分析(DGE)的结果与qPCR的结果的相关性分析，主要涉及DGE标准分析中表达量TPM(Transcripts Per Million clean reads)的计算方式，具体来说：TPM＝每个基因包含的原始Clean Tags数/该样本中总clean Tags数×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(a)示出了样品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.，30May 2008，Nature Methods|Advance Online Publication)如下所示：

$\mathrm{RPKM}=\frac{{10}^{6}C}{\mathrm{NL}/{10}^3}$

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

图3示出了本发明两样品的RNA-Seq分析的结果与qPCR结果的相关性分析结果的示意图。通常来说，如RNA-Seq数据产量为3Mreads，可以从样品测序数据中随机取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/query/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所示。

表11M测序量标准品与重复实验的数据相关性对照表

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

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

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

Claims (10)

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

对测序技术得到的测序片段分别进行数字基因表达谱分析(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)和转录组分析(RNA-Seq)；

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

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

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

9.如权利要求8所述的装置，其特征在于，所述相关性分析模块进一步包括：

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

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

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

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

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

Images