基于Proton的转录组测序数据的比较分析方法和系统的制作方法
【技术领域】
[0001] 本发明涉及生物信息技术领域，尤其涉及一种基于Proton的转录组测序数据的 比较分析方法和系统。
【背景技术】
[0002] 转录组（transcriptome)广义上指某一生理条件下，细胞内所有转录产物的集 合，包括信使RNA (mRNA)、核糖体RNA、转运RNA及非编码RNA ;狭义上指所有信使RNA的集 合。由于转录组代表了生物在某一时刻的基因表达状态，因此，对转录组的研宄具有极大的 生物学意义。
[0003] 基因可以具有产生mRNA的多个不同拼接的转录本，且mRNA也可以由基因组上的 其他区域产生。测序技术是用于调查和理解mRNA表达谱的有力工具。相对于传统的芯片 杂交平台，转录组测序无需预先针对已知序列设计探针，即可对任意物种的整体转录活动 进行检测，提供更精确的数字化信号，更高的检测通量以及更广泛的检测范围，且还允许发 现新的序列变异和转录本的可能性，是目前深入研宄转录组复杂性的强大工具。基于Ion Proton平台的RNA-Seq产品可以实现在2~4个小时完成转录组测序，产出高达IOGb的产 量，读长为100_200bp。测序成本方面，Ion Proton摆脱了昂贵的光学测序的束缚，测序成 本大大降低。较低的测序成本，更快的测序速度、更长的数据片段以及良好的一致准确性使 得Ion Proton平台的转录组测序在生命科学研宄领域，尤其是临床诊断方面的优势十分明 显O
[0004] 由于基于Proton的mRNA测序技术是高吞吐量的并产生相当大量的基因数据， mRNA测序数据可能需要使用各种计算技术以根据所获得的的大量数据对测序的mRNA进行 准确和有意义的推断。此外，来自不同试验（其可以是来自不同组织或不同实验条件下的 样本）的mRNA测序数据可能还需要被比较，分析相当大量的mRNA的数据是相当复杂的过 程。对不同试验之间的mRNA测序数据的分析检测会受到测序平台特征的影响，也受到用 于分析数据的计算技术的影响，需要大量的人工检查和验证后续工作。
【发明内容】
[0005] 针对现有技术中存在的问题，本发明的主要目的在于提供一种基于Proton的转 录组测序数据的比较分析方法，可以提高转录组测序数据比较分析的准确度和可靠性。
[0006] 一方面，本发明提供了一种基于Proton的转录组测序数据的比较分析方法，包 括：
[0007] 通过Proton测序平台获取某一物种至少两个转录组的原始测序数据；
[0008] 过滤所述各组原始测序数据中的不合格数据，获得所述各个转录组的待分析数据 (即clean reads，它是原始测序数据经过过滤处理后留下的数据）；
[0009] 对所述各个转录组的待分析数据进行第一步分析和第二步分析，其中，所述第一 步分析包括：
[0010] (1)将所述各个转录组的待分析数据分别比对到所述物种的参考基因组；
[0011] (2)利用RPKM标准化方法对所述各个转录组的待分析数据进行转录本定量分析；
[0012] (3)对转录本定量分析后的所述各个转录组的待分析数据进行两两差异比较，并 根据设定的阀值筛选出显著差异表达的基因；
[0013] (4)对所述显著差异表达的基因进行功能分析；
[0014] 所述第二步分析包括：
[0015] (1)选取合适的参数，将所述各个转录组的待分析数据分别比对到所述物种的参 考基因组；
[0016] (2)可变剪接分析：得到各个基因的可信的剪接位点，对剪接位点信息进行计算， 获得含有可变剪接形式的基因注释结果；和/或
[0017] 基因结构优化分析：根据比对结果进行基因结构的优化分析；和/或
[0018] 新转录本预测：根据比对结果进行新转录本的预测；和/或
[0019] SNP分析：根据比对结果进行SNP预测；和/或
[0020] 基因融合检测：根据比对结果进行融合基因检测。
[0021] 在本发明的一个实施例中，所述基因结构优化分析的步骤为：将比对上的reads 进行组装，过滤掉组装结果小于阀值的转录本，且过滤掉与所述物种已知的转录本同源性 较高的转录本，再将剩下的转录本与所述物种的参考基因组进行比对，如果某个转录本与 其对应的基因上游或下游的重叠部分（overlap)超过阀值，且与其对应位置的已知转录本 覆盖度大于阀值，则将该转录本作为最终结果，往其已知的注释结果的两端延伸。在本发 明的又一个实施例中，基因结构优化分析中过滤掉组装结果小于阀值的转录本中的阀值为 ISObp ;如果某个转录本与其对应的基因上游或下游的重叠部分（overlap)超过阀值中的 阀值为200bp ;且与其对应位置的已知转录本覆盖度大于阀值中的阀值为0. 9。
[0022] 在本发明的又一个实施例中，所述新转录本预测的步骤为：将比对上的reads进 行组装，过滤掉组装结果小于阀值的转录本，且过滤掉与所述物种已知的转录本同源性较 高的转录本，再将剩下的转录本与所述物种的参考基因组进行比对，如果某个转录本与其 对应的基因上游或下游的重叠部分（overlap)没有达到阀值，且满足一定的测序深度， 则称为新转录本。在本发明的又一个实施例中，新转录本预测中过滤掉组装结果小于阀 值的转录本中的阀值为180bp ;如果某个转录本与其对应的基因上游或下游的重叠部分 (overlap)超过阀值中的阀值为200bp。
[0023] 在本发明的又一个实施例中，不合格数据包括序列长度低于30的reads ;修剪掉 测序接头后序列长度低于30的reads ;修剪掉低质量的碱基序列后序列长度小于30的 reads，其中所述低质量碱基序列是指从3'端开始，以15为窗口统计reads平均质量，若平 均质量值小于9,则继续向5'端滑动，直到窗口平均质量大于9时该窗口 3'端到该reads 3'端的所有喊基序列。
[0024] 在本发明的又一个实施例中，第一步分析中的比对是使用Tmap软件进行的，所述 第一步分析中的所述两两差异比对是使用泊松分布模型进行的。
[0025] 在本发明的又一个实施例中，利用所述Tmap软件比对时，参数的具体设置如下： 设置比对结果输出类型a为2,即为所有最好的位点；设置程序运行时的线程数为1~16 ; 设置输出文件类型为压缩的bam格式；设置比对策略为map4,即为bwa快速比对。
[0026] 在本发明的又一个实施例中，第一步分析中的功能分析包括聚类分析、GO富集分 析和/或Pathway分析。
[0027] 在本发明的又一个实施例中，第二步分析中的所有比对是使用tophat软件进行 的。
[0028] 在本发明的又一个实施例中，利用所述tophat软件比对时，参数的具体设置如 下：设置比对reads的错配数为0~5 ;设置reads开gap最大的缺口数为0~20 ;设置 reads的最大的编辑距离为5~20 ;设置reads剪接比对时最大错配数为0~4 ;设置reads 比对最大插入长度为1~10 ;设置reads比对最大缺失长度为1~10 ;设置片段的比对最 大错配数为0~5 ;设置Bowtie2比对内设参数为1~3 ;设置线程数为1~16。
[0029] 在本发明的又一个实施例中，利用所述tophat软件比对时，参数的具体设置如 下：设置比对reads的错配数为3~5 ;设置reads开gap最大的缺口数为6~15,设置 reads的最大的编辑距离为10~20,设置reads剪接比对时最大错配数为2~4 ;设置 reads比对最大插入长度为5~10 ;设置reads比对最大缺失长度为5~10 ;设置片段的 比对最大错配数为2~5 ;设置Bowtie2比对内设参数为1~3 ;设置线程数为6~16。
[0030] 在本发明的又一个实施例中，可变剪接分析中的对剪接位点信息进行计算包括以 下步骤：利用自己编写的perl脚本对剪接位点信息进行计算。
[0031] 在本发明的又一个实施例中，基因结构优化分析中和新转录本预测中reads组装 是使用cufflink软件进行的；所述SNP分析中SNP检测是使用pyrohmmvar软件进行的；所 述基因融合检测是使用tophat-fusion软件进行的。
[0032] 在本发明的又一个实施例中，使用所述cufflink软件组装时，参数的具体设置如 下：设置参数u，对多reads进行矫正；设置线程数p为1~16。
[0033] 在本发明的又一个实施例中，使用所述pyrohmmvar软件进行SNP检测时，参数的 具体设置如下：滑动窗口的长度设为3~30 ;最小的变异质量值分数设为20~50 ;最小碱 基质量值设为5~20 ;最小的reads长度设为15~30 ;设置忽略indels I ;设置最小比对 打分为〇 ;设置最小的比对一致性为〇. 50~0. 90 ;设置用于计算SNP最小支持数reads的 系数为〇.