SPOT-RNA:使用二维深度神经网络和迁移学习的集合进行 RNA 二级结构预测

摘要：我们的大部分人类基因组都会转录成结构和功能未知的非编码 RNA。获得非编码 RNA 的功能线索需要准确的碱基配对或二级结构预测。然而，当前基于折叠的算法的此类预测性能十多年来一直停滞不前。在这里，我们建议使用深度上下文学习进行碱基对预测，包括那些通过三级相互作用稳定的非规范和非嵌套（伪结）碱基对。由于只有 250 个非冗余的高分辨率 RNA 结构可用于模型训练，因此我们利用最初使用 10,000 个非冗余的最新高质量 bpRNA 数据集训练的模型进行迁移学习通过比较分析获得RNA。由此产生的方法在预测所有碱基对，特别是非规范和非嵌套碱基对方面取得了巨大的、统计上显着的改进。所提出的方法 (SPOT-RNA) 具有免费提供的服务器和独立软件，应该有助于改进 RNA 结构建模、序列比对和功能注释。

1.简介

        RNA 二级结构由其核苷酸序列中通过氢键配对的核苷酸碱基列表表示。这些碱基对的堆叠形成了驱动RNA三维结构折叠的支架 1 。因此，RNA 二级结构的知识对于 RNA 结构建模和理解其功能机制至关重要。因此，已经开发了许多实验方法来通过使用一维或多维探针来推断配对碱基，例如酶、化学物质、突变和与下一代测序相结合的交联技术 2,3 。然而，单碱基对分辨率下的精确碱基配对信息仍然需要通过 X 射线晶体学、核磁共振 (NMR) 或低温电子显微镜确定高分辨率、三维 RNA 结构。 RNAcentral 4 中收集的 1400 万个非编码 RNA 具有实验确定的结构 5 ，非常需要开发准确且具有成本效益的计算方法来直接从序列预测RNA二级结构。

        目前的RNA二级结构预测方法可分为比较序列分析和具有热力学、统计或概率评分方案的折叠算法6。比较序列分析确定同源序列之间保守的碱基对。如果有大量同源序列可用并且这些序列与专家知识手动比对，这些方法是高度准确的。然而，Rfam8 中已知的 RNA 家族只有几千个。因此，RNA二级结构预测最常用的方法是根据适当的评分函数折叠单个RNA序列。在这种方法中，RNA 结构根据最近邻模型被分为子结构，例如环和茎。然后采用动态规划算法从这些子结构中定位全局最小值或概率结构。每个子结构的评分参数可以通过实验获得（例如，RNAfold11、RNAstruct12 和 RNAshapes13）或通过机器学习（例如，CONTRAfold14、CentroidFold15 和 ContextFold16）。然而，总体精度（所有预测碱基对中正确预测的碱基对的比例）似乎已达到“性能上限”6，约为 80%17,18。这部分是因为所有现有方法都忽略了三级相互作用产生的部分或全部碱基对19。这些碱基对包括单独（非堆叠）、假结（非嵌套）和非规范（非 A-U、G-C 和 G-U）碱基对以及三联体相互作用 19,20。虽然一些方法可以预测带有假结的 RNA 二级结构（例如 pknotsRG21、Probknot22、IPknot23 和 Knotty24），而其他方法可以预测非规范碱基对（例如 MC-Fold25、MCFold-DP26 和 CycleFold27），但它们都不能提供两者的计算预测，更不用说单独的碱基对和碱基三联体了。

       本文提出的工作受到 Raptor-X28 和 SPOT-Contact29 使用残差网络 (ResNets)30 和二维双向长网络等深度学习神经网络算法直接预测蛋白质序列中蛋白质接触图的最新进展的启发。短期记忆细胞（2D-BLSTM）31,32。 SPOT-Contact 将整个蛋白质“图像”视为上下文，并使用 ResNet 超深混合网络与 2D-BLSTM 相结合进行预测。 ResNet 可以从每一层的整个序列“图像”中捕获上下文信息，并映射输入和输出之间的复杂关系。此外，2D-BLSTM 被证明在蛋白质结构预测中传播长程序列依赖性方面非常有效29，因为 LSTM 细胞能够记住训练期间序列位置彼此远离的残基之间的结构关系。与蛋白质接触图类似，RNA 二级结构是二维接触矩阵，尽管其接触定义不同（分别为 RNA 碱基对的氢键和蛋白质接触的距离截止）。然而，与蛋白质不同，蛋白质数据库 (PDB)5 中可用的非冗余 RNA 结构数量较少，使得深度学习方法不适合直接基于单序列的 RNA 二级结构预测。因此，机器学习技术很少得到利用。据我们所知，唯一的例子是 mxfold33，它采用小型机器学习算法（结构化支持向量机）进行 RNA 二级结构预测。与热力学模型相结合后，其性能比基于折叠的技术有了一些改进。然而，mxfold 仅限于规范碱基对，而不考虑假结。

       最近，发布了包含超过 100,000 个 RNA 序列 (bpRNA34) 的大型数据库，该数据库具有二级结构的自动注释。虽然这个数据库足够大，足以让我们采用深度学习技术，但比较分析中带注释的二级结构在单碱基对水平上可能并不可靠。为了克服这个限制，我们首先使用 bpRNA 来训练 ResNet 和 LSTM 网络的集合，类似于我们通过 SPOT-Contact29 用于蛋白质接触图预测的集合。然后，我们使用源自高分辨率 RNA 结构的精确碱基对小型数据库进一步训练大型模型。我们成功地使用这种迁移学习技术来识别蛋白质本质无序区域中的分子识别特征36。由此产生的方法称为 SPOT-RNA，是一种深度学习技术，用于预测所有配对的碱基，无论它们是否与三级相互作用相关。与 62 个高通量独立测试集相比，新方法在非嵌套、非规范和所有碱基对的 F1 分数上分别比次优方法提高了 53%、47% 和 10% 以上。通过 X 射线晶体学解析 RNA 结构。 SPOT-RNA 的性能通过由 NMR 测定的 39 个 RNA 结构和 PDB 中最近发布的 6 个非冗余 RNA 确定的单独测试集进一步证实。

