CN109559780A - 一种高通量测序的rna数据处理方法 - Google Patents

Abstract

本发明公开了一种高通量测序的RNA数据处理方法，所述数据处理方法以高通量测序的RNA读段建立数据库，去除读段中的接头引物和低质量碱基，多重序列对比识别内源污染序列和外源污染序列并清除，获得处理后的RNA数据。本发明不仅达到数据高质量化的目标，同时将更多的重点放在污染识别和去除的方法上，比较有针对性地且使用DNA条形码的方式来识别污染，而且在处理后能够保留大部分数据信息，缩小与真实情况的差距，经过组装并利用生物DNA条形码(barcode gene)建立隐马尔科夫模型(HMM)的方法，搜索识别污染来源，比对去除污染序列，使得处理过后的高质量数据达到研究所需的要求。

Description

一种高通量测序的RNA数据处理方法

技术领域

本发明涉及一种高通量测序的RNA数据处理方法，属于高通量测序质量控制领域。

背景技术

下一代测序技术(NGS)又称高通量测序，以高输出量和高解析度为主要特色，能一次并行对几十万到几百万条DNA分子进行序列读取，在提供丰富的遗传学信息的同时，还可以大大降低测序费用，缩短测序时间的测序技术。随着下一代测序(NGS)技术的发展， 测序数据实现高通量地产出，利用转录组测序(RNA-Seq)技术可以实现将整个微生物群 落中的全部DNA或者RNA片段进行测定，得到微生物宏基因组或宏转录组数据。然而，原 始数据中存在多种质量问题、污染问题，并不能直接用于下游分析，如差异表达分析。由 于高通量测序技术数据处理量大，处理内容繁杂，因此对于测序质量的控制、污染源的确 定与排除成为了一个重要的研究课题。测序质量的影响因素是多方面的，常见的影响因素 多数为操作中的误差，跨越日期和组处理数据的批次效应的主要来源已被确定为实验性的， 如来自DNA提取试剂盒，PCR批次或测序仪器。

在RNA提取的过程中，由于rRNA、tRNA占全部RNA量的很大一部分，在送测的样本中不可避免的会混入rRNA和tRNA；在PCR扩增过程中，如果插入片段较短，则被测片段可 能会被测通，读段中可能会包含接头引物序列；另外，PCR偏差会导致测序数据中包含重 复读段(duplicatedreads)；由于酶的活性、测序仪器以及实验操作等多种因素，测序数据中 序列会含有低质量碱基、GC含量偏差以及未知的污染序列。其中，未知来源的污染序列较 难处理，目前发表的质量控制工具均没有有效的处理方法。

目前，常见的测序数据处理方法大致有以下步骤，先去除测序接头以及引物序列，再 过滤低质量碱基数据、重复序列，然后比对过滤rRNA、tRNA和污染序列。目前，有一些 已经发表的质控软件，它们着重强调高通量测序数据的质量问题，包括碱基质量值低下、 污染序列、接头序列以及碱基组成的偏差等，如FastQC、Trimmomatic、Cutadapt、NGS QCToolkit、FASTX-toolkit、PRINSEQ等，可以根据不同的文件类型、测序方法、质控要求选 择不同的质控工具来高效地完成任务。

尽管在测序技术不断发展的背景下，各种类型的测序数据质控工具层出不穷，然而大 多数工具的用途主要集中在对测序数据进行修剪，却没有足够细节化，没有针对转录组测 序数据里出现的质量问题做出有效的解决。虽然也有特意为转录组测序数据设计的质控工 具，但是它们会有各种限制。比如RSeQC，它主要提供的是读段比对的质控概要统计结果， 在一定程度上依赖于UCSC(the University of California,Santa Cruz)的基因组浏览器 (http://genome.ucsc.edu/)，然而它缺少序列修剪和污染过滤的功能，而且运行速度慢。不仅 如此，测序数据中未知来源的污染序列也没有相应的工具用于识别和去除，这将会严重损 害测序数据中包含的大量信息。因此，从实际情况出发，迫切需要设计一套流程方案，有 针对性的整合各个质控工具，从多方面对转录组测序数据进行质量控制。

发明内容

针对现有技术存在的上述问题，本发明的目的是获得一种高通量测序的RNA数据处理 方法。

为实现上述发明目的，本发明采用的数据处理方法的技术方案如下：

所述数据处理方法以高通量测序的RNA读段建立数据库，去除读段中的接头引物和低 质量碱基，多重序列对比识别外源污染序列并清除，获得处理后的RNA数据。通过对内源 性和外源性的污染物去除，可以获得处理后纯净的RNA数据。

优选的，所述内源污染序列包括测序数据中的：16S rRNA序列、18S rRNA序列、23SrRNA序列和28S rRNA序列。

优选的，所述内源污染序列包括测序数据中来源于动植物的RNA序列。

优选的，所述数据处理方法包括如下步骤：

a)获得高通量测序的RNA-seq测序读段，确定可能的污染源基因组，并下载生物DNA 条形码；污染源是取样过程中可能混入的来源于动植物的RNA；比如采样环境是人的粪便， 则污染源可能是人的结肠细胞RNA序列以及大麦RNA序列；

