CN111584006B - 基于机器学习策略的环形rna识别方法 - Google Patents

Abstract

本发明公开了一种基于机器学习策略的环形RNA识别方法，输入数据，在参考基因组上对每个候选环形RNA进行定位并提取这些环形RNA区域附近的Reads特征；使用提取好的特征训练有监督的机器学习模型；使用训练好的模型对候选环形RNA集合进行真假阳性分类，并输出最终的环形RNA。本发明属于一类机器学习过滤策略，具备机器学习过滤策略的优势，在临床实践中可以显著节省成本、时间等。

Description

基于机器学习策略的环形RNA识别方法

技术领域

本发明属于数据科学技术领域，具体涉及一种基于机器学习策略的环形RNA 识别方法。

背景技术

环形RNA(英文名称：Circular RNA，英文缩写：CircRNA)是非编码RNA (英文名称：non-coding RNA，英文缩写：ncRNA)家族中的重要成员。环形 RNA的定义是：circRNA(circular RNA，环状RNA)是一类具有闭合环状结构的非编码RNA分子，没有5′帽子结构和3′poly(A)结构。它的存在早在上世纪 70年代就被发现了，但由于当时技术和知识水平的限制使得环形RNA一度被认为是剪接错误或转录噪音的结果。近年来随着研究的深入和测序技术的发展，从 2012年首次通过RNA测序(英文名称：RNA sequencing，英文缩写：RNA-seq) 技术发现了大量环形RNA的存在并证实环形RNA是区别于传统线性RNA的一类新型RNA分子，通过上游5’剪接位点共价连接下游3’剪接位点形成闭合环状结构，大量且稳定的存在于真核转录组中，并显示出包含一个或多个外显子，以及内含子、基因间区和UTR区的复杂结构。随着RNA-seq技术的发展，越来越多的环形RNA被鉴定和验证，并且大量研究表明并强调了它们参与多种生理过程，并与多种复杂疾病相关，可作为疾病的生物标注物，例如消化系统癌症、乳腺癌和糖尿病等。其中部分环形RNA通过抑制特定miRNA(英文名称：microRNA，英文缩写：miRNA)的活性，来上调或下调特定基因的表达。最近，许多研究表明环形RNA具有翻译潜力，这意味着一些环形RNA可能在体内被翻译为多肽或蛋白质。因此环形RNA所隐藏的巨大潜力都需要进行深入的研究，这将对疾病的发生，发展和精确诊断和治疗具有深远的意义。

检测环形RNA是任何相关研究的基础且关键的步骤。随着RNA-seq技术的普及，利用计算方法从高通量测序数据中对环形RNA进行识别和分析成为一种广泛使用的方式。目前已提出了多种从RNA-seq数据中识别环形RNA的计算方法，然而这些方法普遍存在识别灵敏度和准确率不高的问题。目标是尽可能减少假阳性和假阴性结果，提高环形RNA数据的精度，以求降低数据噪音对下游研究的干扰，加快环形RNA的研究进程。为此，所有针对环形RNA识别的计算方法在得到检测结果后都需要进一步验证数据的可靠性，正确分类其中的真假阳性环形RNA，以期获得高精度的环形RNA数据。

发明内容

本发明所要解决的技术问题在于针对上述现有技术中的不足，提供一种基于机器学习策略的环形RNA识别方法，通过对候选环形RNA潜在特征进行提取并利用机器学习模型进行训练，二分类结果中真类的结果是最后输出的环形RNA，而假类的结果是要被过滤掉的假阳性环形RNA。

本发明采用以下技术方案：

一种基于机器学习策略的环形RNA识别方法，包括以下步骤：

S1、输入SAM格式和FASTA格式的数据，运行环形RNA识别算法获取候选环形RNA集合作为训练样本集并确定候选环形RNA的断点位置；候选环形 RNA包括环形RNA的左断点brk₁和右断点brk₂，将候选环形RNA集合中的环形 RNA按照左断点brk₁升序进行排序，每个环形RNA用参考基因组序列上的坐标 [brk₁,brk₂]表示；