2.结果

2.1 bpRNA初始训练

       我们通过构建一组非冗余的 RNA 序列来训练我们的 ResNets 和 LSTM 网络模型，这些 RNA 序列具有来自 bpRNA34 的注释二级结构，序列同一性截止值为 80%，这是程序 CD-HIT-EST37 允许的最低序列同一性截止值，并且具有之前曾被许多研究用于相同目的38,39。该数据集包含 13,419 个 RNA，在排除那些 >80% 的序列同一性后，进一步随机分为 10,814 个 RNA 用于训练 (TR0)、1300 个用于验证 (VL0) 和 1,305 个用于独立测试 (TS0)。通过使用 TR0 进行训练，VL0 进行验证，以及单个序列（Lx4 的单热向量）作为唯一输入，我们训练了许多具有 ResNet、BLSTM 的数量和大小组合的二维深度学习模型和 FC 层，布局如图 1 所示。最好的 5 个模型（仅通过 VL0 验证）在 VL0 和 TS0 上的集成性能如表 1 所示。与 Matthews 相关系数 (MCC) 的性能基本相同VL0 为 0.632，TS0 为 0.629，表明所训练的集成的稳健性。验证和测试之间的 F1 分数、精度调和平均值和灵敏度也基本相同（0.629 与 0.626）。补充表 1 进一步比较了单个模型与整体模型的性能。在 TS0 中，MCC 从 0.617（最好的单一模型）提高到 0.629，提高了 2%，证实了集成在消除单个模型中的随机预测误差方面的有用性。

图1 SPOT-RNA的广义模型架构。 SPOT-RNA 的网络布局，其中 L 是目标 RNA Act 的序列长度。表示激活函数 Norm。表示标准化函数，PreT表示在bpRNA数据集上训练的预训练（初始训练）模型。 

表 1 初始训练、迁移学习和直接训练后 SPOT-RNA 在验证和测试集上的性能。 

a Matthews 相关系数 b 精度和灵敏度的调和平均值 c 基于五重交叉验证的标准偏差

 

2.2 RNA 结构的迁移学习 

        从 bpRNA 数据集中获得的模型被转移到源自高分辨率非冗余 RNA 结构的碱基对上进行进一步训练，其中 TR1（训练集）、VL1（验证集）和 TS1（测试集）具有 120、30 和 67 个 RNA ， 分别。 TS1 集独立于训练数据（TR0 和 TR1），因为它是通过首先以允许的最低序列同一性截止值 (80%) 通过 CD-HIT-EST 进行过滤而获得的。为了进一步消除潜在的同源性，我们对训练数据（TR0 和 TR1）使用 BLAST-N40，e 值截止值为 10。为了检查所构建模型的一致性，我们通过组合 TR1 和 VL1 数据集进行 5 倍交叉验证。表 1 显示了对相同前 5 个模型的集成的训练数据 (TR1+VL1) 和未见过的 TS1 的交叉验证结果。5 倍 MCC 的微小波动为 0.701 ± 0.02，F1 为 0.690 ± 0.02，5 倍交叉验证和测试集 TS1 之间的微小差异（MCC 为 0.701 vs. 0.690）表明针对看不见的数据。表1还表明，与迁移学习后的模型相比，直接应用bpRNA训练的模型在TS1上的性能合理但较差。在 TS1 上迁移学习之前（0.650）和之后（0.690）MCC 的改进为 6%。补充表 2 和补充表 3 分别比较了五重交叉验证 (TR1+VL1) 和独立测试集 (TS1) 的模型集合和五个单独模型的结果。在交叉验证和独立测试中，观察到该集成相对于最佳单一模型有显着改进，MCC 提高了 3%。

2.3迁移学习和直接学习的比较。

        为了证明迁移学习的实用性，我们还使用相同的集成网络架构和超参数（层数、层深度、内核大小、膨胀因子和学习率）对 5 个模型进行直接训练在结构化 RNA 训练集 (TR1) 上进行，并由 VL1 验证并由 TS1 测试。在 VL1 和 TS1 上直接学习的五个模型的集成性能如表 1 所示。验证和测试之间的相似性能分别为 MCC = 0.583、0.571，证实了直接学习的鲁棒性。然而，这种性能大大低于迁移学习（MCC 值降低 21%，F1 分数降低 30%）。这证实了使用 TR1 的小型训练数据集进行直接学习的难度，以及使用可以有效利用深度学习网络能力的大型数据集（bpRNA）的必要性。补充表 4 通过直接学习 TR1 进一步比较了各个模型与集成模型的性能。图2a比较了独立测试集TS1上初始训练（SPOTRNA-IT）、直接训练（SPOT-RNA-DT）和迁移学习（SPOT-RNA）给出的精确召回（PR）曲线。结果来自减少的 TS1（62 个 RNA，而不是 67 个），因为同一图中显示的一些其他方法不能预测缺失或无效碱基的序列的二级结构。有趣的是，直接训练以非常低的灵敏度（召回率）下的 100% 精度开始，而初始训练和迁移学习在区分正预测和负预测的最高可能阈值下，在可实现的最低灵敏度下具有较高但 <100% 的精度。这表明 bpRNA 中假阳性的存在“污染”了初始训练。尽管如此，迁移学习在 50% 的召回率下实现了 93.2% 的准确率。这表明 bpRNA 中潜在假阳性的比例很小。