b)查看步骤a)中的高通量测序读段是否有接头序列，评估读段各位置的碱基质量，并 评估读段5端和3端的质量；

c)处理低质量碱基和碱基比例波动，过滤内源性污染物；由于本发明的目的是获得高 纯度微生物宏转录组数据，则测序数据中的16S rRNA序列、18S rRNA序列、23S rRNA序 列以及28S rRNA序列被认为是内源性污染物并且被过滤；

d)采用MUSCLE、Trinity和HMMER识别外源污染序列；由于本发明的目的是高纯 度微生物宏转录组数据，因此测序数据中的动植物RNA序列被认为是外源性污染物并且被 过滤；

e)采用KneadData、SAMtools和Bowtie2去除外源污染序列。

具体的，所述真实数据的保存格式为FASTQ，可以从NCBI下载，也可以是转录组测序数据，然后根据测序物种信息以及可能的污染源信息下载相应的参考基因组序列；之后使用FASTQC和Trimmomatic两个软件分别进行测序数据的质量评估和测序数据质量修正，这步完成后的结果文件是去除了低质量碱基和读段的高质量FASTQ格式文件。然后上述得到的高质量FASTQ格式文件还包含了污染序列，需要进一步识别去除，使用RNA-QC-Chain过滤掉内源污染序列。最后是识别外源污染序列并且将其去除，这一步将动植物DNA条形码序列进行多重序列比对；转录组从头组装；通过动植物中特有的DNA条形码序列，对测 序数据中的来自于动植物的污染序列进行识别。识别出污染源后，需要下载污染源参考基 因组，然后使用KneadData、SAMtools、Bowtie2等工具将之前去除内源污染的序列输入， 得到去除外源污染的最终数据文件。

优选的，所述动物DNA条形码从BOLD数据库中下载，植物物种的DNA条形码序列 集采用宁康教授实验室数据。BOLD数据库收录了绝大部分已知物种的DNA条形码(Barcodegene)序列，在“DATA PORTAL”中搜索关键系“animal”，下载所有搜索结果 的FASTA格式的序列文件作为动物物种的DNA条形码序列集。

优选的，数据评估步骤采用FastQC进行。

优选的，数据评估包括：每个位置的碱基质量情况、读段每个位置的各碱基频率分布图、 读段的GC含量频率分布图和接头含量。

优选的，步骤c)具体为：采用Trimmomatic处理低质量碱基和碱基比例波动，采用RNA-QC-chain过滤内源污染序列。处理低质量碱基和碱基比例波动问题时，数据中的读段均去除其5端的12个碱基，然后以5个碱基为单位窗口设置阈值为20进行滑动，去除窗 口中碱基平均质量值低于阈值的读段，最后控制读段的长度不小于85个碱基。该处理方式 既解决了读段5端碱基比例波动问题，又对读段的质量进行了修正。Trimmomatic的输入文 件是FASTQ格式文件，输出文件是高质量读段数据的FASTQ格式文件。

优选的，RNA-QC-chain采用隐马尔科夫模型(Hidden Markov Model，HMM)过滤内源污染序列。即使用HMMER的隐马尔科夫模型搜索，从输入序列中预测出rRNA片段， 然后提取。

优选的，步骤e)具体为：先使用Bowtie2-build建立污染源参考基因组的索引，然后使 用KneadData同时输入多个参考基因组的索引文件，并且比较快速的得到结果。