S2、在参考基因组上对每个候选环形RNA进行定位并提取环形RNA区域附近的Reads特征；

S3、选择监督学习算法构造初始分类器，使用步骤S2提取好的特征训练有监督的机器学习模型；

S4、使用训练好的模型对需要进行预测的候选环形RNA集合进行真假阳性分类并输出，完成RNA识别。

具体的，步骤S1具体为：

S101、读取现有检测方法的输出结果作为训练样本集，给定候选环形RNA 所在参考基因组号，上下游断点位置以及所检测到上下游断点处的BSJ reads信息；

S102、分别读取经人工标记的初始训练样本集L和经人工标记的初始测试集 Lt，初始训练集L的样本量大于等于100。

具体的，步骤S2具体为：

S201、根据环形RNA上下游两处剪接位点的discord和concord的数量得到discord_brk₁、discord_brk₂、concord_brk₁和concord_brk₂特征；

S202、从给定的SAM文件中分别按照给定的环形RNA的两个断点位置处 reads的CIGAR值提取SM_l，SMS_l，SMS_l，SM_r，MS_r，SMS_r特征

S203、依据BSJ特征分别对环形RNA左断点和右断点的支持数和比对质量进行统计，得出Support_l、Mapping_Quality_l、Support_r、Mapping_Quality_r特征；

S204、分别对环形RNA上两个断点的上游和下游区域的平均测序深度和位点覆盖度进行统计，得到cov_brk1_Up，cov_brk1_Down，cov_brk2_Up， cov_brk2_Down，depth_brk1_Up，depth_brk1_Down，depth_brk2_Up，depth_brk2_Down特征；

S205、从给定的FASTA格式文件中根据给定的环形RNA上下游断点位置提取GTAG特征。

进一步的，步骤S203中，断点处的平均比对质量计算如下：

其中，m_i是BSJ read的比对质量，l_d是支持数。

进一步的，步骤S204中，区域的平均测序深度计算如下：

其中，l_d是区域的长度，d_i是映射在位置i处的测序深度；

每个碱基测序深度的总和计算如下：

其中，c_i是映射在位置i处的碱基覆盖度。

具体的，步骤S3中，对每个候选环形RNA提取步骤S2中所列出的23个特征，使得每一个候选环形RNA得到23维的特征向量，把每个候选环形RNA的 23个特征向量存为一行表示对应的候选环形RNA，使用GBDT对训练样本集中的这些特征进行训练，GBDT使用默认的随机树核函数训练模型，并使用网格搜索和10倍交叉验证查找最优参数设置。

具体的，步骤S4中，输出数据为候选环形RNA集合标注真假阳性的文件，输出数据的格式为每一行包含候选环形RNA信息和真假阳性的标注文件，输出数据为候选环形RNA集合标注真假阳性的文件，输出数据的格式为每一行包含候选环形RNA信息和真假阳性的标注文件。

与现有技术相比，本发明至少具有以下有益效果：

本发明一种针对靶向捕获转录组测序数据的假阳性环形RNA过滤方法，基于环形RNA检测软件输出的环形RNA的检测报告文件，针对其中存在假阳性错误，引用有监督学习模型；模型提取候选环形RNA的数据特征，对特征进行二分类学习训练(区分确实存在的环形RNA和假阳性错误)，当训练完成后，使用经过训练的模型将检测报告文件中的假阳性错误过滤掉。本发明利用了有监督学习模型的优势，解决了强制过滤所造成的灵敏度降低，并进一步提高了检测精确度。