图 2 在测试集 TS1 上使用 PR 曲线和箱线图比较 SPOT-RNA 与其他 12 个预测因子的性能。 a 通过初始训练（SPOT-RNA-IT，绿色虚线）、直接训练（SPOT-RNA-DT，蓝色点虚线）和迁移学习（SPOT-RNA-IT）在独立测试集 TS1 上的精确回忆曲线RNA，洋红色实线）。当前使用的十个预测变量的精度和灵敏度结果也显示为用空心符号标记用于考虑伪结的方法，用实心符号标记用于不考虑伪结的方法。 CONTRAfold 和 CentroidFold 也显示为曲线（金色和黑色），因为它们的方法提供了预测概率。 b 通过标记的各种方法给出的独立测试集 TS1 上各个 RNA 的 F1 分数分布。在每个方框上，中心标记表示中位数，方框的底部和顶部边缘分别表示第 25 个和第 75 个百分位数。使用“+”符号单独绘制异常值。

2.4与其他二级结构预测因子的比较。

图 2a 进一步比较了我们的迁移学习集成模型给出的精度/召回率曲线与独立测试集 TS1 上其他 12 个可用的 RNA 二级结构预测因子。具有概率输出的两个预测变量（CONTRAfold 和 CentroidFold）也由 PR 曲线表示，其余的显示为奇点。大多数现有方法的性能集中在 50% 的灵敏度和 67-83% 的精度附近（表 2）。相比之下，与次佳 mxfold 相比，我们的方法 SPOT-RNA 在 MCC 中提高了 9%，在 F1 分数中提高了 10% 以上。

        图 2a 中显示的结果是碱基对水平的整体性能。图 2b 显示了各个 RNA 中 F1 得分的中位数、第 25 个百分位数和第 75 个百分位数的分布。 SPOT-RNA 具有最高的中位 F1 分数以及表现最差的 RNA 的最高 F1 分数 (0.348)，而所有其他方法的 F1 分数几乎为 0。这凸显了 SPOT-RNA 相对于所有其他基于折叠的技术（包括混合了热力学和机器学习模型的 mxfold）的高度稳定的性能。 SPOT-RNA 与 TS1 上的次佳 mxfold 之间的差异具有统计显着性，通过配对 t 检验获得的 P 值 < 0.006。此外，我们还根据 TS1 上 SPOT-RNA、CONTRAfold 和 CentroidFold 的预测碱基对概率计算了整体缺陷（参见“方法”部分）。集合缺陷度量描述了概率结构集合与其相应的天然 RNA 二级结构的偏差，其中 0 代表完美的预测。 SPOT-RNA 的整体缺陷为 0.19，而 CONTRAfold 和 CentroidFold 的整体缺陷分别为 0.24 和 0.25，这表明与其他两个预测因子相比，SPOT-RNA 预测的结构整体与目标结构更相似。 

        由于硬件限制，我们的方法针对最大长度为 500 个核苷酸的 RNA 进行了训练。确定我们的方法在尺寸依赖性方面的表现是很有趣的。由于 TS1 中的最大序列长度为 189，因此，我们通过将分辨率要求放宽至 4 Å，并包括与其他 RNA 复合的 RNA 链（但忽略 RNA 间碱基对），向 TS1 添加了 32 个序列长度从 298 到 1500 的 RNA。 。之所以将分辨率放宽至4 Å并包含与其他RNA复合的RNA链，是因为PDB中没有太多高分辨率的单链长RNA。补充图 1 将 SPOT-RNA 给出的每个 RNA 的 F1 分数与下一个最佳 mxfold 的 F1 分数作为 RNA 长度的函数进行了比较。正如预期的那样，这两种方法都存在较长 RNA 链性能较低的趋势。在我们的方法训练的 500 个核苷酸范围内，SPOT-RNA 的性能始终优于 mxfold。补充图 1 还显示，mxfold 表现更好，平均 F1 得分为 0.50，而 SPOT-RNA 在 21 个长 RNA 上的平均 F1 得分为 0.35 (L>1000)。我们发现SPOTRNA性能不佳主要是因为SPOT-RNA无法捕获序列间隔>300的超长距离对。这种失败是由于训练中长RNA数据有限造成的。相比之下，无论碱基对序列位置之间的距离如何，mxfold 中的热力学算法都可以找到全局最小值。

        补充图 1. 对于大多数 RNA，除了那些长 RNA (>1000) 外，SPOT-RNA 的性能都提高了 mx 倍。 SPOT-RNA 与 TS1 加 32 长序列测试集中单个 RNA 的第二最佳预测因子 (mxfold) 的比较。绿色表示 SPOT-RNA 对 mxfold 的改进，而洋红色表示对 mxfold 缺乏改进。 

        上述比较可能对我们的方法有偏差，因为几乎所有其他比较方法只能预测规范碱基对，其中包括 Watson-Crick (A-U 和 G-C) 对和 Wobble 对 (G-U)。为了解决这种潜在的偏差，表 2 进一步比较了 SPOT-RNA 与其他在 TS1 上分别在规范对、Watson-Crick 对（A-U 和 G-C 对）和 Wobble 对 (G-U) 上的性能。事实上，当非规范对被排除在性能测量之外时，所有方法的性能都会得到提升。 SPOT-RNA 继续保持最佳性能，规范对和 Watson-Crick 对的 F1 分数比次佳 mxfold 提高了 6%，Wobble 对比次佳 ContextFold 提高了 7%。 mxfold 在预测摆动对方面表现不佳，仅排名第四。

        与假结相关的碱基对对于基于折叠和基于机器学习的方法都具有挑战性，因为它们通常与难以预测的三级相互作用相关。为了直接比较预测假结中碱基对的能力，我们将假结对定义为可以去除以产生无假结二级结构的碱基对的最小数量。程序 bpRNA34（可在 https://github.com/hendrixlab/bpRNA 获取）用于从天然和预测的二级结构中获取假结中的碱基对。表 3 比较了 SPOT-RNA 与所有其他 12 种方法的性能，无论它们是否可以处理独立测试 TS1 中至少有一个假结的 40 个 RNA。由于其他方法都无法预测多重峰，因此我们在分析中忽略与多重峰相关的碱基对。 mxfold 仍然排名第二，仅次于 SPOT-RNA，尽管它无法预测假结，因为假结中的碱基对数量仅占所有碱基对的 10%（见补充表 7）。表 3 显示，对于与假结相关的碱基对，所有方法都表现不佳，F1 分数 < 0.3。尽管这个问题具有挑战性，但 SPOT-RNA 在 F1 分数上比次优 (pkiss) 显着提高了 52%。