识别外源污染序列时，先进行动植物DNA条形码序列的多重序列比对步骤，作为下一 步要建立的HMM数据库需要该多重序列比对文件作为输入；分别得到多重序列比对文件后，使用HMMER中的hmmbuild进行HMM数据库的建立，并且将这两个HMM格式文 件合并到一起，便于序列在其中搜索。

优选的，所述数据处理方法还包括差异性表达分析步骤，以保证高通量测序数据不受污 染物影响。

与现有技术相比，本发明有针对性的整合各个质控工具，从多方面对转录组测序数据进 行质量控制，对质量控制流程处理前后的数据在所有读段的平均质量情况、rRNA过滤情况、 各步骤数据中的读段对数量情况以及差异表达结果的差异情况等多方面进行了评估，得到 以下结论：质量预处理步骤可以有效提高数据质量。污染识别与去除步骤在误差允许范围 内，可以识别出来自动植物的污染序列并去除，且去除污染后的数据中残余的污染序列不 影响差异表达分析。本发明相比于一般的质控流程，不仅达到了数据高质量化的目标，同 时将更多的重点放在污染识别和去除的方法上，比较有针对性地且使用DNA条形码的方式 来识别污染，而且在处理后能够保留大部分数据信息，缩小与真实情况的差距，经过组装 并利用生物DNA条形码(barcode gene)建立隐马尔科夫模型(HMM)的方法，搜索识别 污染来源，比对去除污染序列，使得处理过后的高质量数据达到研究所需的要求。差异表 达分析结果表明研究中的质量控制流程可以作为一种比较合理的、创新的数据质量处理方 案。

附图说明

图1为本发明数据处理方法的流程图；

图2为本发明数据处理方法的示意图；

图3为本发明数据处理方法的HQRD数据集数据质量图；

图4为本发明数据处理方法的HQACD数据集数据质量图；

图5为本发明数据处理方法的各组数据每步处理后污染占比柱状图。

具体实施方式

下面结合实施例对本发明提供的一种高通量测序的RNA数据处理方法作进一步详细、 完整地说明。下面描述的实施例是示例性的，仅用于解释本发明，而不能理解为对本发明 的限制。

下述实施例中的实验方法，如无特殊说明，均为常规方法。下述实施例中所用的实验 材料如无特殊说明，均为市场购买得到。