进一步的，步骤S1中运行现有环形RNA检测方法获取候选集，选取大量样本进行人工标记和验证是使用有监督机器学习模型的前提条件，直接影响模型的精度和泛化性等主要性能指标。设置初始训练集L的样本量大小不低于100，有监督学习模型需要输入适当大小的训练集，训练集的大小直接影响模型的精度和泛化性等主要性能指标；根据理论推导和数据测试获得的规模的一个下界，当规模小于本设置时，模型的精度等将出现显著下降。

进一步的，特征提取中所需要的文件包括FASTA格式的参考基因组文件以及已有检测方法所依赖的SAM文件，SAM是一种序列比对格式的标准，主要应用于测序序列mapping到基因组上的结果表示，SAM格式中的多列信息是识别环形RNA的重要依据；FASTA文件存储的碱基序列可用来分析剪切位点周围是否存在GT-AG信号来判断是否为潜在的环形RNA。SAM文件和FASTA文件是目前环形RNA识别算法的常用输入文件。

进一步的，在仿真实验中，可用环形RNA仿真器生成的标注文件对环形RNA 识别算法检测到的候选环形RNA进行样本标注，线性RNA的label标注为0，环形RNA的label标注为1。另外，对于真实数据实验，在一对没有经过RNase R 处理的RNA-seq数据集和经过RNaseR处理的RNA-seq数据集中，如果环形RNA 识别算法在没有经过RNase R处理的RNA-seq数据中识别到的环形RNA对比在经过RNase R处理的数据集中呈现出明显的富集，则标注为环形RNA，label值为1；若没有呈现出富集则标注为线性RNA，label值为0。

进一步的，使用大量的已标记数据进行机器学习，为了达到设计精度，必须进行反复学习，才能降低分类错误，获得符合性能要求的模型，继而可以使用模型对假阳性错误进行过滤。

进一步的，在选择机器学习模型中尽量选择使用支持非线性数据，并对异常值不敏感的，可以增加分类效果。

综上所述，本发明属于一类机器学习过滤策略，设计和使用了一种有监督学习模型，具备机器学习过滤策略的优势。同时，针对现有环形RNA检测方法在过滤环节中的不足——其一，除了提取现有环形RNA通常使用的过滤特征如BSJ 的属性，还另外对测序数据的一致性、完整性、测序深度等特征进行了提取；其二，使用了机器学习有效的将这些特征进行组合分析，这是环形RNA检测中首次使用机器学习技术来识别环形RNA。实验数据证明，本发明模型性能显著优于已有方法，在临床实践中可以显著节省成本、时间等。

下面通过附图和实施例，对本发明的技术方案做进一步的详细描述。

附图说明

图1为本发明流程示意图；

图2为不同测序覆盖下识别灵敏度示意图；

图3为不同测序覆盖下识别精确度示意图；

图4为不同测序覆盖下识别F1值示意图；

图5为CIRCPlus2和CIRI2的在识别环形RNA上的灵敏度、精确度和F1 值示意图。

具体实施方式

请参阅图1，本发明一种基于机器学习策略的环形RNA识别方法，包括以下步骤：

S1、输入数据

RNA-seq生成的SAM文件和基因注释FASTA文件，输入数据的格式要求是： SAM格式和FASTA格式；运行现有环形RNA识别算法来获取候选环形RNA集合，并确定候选环形RNA的断点位置；

运行现有的环形RNA识别算法输出候选环形RNA集合，这里用参考基因组号，环形RNA的5’剪接位点的位置定义为左断点，用brk₁表示；环形RNA的3’剪接位点的位置定义为右断点，用brk₂表示；两者可以共同表示一个候选环形 RNA；将候选环形RNA集合中的环形RNA按照左断点升序进行排序；此时，每个环形RNA可用参考基因组序列上的一对坐标来表示，即[brk₁,brk₂]，该坐标对分别表示环形RNA的上游和下游成环剪接位点的位置。

S101、读取现有检测方法的输出结果，其中需要给定环形RNA的所在参考基因组号，上下游断点位置以及所检测到上下游断点处的BSJ reads的具体信息 (每条read名称，比对质量)；