图 3 测试集 TS1 中 40 个假结 RNA 上所有预测因子的性能。

        非规范对、三联体和孤碱基对也与假结以外的三级相互作用相关。这里，孤立碱基对是指没有相邻碱基对的单个碱基对（即，在不存在[i−1,j+1]和[i+1,j−1]的情况下为[i,j]）。三联体是指一个碱基与另外两个碱基形成碱基对的罕见情况。如补充表 5 所示，与 CycleFold 相比，SPOT-RNA 在预测非规范碱基对方面的 F1 得分提高了 47%。 SPOT-RNA预测的灵敏度虽然较低（15.4%），但精度却很高，达到73.2%。观察到三胞胎和孤对（F1 分数 < 0.2）的性能非常低。

        RNA 的二级结构以其布局中的结构基序为特征。对于每个天然或预测的二级结构，二级结构基序由程序 bpRNA34 进行分类。不同方法预测不同二级结构基序中碱基的性能如表 4 所示。根据 F1 评分，SPOT-RNA 在茎碱基对（较次优提高 6%）、发夹环核苷酸（ 8% 改进）和凸出核苷酸（11% 改进），尽管它在多环中的表现比 CONTRAfold 稍差（2%）。 mxfold 最适合内部循环预测，比第二好的预测器 Knotty 高出 18%。为了证明 SPOT-RNA 预测三级相互作用以及规范碱基对的能力，补充图图 2 和图 3 显示了来自 TS1 的两个示例（核糖开关41 和 t-RNA42），分别具有高性能和平均性能。对于这两个示例，SPOT-RNA 能够预测非规范碱基对（绿色）、假结碱基对和孤对碱基对（蓝色），而 mxfold 和 IPknot 仍然无法成功预测非规范碱基对和假结碱基对。

        为了进一步确认 SPOT-RNA 的性能，我们编译了另一个测试集（TS2），其中包含 NMR 解析的 39 个 RNA 结构。与 TS1 一样，TS2 通过使用 CD-HIT-EST 和 BLAST-N 对我们的训练数据变得非冗余。图 3a 比较了 SPOT-RNA 给出的精确回忆曲线与测试集 TS2 上的 12 个其他 RNA 二级结构预测因子。 SPOT-RNA 在此测试集上的表现优于所有其他预测因子（补充表 6）。此外，图 3b 显示了各个 RNA 中 F1 得分的中位数、第 25 个百分位数和第 75 个百分位数的分布。 SPOT-RNA 以最小的波动实现了最高的中位 F1 分数，尽管 SPOT-RNA 与单个 RNA 上的次佳（这次是 Knotty）之间的差异（如补充图 4 所示）并不显着，获得的 P 值 < 0.16通过配对t检验。 SPOT-RNA 的 TS2 上的整体缺陷最小（SPOT-RNA 为 0.14，而 CentroidFold 和 CONTRAfold 分别为 0.18 和 0.19）。在这里，我们没有比较伪结中的性能，因为该数据集中伪结中的碱基对数量（总共 21 个）太少，无法进行具有统计意义的比较 

        此外，根据 CD-HIT-EST 和 BLAST-N，我们发现总共 6 个具有最近解决的结构（2019 年 3 月 9 日之后）的 RNA 对于我们的训练集（TR0 和 TR1）和测试集（TS1和 TS2)。对合成构建体 RNA（2019 年 6 月 26 日发布，PDB ID 6dvk 中的 H 链）43 的预测与图 4a 中的天然结构进行了比较。对于这种合成 RNA，SPOT-RNA 产生的结构拓扑与天然二级结构非常相似，F1 得分为 0.85，精度为 97%，灵敏度为 77%。特别是，SPOTRNA 正确捕获了 G46 和 A49 之间的一个非规范碱基对，但错过了假结中的其他碱基对。谷氨酰胺 II 核糖开关（PDB ID 6qn3 中的链 A，于 2019 年 6 月 12 日发布）44 和合成构建 Hatchet 核酶（PDB ID 6jq6 中的链 U，于 2019 年 6 月 12 日发布）45 的 SPOT-RNA 预测与其各自的天然二级结构分别如图4b、c所示。对于这两种 RNA，实验证据表明二聚化中存在链交换 44,45。因此，它们的单体天然结构是通过用其原始立场替换交换的立场来获得的。 SPOT-RNA 能够预测茎和假结（蓝色），Glutamine II Riboswitch 的总体 F1 得分为 0.90。对于 Hatchet Ribozyme，SPOT-RNA 能够预测类似天然的结构，F1 得分为 0.74，尽管它遗漏了非规范碱基对和假结碱基对。

        其他三种 RNA 是 Pistol Ribozyme（PDB ID 6r47 中的链 A 和 B，于 2019 年 7 月 3 日发布）46、Mango Aptamer（PDB ID 6e8u 中的链 B，于 2019 年 4 月 17 日发布）47 和腺病毒病毒相关 RNA （PDB ID 6ol3 中的链 C，2019 年 7 月 3 日发布）48。 SPOT-RNA 在手枪核酶、芒果适体和腺病毒病毒相关 RNA 上的 F1 分数分别为 0.57、0.41 和 0.63。对于这种性能水平，显示 RNA 二级结构的一维表示更具说明性（图 5a-c）。图中显示，Pistol Ribozyme 和 Mango Aptamer RNA 的性能相对较差，部分原因是大量非规范碱基对的罕见存在（绿色）。对于腺病毒病毒相关RNA（VA-I），SPOT-RNA的预测效果很差。它包含三个假阳性茎，带有错误预测的假结（图 5c）。

        这 6 种 RNA 与其他 12 个二级结构预测因子的性能比较如图 6 所示。SPOT-RNA 在合成构建体 RNA（图 6a）、谷氨酰胺 II 核糖开关（图 6b）和手枪核酶（图 6b）方面优于所有其他预测因子。图6c)。它在 Mango Aptamer 中并列第一（与 mxfold 相同）（图 6e），在 Hatchet Ribozyme 中排名第二（仅落后于 mxfold）（图 6d）。然而，与其他方法相比，它在腺病毒相关RNA（RNA puzzle2017 的一部分）上表现不佳（图6f）。与其他方法相比，这种预测效果很差，可能是因为这种没有假结的密集接触的碱基配对网络（除了非规范碱基对引起的网络）最适合基于折叠的算法，可以最大化堆叠的规范碱基对的数量。