本实施例以从美国国家生物技术信息中心(National Center forBiotechnology Information，NCBI，https://www.ncbi.nlm.nih.gov/)下载的来自人类粪便的人类肠道微生物 宏转录组数据以及人结肠细胞全部RNA数据和大麦转录组数据为例进行质量控制，具体操 作步骤如下：

一、数据集准备

1.下载数据集

本实验方案中采用的测试数据均为转录组测序数据，从美国国家生物技术信息中心 (National Center for Biotechnology Information，NCBI，https:// www.ncbi.nlm.nih.gov/)下载， 得到5个样本数据(如表1-1所示)，称为真实样本数据，其中3个原始的真实数据(Real data raw，RDR)分别为SRR769409、SRR769424、SRR769439，均来自人类粪便，是人类肠道 微生物宏转录组数据；其余2个样本(SRR1756570、SRR1722292)分别是人类结肠细胞的 全部RNA数据和大麦转录组数据，该2个样本数据被作为可能的动植物源污染。

使用NCBI的SRA数据处理软件包sratoolkit，运行SRA数据下载软件prefect，参数为 样本在SRA(Sequence Read Archive)数据库中的测序结果集(Runs)的序号，如prefectSRR769409，即得到样本的SRA格式文件(如SRR769409.sra)；然后运行格式转换软件fastq-dump，将SRA格式文件转换为FASTQ格式或者FASTA格式(如SRR769409_1.fastq,SRR769409_2.fastq)

表1-1样本(真实)数据信息表

2.测试数据模拟

为了便于检验实验方案的可行性，采取将原始的人为污染数据(Artificialcontaminated data raw，ACDR)作为测试数据的方式进行数据处理流程。人为污染数据一共有3组，即 3个样本数据，6个FASTQ格式文件(如表2-2所示)，由5个原始的真实数据按照不同的 情况混合。混合方式如下：(1)ACDR769409，将SRR769409和SRR1756570的读段对分别混合在一起，其中SRR769409中的读段作为目标读段(target reads)，SRR1756570中的读段作为来自人类的污染序列(contamination)；(2)ACDR769424，将SRR769424和SRR1756570的读段对分别混合在一起，其中SRR769424中的读段作为目标读段，SRR1756570中的读段作为来自人类的污染序列；(3)ACDR769439，将SRR769439，SRR1756570和SRR1722292的读段对分别混合在一起，其中SRR769439中的读段作为目 标读段，SRR1756570和SRR1722292中的读段分别作为来自人类(动物)和大麦(植物) 的污染序列。同时，对人为污染数据中的目标读段数量和污染读段数量之比进行统计，定 义污染率(Contamination rate，CR)为：

ACDR769409、ACDR769424、ACDR769439的污染率分别为85.6％、83.3％以及93.4％， 和真实的数据污染情况比较相近。

表1-2本文所用数据集概况

3.参考基因组

根据模拟数据的信息，设定人类样本的转录组数据和大麦样本的转录组数据作为来自动 物和植物的污染序列，并且实验方案采取的是读段比对到参考基因组的方式去除污染，所 以需要下载相关的参考基因组数据。

从NCBI下载人类参考基因组序列，版本为GRCh38，格式为FASTA，下载地址为ftp://ftp.ncbi.nlm.nih.gov/refseq/H_sapiens/annotation/GRCh38_latest/refseq_ identifiers/GRCh38_latest_genomic.fna.gz，压缩文件GRCh38_latest_genomic.fna.gz大小为916M，解压后文件 GRCh38_latest_genomic.fna大小为3.1G；从NCBI下载大麦参考基因组序列，版本为v1.3， 格式为FASTA，下载地址为 ftp://ftp.ncbi.nlm.nih.gov/ genomes/all/GCA/000/182/155/GCA_000182155.3_O.barthii_v1.3/GCA_000182155.3_ O.barthii_v1.3_genomic.fna.gz，压缩文件GCA_000182155.3_O.barthii_v1.3_genomic.fna.gz 大小为93.4M，解压后文件GCA_000182155.3_O.barthii_v1.3_genomic.fna大小为598M。

模拟数据在经过质控流程处理后，还要进行差异表达分析，该过程需要微生物参考基 因组，考虑到模拟数据中微生物来源于肠道，则下载已经整合的人类肠道微生物参考基因 组，下载地址为ftp://climb.genomics.cn/pub/10.5524/100001_101000/100064/ 1.GeneCatalogs/IGC.fa.gz， 压缩文件IGC.fa.gz大小为2.3G，解压后文件IGC.fa大小为7.7G。

4.生物DNA条形码下载

BOLD(http://www.boldsystems.org)是加拿大生物多样性基因组中心开发的基于云的 数据存储和分析平台。它由四个主要模块组成：一个数据入口、一个教育入口、一个BINS 登记(新兴物种)以及一个数据收集和分析工作台。BOLD的数据库收录了绝大部分已知 物种的DNA条形码(Barcode gene)序列，在“DATA PORTAL”中搜索关键系“animal”， 下载所有搜索结果的FASTA格式的序列文件作为动物物种的DNA条形码序列集，重命名 为COI.fa，文件大小为2.5M。植物物种的DNA条形码序列集使用宁康教授实验室已经整 个好的数据，重命名为ITS2.fa，文件大小为36M。

二、数据分析

1.数据质量评估

转录组测序数据在进行处理之前需要进行各项质量信息的获取，第一步使用测序数据质 量控制软件FastQC，从3组原始的人为污染数据和3组原始的真实数据中分别获取到基本 的统计信息、每个位置的碱基质量情况、读段每个位置的各碱基频率分布图、读段的GC 含量频率分布图、接头含量等。该步骤主要的内容有：

根据六组数据各自的接头含量图得知，6组数据均不含接头序列，故不用做接头序列处 理；三组原始真实数据的平均质量值基本都在Q30以上，属于高质量数据，但是在读段靠 近3端处碱基的质量值开始下滑，有部分读段的质量值小于Q10，需要进行修剪。

从原始真实数据的读段每个位置的碱基频率分布图分析(如图3所示)，该三组数据均 存在相同的问题，即读段前10个碱基位置的ATGC碱基比例波动较大，可能是由于测序反 应开始时脱氧核糖核酸连接酶和底物结合不稳定导致的。

查看另外3组原始人为污染数据读段每个位置碱基的平均质量图，发现有部分读段的质 量值低于20，需要适当修剪。根据原始人为污染数据的读段每个位置的碱基频率分布图可 知，三组数据在5端测序起始的10至12个碱基处碱基比例波动较大，可能同样是由于测 序反应开始时脱氧核糖核酸连接酶和底物结合不稳定导致的，需要进行处理。

2.数据质量修正

该步骤主要是针对上一步数据评估中提及的质量问题进行处理，主要用到的质量控制软 件是Trimmomatic。六组数据中的读段均去除其5端的12个碱基，然后以5个碱基为单位 窗口设置阈值为20进行滑动，去除窗口中碱基平均质量值低于阈值的读段，最后控制读段 的长度不小于85个碱基。该处理方式既解决了读段5端碱基比例波动问题，又对读段的质 量进行了修正。

3.内源污染序列过滤

转录组测序数据在采集的样本提纯后，仍然会包含大量的rRNA和tRNA，肠道微生物 的转录组数据中可能还会含有来自人类的rRNA和tRNA污染序列。RNA-QC-Chain中的RQC-rRNA-filter原理是用HMMER的隐马尔科夫模型(Hidden Markov Model,HMM)搜 索，从输入序列中预测出rRNA片段，然后把它们提取出来。rRNA模式模型是通过来自 SILVA数据库(123版本)的16S/18S/23S/28S rRNA片段建立的。由于HMM算法不依赖 于rRNA的来源基因组的数据，而是依赖于rRNA序列的模式，因此RNA-QC-Chain使得 rRNA片段的去除不依赖于比对和注释。该步骤最后得到高质量以及内源污染序列被过滤的 预处理数据，统计各组数据中的目标读段和污染读段情况。

该步骤完成后，三组人为污染数据需要继续进行后续处理，其余三组真实数据则质量控 制处理完成(Clean real data，RDC)，可用于差异基因表达分析。

三、外源污染序列识别去除

1.测序数据外源污染序列识别

首先要进行的是动植物DNA条形码序列的多重序列比对步骤，因为下一步要建立的 HMM数据库需要该多重序列比对文件作为输入。动物的barcode gene选用细胞色素c氧化 酶I基因(COI)，植物的barcode gene选用基因ITS2，分别得到多重序列比对文件后，使用HMMER中的hmmbuild进行HMM数据库的建立，并且将这两个HMM格式文件合并 到一起，便于序列在其中搜索。

然后对预处理数据进行从头组装，目的是为了得到较长的contig或者scaffold序列，用 于在HMM数据库中比对搜索。该步骤需要组装动植物污染源的读段，动植物中的基因存 在可变剪切，而Trinity对于可变剪切问题处理较好，因此它比较适合用于该步骤的从头组 装。接着使用HMMER中的hmmsearch，输入HMM格式文件和组装完成的长片段序列的FASTA格式文件，搜索得到结果。通过结果文件中的信息可以得知在数据库中比对上的contig/scaffold序列的匹配度、比对上的基因的信息等，由此得知预处理数据中外源污染情 况。

2.测序数据外源污染序列去除

首先使用Bowtie2-build建立污染源参考基因组的索引，这一步和一般的比对过程中建 立参考基因组的索引没有差别。然后使用KneadData，它可以同时输入多个参考基因组的索 引文件，并且比较快速的得到结果，它的核心比对方法也是Bowtie。

三、结果评估

1.在测序数据质量修正方面

如图3和图4所示，六组数据中的读段均去除其5端的12个碱基，然后以5个碱基为单位窗口设置阈值为20进行滑动，去除窗口中碱基平均质量值低于阈值的读段，最后控制读段的长度不小于85个碱基。质量修正处理过后，读段的质量得到了大幅度提升。该处理方式既解决了读段5端碱基比例波动问题，又对读段的质量进行了修正。

2.在外源污染序列识别去除方面

表3-2所示为HQADC数据集的转录本从头组装出的contig的统计表，发现三组数据的N50长度有差异，最长的为1798碱基，N50长度已经可以达到覆盖DNA条形码序列的 目的，且每组数据中有80％左右为污染序列，所以组装出来的contig也会有很大几率属于 污染源生物，但是具体多大概率不得而知，目前能确定的是该方法有一定鉴别能力；三组 数据的contig长度中值均不太理想，片段过短，如果比对打分较高的contig序列长度太短， 依然不能作为污染识别的依据。

组装完成的contig序列在建立好的HMM搜索库中进行比对搜索(如表3-3所示)，可以识别出某些contig来源于人类(homo sapiens)或者大麦(oryza barthii)，但是同样会有 其它比对打分较高的contig，这类contig并不来自于理想的污染源，如在HQACD769439中识别出双形蓟马(Stenchaetothrips biformis)，这种误差肯定存在，但是不影响整个质量控 制流程。因此可以将误差污染源作为一种真实存在的污染，下载其参考基因组数据，同样 将读段去做比对。

表3-2转录本从头组装contig信息统计表

表3-3污染序列识别情况

3.在效率方面

对质量控制流程中每一步处理得到的数据进行了组分(目标读段和污染读段)的统计(如 表3-4所示)，并从数据组成层面上比较人为污染数据和真实数据在每一步的差异情况，以 此来评估质量控制流程在每一步的可信度。从表中可以得知，ACD769409、ACD769424、 ACD769439三组数据中的污染序列的占比均大幅下降，分别为从85.6％降至10.9％，从 82.5％降至9.7％，从93.4％降至41.5％，该数据说明了测序数据经过设计的质量控制流程的 处理，可以在一定程度上降低数据的污染率。

表3-4各步骤数据组分统计表

各组数据每步处理后，数据中污染占比的具体信息如图5所示，三组数据中均发现在污 染序列识别与去除步骤中，污染序列占比被降低最显著，说明该步骤在整个质量控制流程 中达到了污染去除的目的。图5中横坐标表示人为污染数据集ACD中三组数据样本，纵坐 标表示污染读段数占总读段数的百分比，不同颜色的柱子表示各步骤处理后的数据。从图5 可以看出经过处理后，污染读段显著减少。

最后有必要在此说明的是：以上实施例只用于对本发明的技术方案作进一步详细地说 明，不能理解为对本发明保护范围的限制，本领域的技术人员根据本发明的上述内容作出 的一些非本质的改进和调整均属于本发明的保护范围。

Claims (10)

1.一种高通量测序的RNA数据处理方法，其特征在于：所述数据处理方法以高通量测序的RNA读段建立数据库，去除读段中的接头引物和低质量碱基，多重序列对比识别内源污染序列和外源污染序列并清除，获得处理后的RNA数据。

2.根据权利要求1所述的高通量测序的RNA数据处理方法，其特征在于，所述内源污染序列包括测序数据中的：16S rRNA序列、18S rRNA序列、23S rRNA序列和28S rRNA序列。

3.根据权利要求1所述的高通量测序的RNA数据处理方法，其特征在于，所述内源污染序列包括测序数据中来源于动植物的RNA序列。

4.根据权利要求1所述的高通量测序的RNA数据处理方法，其特征在于，所述数据处理方法包括如下步骤：

a)获得高通量测序的RNA-seq测序读段，根据取样环境确定可能的污染源基因组，并下载生物DNA条形码；

b)查看步骤a)中的高通量测序读段是否有接头序列，评估读段各位置的碱基质量，