S102、分别读取经人工标记的初始训练样本集L和经人工标记的初始测试集 Lt，可以取消现有环形RNA检测算法中的过滤步骤允许增加训练集中的负样本。

初始训练集L的样本量大小不低于100。

训练样本集和测试集中的样本标注方法如下：

在仿真实验中，可用环形RNA仿真器生成的标注文件对环形RNA识别算法检测到的候选环形RNA进行样本标注，线性RNA的label标注为0，环形RNA 的label标注为1；

对于真实数据实验，在一对没有经过RNase R处理的RNA-seq数据集和经过RNaseR处理的RNA-seq数据集中，如果环形RNA识别算法在没有经过RNase R处理的RNA-seq数据中识别到的环形RNA对比在经过RNase R处理的数据集中呈现出明显的富集，则标注为环形RNA，label值为1；若没有呈现出富集则标注为线性RNA，label值为0。

S2、在参考基因组上对每个候选环形RNA进行定位并提取这些环形RNA区域附近的Reads特征；

每个环形RNA在参考基因组有上下游两个断点，环形RNA断点周围的Reads 与其它区域的Reads映射到参考基因组上会呈现出不同的特点。其中Reads映射后的一致性、完整性、测序深度这三种属性所表现出的潜在特征可以用于揭示环形RNA的存在。针对这三种属性所提取的特征如下：

S201、按照一致性提取特征

对于一个给定的环形RNA，用[brk₁,brk₂]来表示，其中brk₁,brk₂分别为该环形RNA的左断点和右断点。针对环形RNA的范围，算法分别向左右各扩展一个 range长度，此处的range＝μ+3σ，其中，μ、σ分别是文库的均值和标准差。收集区域[brk₁-range,brk₁+range]以及[brk₂-range,brk₂+range]之间的Paired-end Reads。如果Paired-end Reads映射在参考基因组后插入距离异常或方向异常可能预示环形RNA的存在，而插入距离正常的Paired-end Reads可能被认为是环形RNA不存在的表征。因此，在这里对于给定的候选物，收集跨越环形RNA断点的所有Paired-end Reads，即[brk₁-range,brk₁+range]以及[brk₂-range,brk₂+range]区域内不一致和一致的读对，并将这两种读对分别记为discord和concord。因此，可根据环形RNA上下游两处剪接位点的discord和concord的数量得到4个特征，discord_brk₁、discord_brk₂、concord_brk₁和concord_brk₂。

具体的，需要从给定的SAM文件中分别对环形RNA的两个断点位置按以上要求提取4个特征值。

S202、根据完整性提取特征

完整性描述了Paired-end Reads在参考基因组上比对的具体情况，是环形 RNA识别的重要标志。完整性主要分为这三种比对情况：

1)完全映射(Mapped Reads)，完全映射意味着Paired-end Reads没有发生剪切，可以整体直接映射在参考基因组上，表明了在映射位置处不存在环形RNA；

2)部分映射(Split Reads)，部分映射意味着Reads不能完整的映射在参考基因组上，需要拆分后才能比对上，这一特性支持了环形RNA的存在，剪接映射的读段通常有助于找到环形RNA的确切断点；

3)未映射(Unmapped Reads)，指的是无论怎样拆分Reads都无法在参考基因组上找到合适的映射位置，在环形RNA的识别中可忽略。

因此，在环形RNA的识别中只考虑Split Reads，具体指映射在环形RNA断点周围的Reads。

完整性中只有部分映射是环形RNA存在的标志。部分映射比对到参考基因组上的映射形式不同，在环形RNA识别过程中起到的作用大小也不同。部分映射是指Reads无法完整比对在参考基因组上，需要拆分成多段来进行映射，即 SAM文件中的CIGAR值会出现除了‘M’类型以外的其它映射状态，包含‘I’、‘D’、‘N’、‘P’、‘S’、‘H’这些映射形式。目前在环形RNA的识别中，只有‘S’和‘H’可以有效表征剪接事件，而其它的映射状态并无正面作用可以被忽略。