3.讨论

        这项工作开发了纯粹基于深度神经网络学习单个 RNA 序列的 RNA 二级结构预测方法。由于只有少量高分辨率的RNA结构可用，深度学习模型必须首先使用根据比较分析注释的RNA二级结构（bpRNA）的大型数据库进行训练，然后对源自RNA二级结构的精确二级结构进行迁移学习。 3D 结构。尽管 bpRNA 中的轻微噪声数据导致精度提高了约 96%（图 2a），但通过迁移学习生成的模型比基于直接学习 TS1 的模型有了显着的改进（F1 分数提高了 30%）。无需基于折叠的优化，迁移学习模型产生的方法不仅可以预测规范碱基对，还可以预测那些通常与三级相互作用相关的碱基对，包括假结、孤碱基对和非规范碱基对。通过对 62 个 RNA 高分辨率 X 射线结构的独立测试，与现有的 12 种二级结构预测技术进行比较，该方法（SPOT-RNA）的精度达到 93%，比第二种方法提高了 13%。当 SPOT-RNA 的灵敏度设置为 50.8%（如 mxfold 中所示）时，最佳方法 mxfold。

        纯机器学习方法的优点之一是所有碱基对都可以训练和预测，无论它是否与局部或非局部（三级）相互作用相关。相比之下，基于折叠的方法必须具有准确的能量参数来捕获非规范碱基对，以及用于全局最小搜索的复杂算法以解决伪结。 SPOT-RNA 代表了预测非规范碱基对方面的重大进步。其 F1 分数比 CycleFold 提高了 47%，从 17% 提高到 26%，尽管两种方法的灵敏度较低，约为 16%（补充表 5）。 SPOTRNA 还可以实现假结中碱基对的最佳预测，尽管所有方法的性能仍然较低，SPOT-RNA 的 F1 分数为 0.239，次优的 F1 分数为 0.157（pkiss，表 3）。这主要是因为结构数据集中假结中的碱基对数量较少（TS1 中每个假结 RNA 平均有 3-4 个碱基对，参见补充表 7）。此外，许多堆叠碱基对的长茎比假结中的一些非局部碱基对更容易学习和预测。作为未来方法开发的参考，我们还检查了 SPOT-RNA 捕获三重相互作用的能力：一个碱基与其他两个碱基配对。精度和灵敏度均较低（分别为 12% 和 7%，补充表 5）。这主要是因为预训练的bpRNA中缺乏碱基三联体的数据，并且结构训练集TR1中三联体和四联体的数量都只有1194个。        

        为了进一步确认性能，将 SPOT-RNA 应用于通过 NMR (TS2) 确定的 39 个 RNA 结构。与 X 射线结构不同，核磁共振确定的结构是通过最小化基于实验距离的约束而产生的。这 39 个 NMR 结构较小，平均长度为 51 个核苷酸，假结中总共只有 21 个碱基对。因此，所有方法的预测都更容易（TS1 的除 SPOT-RNA 之外的 MCC < 0.7，但 TS2 中的大多数方法 > 0.74）。尽管如此，与其他 12 个预测因子相比，SPOT-RNA 仍然具有最佳性能（图 3、补充表 6 和补充图 4）。此外，SPOT-RNA 的性能在 PDB 中最近发布的 6 个非冗余（TR0 和 TR1）RNA 上进行了测试。 SPOT-RNA 在 6 个 RNA 中的 4 个中表现最好或与最好相同，在 1 个中表现第二好（图 6）。

        SPOT-RNA 的局限性之一是，由于我们的硬件限制，它是由短于 500 个核苷酸的 RNA 训练的。在 500 个核苷酸内，SPOT-RNA 比现有技术提供了一致的改进（补充图 1）。然而，对于非常长的 RNA 链（>1000），纯粹基于机器学习的技术并不像一些基于折叠算法的方法（例如 mxfold）那么准确，如补充图 1 所示。缺乏对长 RNA 的训练是主要原因。目前，即使没有硬件限制，PDB结构中>500个核苷酸的高分辨率RNA结构的数量也太少，无法提供足够的训练。因此，现阶段，SPOT-RNA最适合RNA长度<500。

        除了预测精度之外，RNA 二级结构预测还需要高计算效率，因为通常需要基因组规模的研究。我们发现，在 32 核 Intel Xenon(R) E5-2630v4 CPU 的单线程上预测测试集 TS1 中所有 62 个 RNA 的 CPU 时间为 540 秒，比 Knotty（2800 秒）快，但比 IPknot 慢(1.2 秒)、ProbKnot (13 秒) 和 pkiss (112 秒)。然而，我们的分布式版本可以轻松地在多个 CPU 线程或 GPU 上运行。例如，通过在单个 Nvidia GTX TITAN X GPU 上运行 SPOT-RNA，预测所有 62 个 RNA 的计算时间将减少到 39 秒。因此，SPOT-RNA可以切实用于基因组规模的研究

        这项工作使用单个 RNA 序列作为唯一的输入。值得注意的是，在二级结构预测中，仅依靠单个序列就可以获得比现有折叠方法更准确的方法。对于蛋白质接触图预测，PSIBLAST40 和 HHblits49 生成的进化图谱以及同源序列 50 之间的直接耦合分析是关键输入向量，负责最近高精度预测的改进。因此，人们预计，从 BLAST-N 和直接/进化耦合分析生成的类似进化衍生序列谱将进一步改善长 RNA 中非局部碱基对的二级结构预测。事实上，最近我们已经证明，使用进化衍生的序列图谱可以显着提高预测 RNA 溶剂可及性和灵活性的准确性38,39。例如，如果单个序列被 BLAST-N38 的序列图谱替换，则预测溶剂可及性与实际溶剂可及性之间的相关系数从 0.54 增加到 0.63。然而，序列概况和进化耦合的生成在计算上非常耗时。由此产生的改进（或缺乏改进）很大程度上取决于当前 RNA 序列数据库中可用的同源序列的数量。如果同源序列的数量太少（对于大多数 RNA 来说都是如此），它可能会引入比预测信号更多的噪声，如蛋白质二级结构和内在无序预测所示51,52。此外，合成的RNA不会有任何同源序列。因此，我们在本研究中提出了以单序列信息作为输入的方法。使用序列图谱和进化耦合作为 RNA 二级结构预测的输入正在进行中。

        进一步改进 SPOT-RNA 的另一种可能方法是利用预测概率作为具有适当评分函数的折叠约束。这种双重方法可能会改进 SPOT-RNA，因为折叠优化可能具有更好的能力来捕获长 RNA 的 WC 对之间的非局部相互作用，特别是如补充图 1 所示。然而，简单的集成可能不会产生对于较短的链（<500）有很大的改进。在 mxfold 中，将机器学习和热力学模型相结合，在一个测试集中导致 0.6%，在另一测试集中导致 5%33。此外，大多数热力学方法只是忽略非规范碱基对，许多甚至不考虑假结。例如，mxfold 采用无假结热力学方法与其机器学习模型相结合。因此，平衡规范、非规范和伪结的性能需要仔细选择适当的评分方案。简单的集成可能会导致一种类型的碱基对具有高性能，但会牺牲其他类型的碱基对。尽管如此，我们发现，如果我们简单地保留预测三重相互作用中预测概率最高的碱基对，SPOT-RNA 的 F1 分数将再提高 3%（从 TS1 中的 0.69 提高到 0.71），这证实了存在一些改进的余地。我们将把这个推迟到未来的研究中。

        二级结构预测性能的显着提高应该可以大大改善 RNA 3D 结构的建模。这是因为该方法不仅预测规范碱基对，而且还提供非规范和非嵌套碱基对的重要三级接触。因此，它可以作为一个更准确的准三维框架，以实现正确折叠成正确的 RNA 三级结构。 RNA Puzzles（盲 RNA 结构预测）53 证明了 2D 结构预测对于 3D 结构建模的有用性。此外，预测二级结构基序（茎、环和凸起，参见表 4）的改进将允许更好的功能推理 54,55、序列比对 56 和 RNA 抑制剂设计 57。该方法数据集可作为服务器和独立软件公开在http://sparks-lab.org/jaswinder/server/SPOT-RNA/ 和 https://github.com/jaswindersingh2/SPOT-RNA/。