并评估读段5端和3端的质量；

c)处理低质量碱基和碱基比例波动，过滤内源性污染物；

d)采用MUSCLE、Trinity和HMMER识别外源污染序列；

e)采用KneadData、SAMtools和Bowtie2去除外源污染序列。

5.根据权利要求1所述的高通量测序的RNA数据处理方法，其特征在于：所述真实数据的保存格式为FASTQ。

6.根据权利要求1所述的高通量测序的RNA数据处理方法，其特征在于：所述动物DNA条形码从BOLD数据库中下载，植物物种的DNA条形码序列集采用宁康教授实验室数据。

7.根据权利要求1所述的高通量测序的RNA数据处理方法，其特征在于：数据评估步骤采用FastQC进行。

8.根据权利要求5所述的高通量测序的RNA数据处理方法，其特征在于，数据评估包括：每个位置的碱基质量情况、读段每个位置的各碱基频率分布图、读段的GC含量频率分布图和接头含量。

9.根据权利要求1所述的高通量测序的RNA数据处理方法，其特征在于：步骤c)具体为：采用Trimmomatic处理低质量碱基和碱基比例波动，采用RNA-QC-chain过滤内源污染序列。

10.根据权利要求7所述的高通量测序的RNA数据处理方法，其特征在于：RNA-QC-chain采用隐马尔科夫模型过滤内源污染序列。