Split Reads的映射状态只考虑‘S’，‘H’和‘M’的话，其中‘H’只出现在一条Read的前端或末端，但不会出现在中间，因此这三种类型可能会组合出SM，HM， MS，MH，SMS，HMH，HMS，SMH这8种CIGAR值的组合形式。具体来说， Split Reads会以两种形式比对在参考基因组上。

一方面，典型的剪接会呈现两段式映射的结果，位于环形RNA的5’端剪接位点处的BSJ reads的CIGAR值形式为xS/HyM，相应的在3’剪接位点处会显示 xMyS/H的形式；

另一方面，一些环形RNA具有更复杂的比对特征，如长度较短的环形RNA 以及侧接短外显子的环形RNA，这些环形RNA的BSJ reads可能以三段式不连续地映射在参考基因组上，从而在3’和5’剪接位点处显示出xS/HyMzS/H的比对结果。

由于‘S’和‘H’所表征的都是剪接形式，通常会成对出现，因此将它们归类为同一种形式，因此可将上述组合后的8种CIGAR类型总结为三种映射状态SM， MS和SMS。依据完整性可对环形RNA左右断点附近的所有的Split Reads进行提取，左断点记为l，右断点记为r，分别对两个断点计算出具有相应CIGAR类型的Split Reads的数量，可得到6个特征：SM_l，MS_l，SMS_l，SM_r，MS_r， SMS_r。

具体来说，从给定的SAM文件中分别按照给定的环形RNA的两个断点位置处reads的CIGAR值提取出6个特征。

S203、依据BSJ提取特征

支持数(Support)、比对质量(Mapping_Quality)都有助于区分候选物，这些都是基于环形RNA识别算法所识别出的BSJ的特征。在已知环形RNA断点处的测序深度的条件下，可以根据检测到的BSJ reads的数量即支持数来推定环形RNA的表达丰度，从而判断环形RNA的存在或不存在。需要指出的是，环形RNA的左断点和右断点应分开处理，根据两个断点间的表达丰度的差异有助于识别具有不准确断点位置的假阳性环形RNA。

另外，在高通量测序数据中，比对质量可以衡量Read在参考基因组上的比对情况，分值越高说明该Read在参考基因组上的比对位置越唯一，出现错误映射或多处映射的概率更小。尤其对于BSJ reads这样的部分映射读段，因为参考基因组上重复区域的存在，短片段在参考基因组上准确映射的难度更大，如果片段长度小于19bp还会被比对软件忽略，片段的映射结果越可靠对于环形RNA的判定更有利。因此左右断点处的BSJ reads均需要达到一定的比对质量来保证环形RNA的可靠性。

下面的公式用来计算一个断点处的平均比对质量，具体为：

其中，m_i是BSJ read的比对质量，l_d是支持数。

依据BSJ特征分别对环形RNA左断点和右断点的支持数和比对质量进行统计，可得出以下4个特征：Support_l、Mapping_Quality_l、Support_r、 Mapping_Quality_r。

具体来说，这里根据现有检测算法所输出的BSJ reads所给定的信息给出4 个特征。

S204、依据测序深度提取特征

测序深度(Read Depth)是测序数据中广泛使用的特征，其是指映射到特定位点或基因组区域的Reads的数量。环形RNA通常以不同的方式增加或降低Read Depth，因此在参考基因组上具有与预期不同的平均映射读段深度的区域可以支持或反对环形RNA的存在。

用下面的公式来计算一个区域的平均测序深度，具体为：

其中，l_d是区域的长度，d_i是映射在位置i处的测序深度。每个位置或区域的测序深度值d_i可以由SAMtools计算出来。