4.方法

4.1 数据集

        初始训练的数据集来自 bpRNA-1m（版本 1.0）34，它由 102,348 个带有注释二级结构的 RNA 序列组成。使用CD-HIT-EST37去除序列相似性超过80%的序列。大约 80% 的序列同一性截止值是 CD-HIT-EST 允许的最低截止值，之前已被用作 RNA 非冗余截止值38,39。去除序列相似性后，剩下 14,565 个序列。由于我们仅基于具有 PDB 结构的 RNA 准备了单独的数据集，因此该数据集中可用的 PDB5 中具有 RNA 结构的 RNA 序列也被删除了。此外，由于长序列训练的硬件限制，最大序列长度被限制为500。经过预处理，该数据集包含13,419个序列。这些序列被随机分成 10,814 个 RNA 用于训练 (TR0)、1300 个用于验证 (VL0) 和 1,305 个用于独立测试 (TS0)。补充表 7 显示了 RNA 序列的数量及其 Watson-Crick（A-U 和 G-C）、Wobble（G-U）和非规范碱基对计数以及与假结相关的碱基对数量。 TR0、VL0 和 TS0 中的平均序列长度均约为 130。此处，与假结相关的碱基对定义为可以去除以产生无假结二级结构的碱基对的最小数量。假结标签是使用软件 bpRNA34（可在 https://github.com/hendrixlab/bpRNA 获取）生成的。

        用于迁移学习的数据集是通过20195年3月2日从PDB下载高分辨率（<3.5 Å）RNA获得的。这些序列中相似度超过80%的序列被CD-HIT-EST37删除。去除序列相似性后，只剩下 226 个序列。这些序列被随机分成 120、30 和 76 个 RNA，分别用于训练 (TR1)、验证 (VL1) 和独立测试 (TS1)。此外，TS1中与TR0序列相似性超过80%的任何序列也被删除，这使得TS1的RNA减少到69个。由于 CD-HIT-EST 只能去除相似度超过 80% 的序列，因此我们采用 BLAST-N40 进一步去除与 e 值截止值为 10 的训练数据潜在的序列同源性。该过程进一步将 TS1 从 69 降低到 67 RNA。

        为了进一步对 RNA 二级结构预测进行基准测试，我们采用了 NMR 解析的 641 个 RNA 结构。使用同一性截止值为 80% 的 CD-HIT-EST，然后针对 TR0、TR1 和 TS1 使用 e 值截止值为 10 的 BLAST-N，我们获得了 39 个 NMR 解析的结构作为 TS2。

        所有PDB集的二级结构都是使用DSSR58软件从它们各自的结构推导出来的。对于 NMR 解析的结构，使用模型 1 结构，因为它被认为是所有结构中最可靠的结构。所有集合的规范、非规范和假结碱基对以及碱基多重态（三重态和四重态）的数量列于补充表 7 中。这些数据集以及带注释的二级结构可在 http://sparks-lab 上公开获取。 org/jaswinder/server/SPOT-RNA/ 和 https://github.com/jaswindersingh2/SPOT-RNA。

 4.2 RNA 二级结构类型。

        对于不同 RNA 二级结构类型的分类，我们使用了 bpRNA34 之前使用的相同定义。茎被定义为不间断碱基对的区域，没有中间的环或凸起。发夹环是不配对核苷酸的序列，其两端在茎区的两条链处相遇。内部环被定义为两侧都有闭合碱基对的两条不成对的链。凸出是内部环的一种特殊情况，其中一根线的长度为零。多环由两个以上不成对的链组成，通过茎连接。 TR1、VL1和TS1中不同二级结构类型的分布（不包括多重碱基对）如补充表8所示。这些二级结构分类是通过使用二级结构分析程序bpRNA34获得的。