此外，还使用SAMtools计算了在提供的BED文件中指定的每个基因组区域的总碱基数量(即每个碱基测序深度的总和)，公式中l_d是区域的长度，c_i是映射在位置i处的碱基覆盖度。

这里分别对环形RNA上两个断点的上游和下游区域的平均测序深度和位点覆盖度(bed coverage)进行统计，可以得到8个特征：cov_brk1_Up， cov_brk1_Down，cov_brk2_Up，cov_brk2_Down，depth_brk1_Up， depth_brk1_Down，depth_brk2_Up，depth_brk2_Down。

S205、依据GT-AG信号提取特征

根据GT-AG信号对候选环形RNA提取GTAG特征，GT-AG信号是真核转录物中的主要剪切信号，在多个环形RNA的识别算法中用于最后验证BSJ reads 的可靠性。GT-AG规则指的是mRNA前体(pre-mRNA)在参与内含子剪接的两个特殊位点，每个内含子的5’端开始都邻接着两个核苷酸GT，3’端末尾都邻接着两个核苷酸AG。通过分析Junction Read的两个区段在参考基因组上的映射情况，检查它们是否侧接了AG和GT双核苷酸(或反向互补的双核苷酸CT和 AC)来推断BSJ reads是否可靠。

根据候选环形RNA指定的上下游两处剪接位点[brk₁,brk₂]在参考基因组上进行定位，检查两处剪接位点旁是否侧接了AG和GT双核苷酸(或反向互补的双核苷酸CT和AC)，对于边界模糊而使AG和GT双核苷酸没有紧邻剪接位点的情况，如果AG和GT双核苷酸在参考基因组上与给定的剪接位点偏离的距离和方向一致的话，也可认定满足GT-AG信号，因此这里可得到一个特征GTAG。

具体来说，需要从给定的FASTA格式的文件(参考基因信息)中根据给定的环形RNA上下游断点位置提取GTAG特征。

综上所述，对于每个环形RNA，共提取出了23个特征。随着对环形RNA 研究的深入，可能会发现更多可有效表征环形RNA的特征，均可添加到本发明创建的特征向量中并用于机器学习训练，即特征提取的这一步骤是可扩展的。

S3、使用提取好的特征训练有监督的机器学习模型；

选择监督学习算法构造初始分类器，即选用1个监督学习自动机并基于初始训练集生成的学习器；

在获得候选环形RNA集合后，接着对每个候选环形RNA提取步骤S2中所列出的23个特征，使得每一个候选环形RNA得到23维的特征向量。这里用参考基因组号，环形RNA的5’剪接位点的位置(左断点，brk₁)，环形RNA的3’剪接位点的位置(右断点，brk₂)来表示一个候选环形RNA，并将候选环形RNA 集合中的环形RNA按照左断点升序进行排序。接着依次对候选环形RNA集合中的每个环形RNA提取特征，对每一个候选环形RNA而言，先用[brk₁,brk₂]在参考基因组上定位，然后收集指定区域的Reads，进行特征值的计算。把每个候选环形RNA的23个特征向量存为一行，用来表示其对应的候选环形RNA，换言之，每一行可以看作一个实例。这些实例都对应着一个标签值(label)，label是类别号，取值为0，1，其中1表示真阳性环形RNA类，0表示假阳性环形RNA类。

对训练集中每个候选环形RNA进行特征提取后，可以使用通过GBDT对训练样本集中的这些特征进行训练，也可换用其它预测性能更优的机器学习模型来训练。GBDT使用默认的随机树(Random Trees)核函数来训练模型，并使用网格搜索和10倍交叉验证来查找最优参数设置。

S4、使用训练好的模型对候选环形RNA集合进行真假阳性分类，并输出最

输出数据是：给候选环形RNA集合标注真假阳性的文件，输出数据的格式要求是：每一行包含候选环形RNA信息和真假阳性的标注文件。

使用训练好的模型来对测试数据集进行真假阳性环形RNA分类，被归类为真阳性的环形RNA作为最终输出结果。