4.3深度神经网络

        我们采用了深度学习神经网络集合进行预训练。该集合由 5 个基于其在 VL0 上的表现而排名靠前的模型组成，其架构如图 1 所示，与之前 SPOT-Contact29 中用于蛋白质接触预测的架构类似。

        每个模型的架构均由 ResNet 块、后跟 2DBLSTM 层和全连接 (FC) 块组成。正如 He 等人 30 中提出的，在我们的 ResNet 块之前使用了用于预激活的初始卷积层。初始卷积层后面是 NA ResNet 块（图 1 中的块 A）。每个 ResNet 块由两个卷积层组成，kernel size大小分别为 3 x3 和 5 x 5，深度为 DRES。使用指数线性单元 (ELU)59 激活函数和层归一化技术60。在每个卷积层之前使用 25% 的丢弃率，以避免训练期间过度拟合。在某些模型中，我们使用了扩张卷积，据报道可以更好地学习较长范围的依赖性62。对于扩张的卷积层，扩张因子设置为2i%n，其中i是卷积层的深度，n是固定标量，%是模运算符。

        该架构中的下一个模块是 2D-BLSTM31,32。最终 ResNet 块的输出在作为 2D-BLSTM 的输入之前被激活（使用 ELU）和归一化（使用层归一化）。每个LSTM方向单元的节点数为DBL。在 2D-BLSTM 之后，使用具有 DFC 节点的 NB FC 层，如图 1 中的块 B 所示。每个 FC 层的输出均使用 ELU 函数激活，并使用层归一化技术进行归一化。隐藏 FC 层采用 50% 的 dropout 率，以避免过度训练。该架构的最后阶段由一个带有一个节点的输出 FC 层和一个 sigmoidal 激活函数组成。 sigmoid 函数将输出转换为每个核苷酸与其他核苷酸配对的概率。输出的数量等于大小为 Ĺ L 的上三角矩阵中的元素数量，其中 L 是序列的长度。

        每个模型均在 Google 的 Tensorflow 框架 (v1.12)63 中实现，并使用具有默认参数的 ADAM 优化算法 64 进行训练。所有模型均在 Nvidia GTX TITAN X 图形处理单元 (GPU) 上进行训练，以加快训练速度65。我们基于图 1 所示的架构，通过在 NA、DRES、DBL、NB 和 DFC 上执行超参数网格搜索，在 TR0 上训练了多个深度学习模型。 NA、DRES、DBL、NB、DFC 分别从 16 至 32、32 至 72、128 至 256、0 至 4 和 256 至 512 进行搜索。这些模型在 VL0 上进行了优化并在 TS0 上进行了测试。然后使用迁移学习在 TR1 上进一步训练这些模型。在迁移学习过程中，VL1用作验证集，TS1用作独立测试集。

4.4迁移学习

        迁移学习35涉及进一步训练一个大型模型，该模型是在大型数据集上针对特定任务进行训练的，并将其训练到具有有限数据的其他相关任务。在这个项目中，我们使用大型数据集 bpRNA 进行初始训练，然后使用小型 PDB 数据集进行迁移学习，如图 1 所示。在 TR0 上学习的所有权重/参数都经过重新训练，以便在 TR1 上进行进一步训练。 在迁移学习期间，训练和验证标签的格式与初始训练完全相同，为二维 (2D) ĹxL 上三角矩阵，其中 L 是 RNA 序列的长度。迁移学习期间使用的所有标签均源自 PDB 中的高分辨率 X 射线结构。迁移学习中的一些方法会冻结特定层的权重并训练其他层。在这里，我们训练了模型的所有权重，而不冻结任何层，因为这提供了更好的结果。之前关于蛋白质分子识别特征 (MoRF) 预测的工作36还表明，通过重新训练所有权重来使用迁移学习可以提供比在重新训练期间冻结某些层更好的结果

        在 TS1 上的迁移学习过程中，我们使用了与 TS0 训练模型相同的超参数（层数、层深度、内核大小、膨胀因子和学习率）。所有模型都针对 VL1 进行了验证，并根据这些模型在 VL1 上的性能，为集成选择了 5 个最佳模型。这些模型的参数如补充表9所示。

4.5输入

        SPOT-RNA 的输入是由大小为 Ĺ 4 的二元 onehot 向量表示的 RNA 序列，其中 L 是 RNA 序列的长度，4 对应于碱基类型 (A、U、C、G) 的数量。在 one-hot 编码中，值 1 分配给向量中相应的基本类型位置，其他位置分配为 0。在 onehot 编码向量中分配了残差值为 -1 的缺失或无效序列。

        该一维 (L x 4) 输入特征通过 RaptorX-Contact28 中所述的外部串联函数转换为二维 (L ́xL ́x 8)。在输入模型之前，输入被标准化为具有零均值和单位方差（根据训练数据）。

4.6输出

        我们模型的输出是一个二维 (2D) L x L 上三角矩阵，其中 L 是 RNA 序列的长度。该上三角矩阵表示每个核苷酸与序列中任何其他核苷酸配对的可能性。使用单个阈值来决定核苷酸是否与任何其他核苷酸配对。阈值的选择方式可以优化验证集的性能。

 4.7性能指标。

        RNA二级结构预测是一个二元分类问题。我们使用灵敏度、精度和 F1 分数来衡量性能，其中灵敏度是所有本地碱基对中预测的碱基对的分数 (SN = TP/(TP + FN))，精度是正确预测的碱基对的分数 (PR= TP/(TP + FP)，F1 分数是其调和平均值 (F1 = 2(PR *SN)/(PR + SN))。这里，TP、FN和FP分别表示真阳性、假阴性和假阳性。除了上述强调积极因素的指标外，还使用了平衡指标马修斯相关系数（MCC）66。 MCC 计算为

        

        其中 TN 表示真阴性。 MCC 衡量预期类别与获得类别之间的相关性。此外，使用精确回忆（灵敏度）曲线将我们的模型与当前可用的 RNA 二级结构预测因子进行比较。为了显示 SPOT-RNA 相对于第二最佳预测因子的改进的统计显着性，对 F1 评分使用配对 t 检验以获得 P 值 67。 P 值越小，两个预测变量之间的差异越显着。由于 SPOT-RNA 的输出是碱基对概率，因此我们可以使用集成缺陷作为附加性能指标。集合缺陷描述了预测的碱基对概率和目标结构之间的相似性68。可以通过向未配对碱基的预测概率矩阵和目标矩阵添加额外的列来计算。

        如果 P 和 S 分别是预测结构和目标结构，并且在附加额外列之后 P' 和 S' 是预测结构和目标结构，则整体缺陷 (ED) 由下式给出

        

其中 L 是序列的长度。 ED的值越小，预测的碱基对概率与目标结构之间的结构相似性越高。 

4.8方法比较

        我们将 SPOT-RNA 与 12 个最佳可用预测因子进行了比较。我们下载了单机版的mxfold33（可在https://github.com/keio-bioinformatics/mxfold获取）、ContextFold16（可在https://www.cs.bgu.ac.il/negevcb/contextfold/获取） 、CONTRAfold14（可在 http://contra.stanford.edu/contrafold/ 获取）、Knotty24（可在 https://github.com/HosnaJabbari/Knotty 获取）、IPknot23（可在 http://rtips.dna.bio 获取）。 keio.ac.jp/ipknot/）、RNAfold11（可在 https://www.tbi.univie.ac.at/RNA/ 获取）、ProbKnot22（可在 http://rna.urmc.rochester.edu/RNAstruct 获取）。 html）、CentroidFold15（可在 https://github.com/satoken/centroid-rna-package 获取）、RNAstruct12（可在 http://rna.urmc.rochester.edu/RNAstruct.html 获取）、RNAshapes13（可在 https://rna.urmc.rochester.edu/RNAstruct.html 获取） //bibiserv.cebitec.uni-bielefeld.de/rnashapes）、pkiss13（可在 https://bibiserv.cebitec.uni-bielefeld.de/pkiss 获取）和 CycleFold27（可在 http://rna.urmc 获取）。 rochester.edu/RNAstruct.html）。在大多数情况下，除了 pkiss 之外，我们都使用默认参数进行二级结构预测。在 pkiss 中，我们使用的策略 C 较慢但彻底，而策略 A 和 B 则较快但不太准确。对于 CONTRAfold 和 CentroidFold，它们的性能指标源自预测的碱基对概率以及最大化 MCC 的阈值。

4.9报告摘要

有关研究设计的更多信息，请参阅本文链接的《自然研究报告摘要》。