步骤S1、S2、S3中，选择使用其它环形RNA识别算法获取候选环形RNA 集合，添加其它潜在的可表征环形RNA的特征进行训练，也可选择使用其它机器学习模型进行最终的分类。

为使本发明实施例的目的、技术方案和优点更加清楚，下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例是本发明一部分实施例，而不是全部的实施例。通常在此处附图中的描述和所示的本发明实施例的组件可以通过各种不同的配置来布置和设计。因此，以下对在附图中提供的本发明的实施例的详细描述并非旨在限制要求保护的本发明的范围，而是仅仅表示本发明的选定实施例。基于本发明中的实施例，本领域普通技术人员在没有作出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。

(1)仿真实验

首先以hg19参考基因组中的1号染色体作为参考基因组。选用CIRCPlus算法输出候选环形RNA测试集，采用了灵敏度(图2)和精确度(图3)指标来进行算法评估，考虑到精确度和灵敏度相互制约影响的事实，为了评估灵敏度和精确度的平衡性能，还采用了F1值(图4)，衡量机器学习过滤算法对测试集中真假阳性的分类结果。其中CIRCPlus_initial为CIRCPlus算法取消过滤步骤所得到的输出结果，CIRCPlus2为所提出的专利方法的输出结果。可以看出，CIRCPlus2 在不降低灵敏度的同时，可以将现有算法的输出结果中的大量假阳性进行过滤，并将检测精确度提高到98％左右，并显示出更优的F1值。

(2)真实实验

在真实数据中，评估CIRCPlus2与CIRI2在真实数据中识别性能上的差异，将这两种工具应用于HEK293细胞系的一对没有经过RNase R处理的RNA-seq 数据集和经过RNaseR处理的RNA-seq数据集中，算法在没有经过RNase R处理的RNA-seq数据中识别到的环形RNA需要根据是否在经过RNase R处理的数据集中呈现出明显的富集，来判定识别是否正确，并用来计算各自算法的精确度，如果呈现出明显的富集即为真阳性环形RNA则在CIRCPlus2训练样本集的标签值为1，若没有呈现出富集即为假阳性环形RNA则在CIRCPlus2训练样本集的标签值为0。两个算法最终识别到的环形RNA作为总环形RNA来计算每个算法的识别灵敏度。

CIRCPlus2具体使用HEK293数据集中的第1到22号染色体的测序数据来识别候选环形 RNA(不考虑性染色体)，并且将CIRCPlus2初步识别到候选环形RNA的70％(2943中的2100 个)用于训练模型，其余的30％的候选物用于评估算法的性能。CIRCPlus2在使用分类模型后总共预测出的真阳性环形RNA为1465个，其中经过RNase R抗性验证后的真阳性环形RNA 为1119个。将CIRCPlus2和CIRI2识别出的所有真阳性环形RNA作为总环形RNA(共计2013 个环形RNA)。图5展示了CIRCPlus2和CIRI2的在识别环形RNA上的灵敏度、精确度和F1值，可以看出CIRI2和CIRCPlus2的灵敏度几乎相同，但CIRCPlus2在精确度和F1得分上高于CIRI2，证明了分类模型在真实数据中也可表现出较优的识别性能。

综上所述，本发明一种基于机器学习策略的环形RNA识别方法，CIRCPlus2 首先收集现有环形RNA检测算法输出的候选环形RNA集合作为训练样本集。接下来，对候选环形RNA样本集中的每个候选环形RNA提取序列特征，包括候选环形RNA剪接位点处序列的一致性、完整性、测序深度这三种属性特征，以及依据BSJ属性和GT-AG信号所提取的特征。然后，使用机器学习方法对所提取的特征向量和带有是否为环形RNA还是线性RNA的标注值进行有监督的模型训练。最后，CIRCPlus2通过训练好的模型可实现对需要进行预测的候选环形RNA 集合进行真假阳性环形RNA的分类。

以上内容仅为说明本发明的技术思想，不能以此限定本发明的保护范围，凡是按照本发明提出的技术思想，在技术方案基础上所做的任何改动，均落入本发明权利要求书的保护范围之内。

Claims (6)

1.一种基于机器学习策略的环形RNA识别方法，其特征在于，包括以下步骤：

S1、输入SAM格式和FASTA格式的数据，运行环形RNA识别算法获取候选环形RNA集合作为训练样本集并确定候选环形RNA的断点位置；候选环形RNA包括环形RNA的左断点brk₁和右断点brk₂，将候选环形RNA集合中的环形RNA按照左断点brk₁升序进行排序，每个环形RNA用参考基因组序列上的坐标[brk₁,brk₂]表示；

S2、在参考基因组上对每个候选环形RNA进行定位并提取环形RNA区域附近的Reads特征，具体为：

S201、根据环形RNA上下游两处剪接位点的discord和concord的数量得到discord_brk₁、discord_brk₂、concord_brk₁和concord_brk₂特征；

S202、从给定的SAM文件中分别按照给定的环形RNA的两个断点位置处reads的CIGAR值提取SM_l，MS_l，SMS_l，SM_r，MS_r，SMS_r特征；

S203、依据BSJ特征分别对环形RNA左断点和右断点的支持数和比对质量进行统计，得出Support_l、Mapping_Quality_l、Support_r、Mapping_Quality_r特征；

S204、分别对环形RNA上两个断点的上游和下游区域的平均测序深度和位点覆盖度进行统计，得到cov_brk1_Up，cov_brk1_Down，cov_brk2_Up，cov_brk2_Down，depth_brk1_Up，depth_brk1_Down，depth_brk2_Up，depth_brk2_Down特征；

S205、从给定的FASTA格式文件中根据给定的环形RNA上下游断点位置提取GTAG特征；

S3、选择监督学习算法构造初始分类器，使用步骤S2提取好的特征训练有监督的机器学习模型；

S4、使用训练好的模型对需要进行预测的候选环形RNA集合进行真假阳性分类并输出，完成RNA识别。

2.根据权利要求1所述的基于机器学习策略的环形RNA识别方法，其特征在于，步骤S1具体为：

S101、读取现有检测方法的输出结果作为训练样本集，给定候选环形RNA所在参考基因组号，上下游断点位置以及所检测到上下游断点处的BSJ reads信息；

S102、分别读取经人工标记的初始训练样本集L和经人工标记的初始测试集Lt，初始训练集L的样本量大于等于100。

3.根据权利要求1所述的基于机器学习策略的环形RNA识别方法，其特征在于，步骤S203中，断点处的平均比对质量计算如下：

其中，m_i是BSJ read的比对质量，l_d是支持数。

4.根据权利要求1所述的基于机器学习策略的环形RNA识别方法，其特征在于，步骤S204中，区域的平均测序深度计算如下：

其中，l_d是区域的长度，d_i是映射在位置i处的测序深度；

每个碱基测序深度的总和计算如下：

其中，c_i是映射在位置i处的碱基覆盖度。

5.根据权利要求1所述的基于机器学习策略的环形RNA识别方法，其特征在于，步骤S3中，对每个候选环形RNA提取步骤S2中所列出的23个特征，使得每一个候选环形RNA得到23维的特征向量，把每个候选环形RNA的23个特征向量存为一行表示对应的候选环形RNA，使用GBDT对训练样本集中的这些特征进行训练，GBDT使用默认的随机树核函数训练模型，并使用网格搜索和10倍交叉验证查找最优参数设置。

6.根据权利要求1所述的基于机器学习策略的环形RNA识别方法，其特征在于，步骤S4中，输出数据为候选环形RNA集合标注真假阳性的文件，输出数据的格式为每一行包含候选环形RNA信息和真假阳性的标注文件。

Images