一种基于深度学习的RBP结合位点预测算法的制作方法

一种基于深度学习的rbp结合位点预测算法
技术领域
1.本发明涉及rna
‑
蛋白质相互作用与深度学习领域，具体涉及一种基于深度学习的rbp结合位点预测算法。


背景技术：

2.rna结合蛋白(rbp)是普遍存在的能与rna结合的蛋白质的总称，rbp通过与rna特异性结合，直接或间接地调控rna的功能，如介导rna的成熟、转运、定位和翻译等。rbp参与了许多生命进程，它们的数量占据了真核蛋白质组的5
‑
10％。在生物体中大部分的rna发挥作用都要与蛋白质相结合，rbp与特异rna结合的过程又称为rna
‑
蛋白质相互作用，此时存在结合位点。鉴定rbp在特异性rna上的结合位点是rna
‑
蛋白质相互作用的典型功能表现形式。在物理度量空间中，蛋白质中特定氨基酸残基和rna中特定碱基之间的最小距离小于时，蛋白质与rna发生结合。蛋白质与rna发生结合的一小段序列被称为rna的motif(模序)，motif是一组rna序列共有的短序列模式，此处也指rna的特征，具体表现形式是四种碱基(或者结构注释)分别出现在短序列各个位置上的概率，一般长度为7bp。
3.蛋白质和rna都是人体中不可或缺的生物大分子，二者结合形成的核糖核蛋白复合体已被确定在许多重要的生物功能中都起着重要作用。如何利用计算方法提高其预测精度成为当今国内外研究的热点，现阶段研究中，由于对蛋白质与rna的结合机制以及结合位置尚未完全了解，传统的机器学习方法特征提取的困难，目前还未出现一种公认的特征提取方法。
4.因此，研究者们纷纷转向研究深度学习方法。一方面，大多数深度学习方法往往仅从rna序列出发，忽略了对应的rna结构信息，而研究表明rna结构的修饰及其空间结构的改变都可以导致与其结合的rbp发生变化，产生不同的生物学功能，因此，蛋白质与rna的结合对特定rna的结构具有一定的倾向性，引入rna二级结构预测具有一定的理论意义。另一方面，多数方法仅考虑结合位点是否存在，忽略了结合发生的序列特异性，研究表明rna结合蛋白的序列特异性对于开发生物系统中调控过程的模型以及确定病因疾病变体至关重要，因此，探明结合位点的序列特异性对于了解结合机制必不可少。另一方面，由于某些rbp数据规模较小，模型预测能力始终无法提升，因此，对于离散rna数据，如何提升数据量较少的模型的预测能力尤为重要。最后，大多数模型训练时参数随机初始化，容易造成局部最优，因此如何避免局部最优也是提高结合位点预测的关键点。
5.研究发现，rbp与许多重要的生物学过程密切相关，例如基因的转录后调控、基因剪接和定位等。rbp中ago蛋白可以与mirna和mrna靶标形成复合物，mirna充当靶标mrna的阻遏物。此外，rbp的失调与突变可能导致多种疾病，某些rbp的缺失会破坏转录后机制引发糖尿病。rbp中fus蛋白和tdp
‑
43蛋白的突变与肌萎缩性侧索硬化症密切相关。通过研究rbp结合特异rna，可以揭示rbp在疾病背后的调控机制。更进一步地，想要探索某些rna的功能，也可以通过与之相关联的蛋白质进行了解。因此，研究蛋白质和rna之间的相互作用已成为后基因组时代的一个有吸引力的主题。


技术实现要素：

6.针对上述问题，本发明的目的在于针对现有技术中的预测算法精度不高的现状提供一种基于深度学习的rbp结合位点预测算法，将深度学习技术与特异rna的领域知识相融合，鉴定特异rna上的rbp结合位点，也对其中的部分问题提出了有效解决方案，例如一些rbp数据规模小、预测模型过拟合与局部最优、以及无法提取序列与结构特异性等。
7.本发明所解决的技术问题可以采用以下技术方案来实现：
8.一种基于深度学习的蛋白质
‑
配体结合位点预测算法，包括如下步骤：
9.步骤1)对于给定的rbp数据集，首先提取序列信息和独立测试集；再根据序列信息，预测出二级结构信息；根据序列和结构信息，分别构建序列和结构的编码矩阵；并构建独立测试集的编码矩阵，对应于每个rbp创建正样本和负样本。
10.步骤2)构建包括生成网络和判别网络的生成对抗网络；在训练过程中，采用强化学习对判别网络优化，使优化过后的判别网络促使生成网络生成更逼真的合成数据，并在步骤1)中所得的序列和结构的编码矩阵上分别进行训练；
11.步骤3)构建卷积自编码器预测模型，将步骤1)中所得到的编码矩阵和步骤2)中所得到的合成数据合并，分别得到数据规模更大的序列和结构矩阵，将序列和结构矩阵分别训练；
12.步骤4)分别提取卷积自编码器中序列和结构编码器中卷积核的参数，分别与rna序列和结构的编码矩阵进行卷积运算，针对每个rbp，将所有滑动窗口中超过规定阈值的短序列拼接起来，对每组短序列进行序列比对，统计序列中每个位置上的一致性，作为motif；
13.步骤5)对每种rbp数据均训练一个对应的预测模型；对于一个待预测的rna序列，首先编码序列，然后将编码信息输入到已训练好的各个预测模型中，最终的结果即为特异rna与各个rbp对应的结合概率。
14.进一步的，所述步骤1)中序列与结构特征矩阵的构造方法如下：
15.步骤1)中构建序列和结构的编码矩阵过程如下：
16.步骤1.1)设定rna序列长度为l，对rna序列和标签进行one
‑
hot编码；对于长度不足l的序列，采用“n”在末尾填充，“n”为等概率填充，“n”的向量表示为[0.25,0.25,0.25,0.25]；对于超过l的序列，将截取长度l作为训练序列；对于碱基a、c、g、u序列编码后编码矩阵维度为l
×
4；
[0017]
步骤1.2)rna二级结构为由focusfold算法预测得到“点
‑
括号”格式，一对左右“括号”表示存在一对碱基配对，而“点”则表示未发生配对，每个符号分别按顺序对应于每个碱基位点，然后根据碱基配对原则，得到每个碱基位点的二级结构注释，所述二级结构注释包括堆叠s、自由端f、连接处j、发卡环h、内环与多环m，对rna结构进行one
‑
hot编码，结构编码后矩阵维度为l
×
5；与序列编码类似，结构编码长度不足l的序列采用“n”等概率向量填充，表示为[0.2,0.2,0.2,0.2,0.2]；长度超过l的序列将截取前l长度作为训练序列。
[0018]
更进一步的，所述步骤2)中生成合成数据的过程如下：
[0019]
步骤2.1)所述生成对抗网络包括生成网络和判别网络，其中，生成网络包括一层循环神经网络和一层全连接网络，全连接网络作为输出单元；判别网络包括一层循环神经网络、一层全连接网络和分类器；
[0020]
步骤2.2)将步骤1.1)和步骤1.2)中所得的序列和结构的编码矩阵分别输入到生
成对抗网络训练，采用强化学习技术优化判别网络；被优化的判别网络的返回值反馈给生成网络，多次训练迭代，生成网络生成更逼真的合成数据；
[0021]
步骤2.3)提取步骤2.2)中已训练好的序列和结构生成网络，分别生成合成的序列和结构矩阵。
[0022]
更进一步的，所述生成网络中，循环神经网络采用双向lstm，包含256个计算单元；全连接网络作为输出层，包含100个神经元；
[0023]
所述判别网络中，循环神经网络采用双向lstm，包含256个计算单元；全连接网络作为中间层，包含100个神经元；分类器采用sigmoid激活函数。
[0024]
更进一步的，所述步骤3)中卷积自编码器构建与训练的过程如下：
[0025]
步骤3.1)所述卷积自编码器包括编码器与解码器，在预训练阶段，编码器与解码器共同工作，序列和结构分别训练，输入数据经过编码器的降维学习，得到中间抽象特征，之后解码器重构输入数据；
[0026]
步骤3.2)在微调阶段，提取编码器模型与权重参数，将序列和结构的编码器并行，增加两层lstm层，最终实现结合位点预测。
[0027]
更进一步的，所述编码器包括两层卷积与两层池化；序列编码器的第一层卷积核大小为(4,10)，数量为16，步长为(1,1)，卷积核的长度为motif的1.5倍；
[0028]
在卷积层后面加入batch normalization，将每一层归一化到同一个正态分布，将其转化为均值为0、方差为1的标准正态分布；批归一化的函数运算如下所示:
[0029][0030]
其中，x
l
是mini
‑
batch中的每个样本，将所有的x
l
归一化为标准正态分布，其中涉及到均值μ和方差δ两个参数；γ和β为可学习重构参数；n为batch
‑
size的大小，即mini
‑
batch中的样本量；y
l
为归一化后的数据，作为下一层的输入数据；
[0031]
所述池化层选用最大池化操作，池化窗口大小为(1,3)，池化步长为(1,2)，最大池化操作保留特征窗口内的最大值；在池化层后添加dropout层，对于保留下来的神经元，得到一个节点较少，规模较小的网络进行训练，以达到避免过拟合的目的；输入层神经元相当于第0层，即a
(0)
，从计算公式上来看，对第a
(z)
层加入dropout，则有：
[0032]
a
(z)
＝r
(z)
*a
(z)
[0033]
其中，r
(z)
表示第z层的神经元随机失活的概率分布，其服从于p＝1
‑
dropout的伯努利分布，a
(z)
表示加入了dropout的层；
[0034]
在编码器中，dropout层后继续添加卷积层与batch normalization，卷积核维度为(1,10)，数量为16，步长与上一层卷积层保持一致为(1,1)；
[0035]
构建的编码器函数如下：
[0036]
e＝encode(x)＝maxpooling(σ(conv(maxpooling(σ(conv(x))))))
[0037]
其中，σ为激活函数，conv(x)为卷积函数，x为输入变量，maxpooling()为最大池化函数；
[0038]
所述解码器包含两层反卷积函数和两层上采样函数，第一层上采样层大小为(1,2)，反卷积层的卷积核维度为(1,12)，数量为16，步长为(1,1)，第二层上采样层大小为(1,2)，接下来反卷积层的卷积核维度为(4,10)，数量为1，步长为(1,1)；
[0039]
构建的解码器函数如下：
[0040]
decode(e)＝σ(deconv(upsampling(σ(deconv(upsampling(e))))))
[0041]
其中，σ为激活函数，deconv(e)为反卷积函数，e为编码器的输出，upsampling()为上采样函数。
[0042]
更进一步的，所述步骤4)中motif的提取与展示过程如下：
[0043]
步骤4.1)在卷积自编码器的微调阶段，模型在多次迭代之后，权重被记录在卷积核中，提取第一层卷积核的参数；
[0044]
步骤4.2)对于起点为i，将每个卷积核与步骤1)得到的序列信息在每个位点进行卷积运算，以长度为10的滑动窗口前进，对于每个碱基位点，得到中间值a
ij
，其中j表示是由第j个卷积核完成的卷积操作；对每个卷积核设定阈值threshold
j
，即该卷积核是所有短序列计算得到的中间值a
ij
的最大值的0.5倍，计算如下：
[0045][0046]
若a
ij
超过该阈值，则认为该起点为i，长度为10的短序列是有效的，将该序列加入motif序列集；通过对每个卷积核滑动运算，得到16组对应的短序列集；
[0047]
步骤4.3)将16组短序列集拼接起来，得到一组含有共性的短序列集；
[0048]
步骤4.4)使用工具weblogo 3将每个rbp的一组短序列集进行序列比对，统计出每个rbp的序列中每个位置上的一致性；
[0049]
步骤4.5)通过步骤4.4)所计算出的一致性，绘制序列标识图。
[0050]
更进一步的，所述序列标识图中每个碱基图形字符的大小与残基在该位置上出现的频率成正比，但纵轴数值不使用频率值，而是计算后的结果；对于每个位置，每种碱基图形字符的高度height具体计算方式如下：
[0051]
height＝f
q
*r
[0052]
其中，q＝{a，u，g，c}，f
q
表示该种碱基的出现频率；r表示特定位点存在的信息量，
[0053]
定义为最大可能熵和碱基分布熵之间的差异，计算如下：
[0054]
r＝log2(m)
‑
(entropy e
k
)
[0055]
其中，entropy表示该位置的整体熵值，m表示碱基种类个数，e
k
表示校正因子，具体为
[0056]
小样本检验修正值，计算如下；
[0057][0058]
e
k
＝(m
‑
1)/(2*ln2*k)
[0059]
其中，k表示用于统计的序列个数。
[0060]
与现有技术相比，本发明的有益效果在于：
[0061]
1.本发明从rna序列信息出发，增加二级结构特征作为输入；
[0062]
2.本发明提出了基于强化学习的生成对抗网络，对rbp进行数据增强，目的是通过增强训练集规模来提高预测模型的预测能力；
[0063]
3.本发明提出了卷积自编码器作为预测模型，预训练过程包括编码器和解码器的训练，提取已训练好的编码器微调，能够更好地提取特征，有效地防止过拟合、局部最优和节约训练时间；能够预测未知结合位点的rna序列；
[0064]
4.本发明提出通过卷积自编码器的卷积核，能够提取与rna发生结合的序列motif与结构motif，探究不同rbp的序列和结构特异性；这些都为最终探明rna和蛋白质相互作用的分子机制，以及从本质上认识相关的细胞生理过程打下了坚实的基础。
附图说明
[0065]
图1为本发明的基于深度学习的rbp结合位点预测算法流程图。
[0066]
图2本发明的生成对抗网络示意图。(a)生成网络模块，(b)判别网络模块。
[0067]
图3为本发明的基于序列的卷积自编码的预训练模块示意图。(a)序列编码器，(b)序列解码器。
[0068]
图4为本发明的基于结构的卷积自编码的预训练模块示意图。(a)结构编码器，(b)结构解码器。
[0069]
图5为本发明的卷积自编码的微调模块示意图。
[0070]
图6为本发明的motif提取模块示意图。
具体实施方式
[0071]
下面将结合本发明实施例中的附图，对本发明实施例中的技术方案做进一步详细说明。本发明通过有效的深度学习方法，对rbp数据集做数据增强，构建卷积自编码器预测模型，通过真实数据与合成数据提高预测模型的预测能力，实现对未知结合位点的rna的预测。
[0072]
一种基于深度学习的rbp结合位点预测方法，包括数据生成模型、训练模型和预测模型，数据生成模型通过对已知的rbp数据集进行特征提取与学习，生成逼真的合成数据；训练模型通过对已知rbp数据集和合成数据集一起进行特征提取，完成神经网络的训练，得到预测所用的模型，即预测模型；训练得到的模型对新的特异rna序列进行相同的预处理后，进行结合位点的预测。
[0073]
参见图1，本发明所述的一种基于深度学习的rbp结合位点预测算法，包括以下步骤：
[0074]
步骤1)对于给定的rbp数据集，提取序列信息以及独立测试集；根据序列信息，利用已有的工具focusfold预测出二级结构信息；根据序列和结构信息，分别构建序列和结构的编码矩阵；构建独立测试集的编码矩阵；对应于每个rbp创建正样本和负样本。
[0075]
建序列和结构的编码矩阵具体实施如下：
[0076]
步骤1.1)神经网络要求输入矩阵中的各个编码序列等长，设定rna序列长度为l，对rna序列和标签进行one
‑
hot编码；对于长度不足l的序列，本领域常用等概率“n”在末尾填充，“n”的向量表示为[0.25,0.25,0.25,0.25]；对于超过l的序列，将截取长度l作为训练
序列；对于碱基a、c、g、u，序列编码后矩阵维度为l
×
4；
[0077]
步骤1.2)rna二级结构由已有的focusfold算法预测得到“点
‑
括号”(dot
‑
bracket)格式，一对左右“括号”表示存在一对碱基配对，而“点”则表示未发生配对，每个符号分别按顺序对应于每个碱基位点，然后根据碱基配对原则，得到每个碱基位点的二级结构注释，本发明共得到5类二级结构注释，包括堆叠(s)、自由端(f)、连接处(j)、发卡环(h)、内环与多环(m)，对rna结构进行one
‑
hot编码，结构编码后矩阵维度为l
×
5；与序列编码类似，结构编码长度不足l的序列采用“n”等概率向量填充，表示为[0.2,0.2,0.2,0.2,0.2]；长度超过l的序列将截取前l长度作为训练序列。
[0078]
步骤2)构建生成对抗网络，包括生成网络和判别网络，在训练过程中，采用强化学习对判别网络优化，优化过后的判别网络促使生成网络生成更逼真的合成数据；并在步骤1)中所得的序列和结构的特征矩阵上分别进行训练。
[0079]
生成合成数据的过程具体如下：
[0080]
步骤2.1)本发明构建的生成对抗网络包括生成网络和判别网络，其中，生成网络包括一层循环神经网络和一层全连接网络，全连接网络作为输出单元；判别网络包括一层循环神经网络、一层全连接网络和分类器，分类器使用sigmoid激活函数。
[0081]
本实施例生成网络结构如图2(a)所示，包括一层循环神经网络和一层全连接网络，其中，循环神经网络采用双向lstm，包含256个计算单元；全连接网络作为输出层，包含100个神经元；
[0082]
本实施例判别网络结构如图2(b)所示，包括一层循环神经网络、一层全连接网络和分类器，其中，循环神经网络采用双向lstm，包含256个计算单元；全连接网络作为中间层，包含100个神经元；最后，分类器使用sigmoid激活函数。
[0083]
步骤2.2)考虑到rna是一种离散的文本数据，在步骤1.1)和步骤1.2)中所得的序列和结构的特征矩阵分别输入到生成对抗网络训练，本发明使用强化学习技术优化判别网络；被优化的判别网络的返回值反馈给生成网络，如图2所示，生成网络与判别网络的网络层参数如图所示，多次训练迭代，生成网络可以生成更逼真的合成数据。
[0084]
步骤2.3)提取步骤2.2)中已训练好的序列和结构生成网络，分别生成合成的序列和结构矩阵。
[0085]
步骤3)构建卷积自编码器预测模型，将步骤1)中所得到的编码矩阵和步骤2)中所得到的合成矩阵合并，分别得到数据规模更大的序列和结构矩阵，将序列和结构矩阵分别训练。
[0086]
卷积自编码器构建与训练的过程如下：
[0087]
步骤3.1)本发明中构建卷积自编码器，包括编码器与解码器，在预训练阶段，编码器与解码器共同工作，序列和结构分别训练，输入数据经过编码器的降维学习，得到中间抽象特征，之后解码器重构输入数据；重构后的数据需要与输入数据尽可能地相似，使得编码器和解码器保留更多的特征；使用步骤2.3)中得到的数据集，训练集与标签相同，训练集与验证集比例是8：2。
[0088]
本实施例的卷积自编码器中编码器包括两层卷积与两层池化，序列和结构编码器分别如图3(a)与图4(a)所示，其中，序列编码器的第一层卷积核大小为(4,10)，数量为16，步长为(1,1)，在这里，卷积核的长度与motif的长度相关，motif是一组rna序列共有的短序
列模式，motif可以看作是与同一种rbp发生结合的一组rna序列的特征，基本表现形式是碱基出现在每个位点上的概率；由于网络中的某些操作可能会导致碱基信息的部分丢失，因此，卷积核的长度被设定为motif的1.5倍。现有研究表明，卷积核的数量与motif的种类相关，因此此处设置为16。
[0089]
其次，在卷积层后面加入batch normalization，也就是批标准化，本质上是一个归一化网络层；当网络层的输入值未进行归一化，与权重矩阵运算后，容易产生一些偏离较大的差异值；在反向传播中，这些现象会导致梯度发散，需要更多的训练继续学习那些被发散掉的数据。batch normalization可以将每一层归一化到同一个正态分布，将其转化为均值为0、方差为1的标准正态分布，不改变数据维度，缓解了训练中的梯度消失问题，加快了模型的收敛速度，批归一化的函数运算如下所示：
[0090][0091]
其中，x
l
是mini
‑
batch中的每条数据，将所有的x
l
归一化为标准正态分布，其中涉及到均值μ和方差δ两个参数。γ和β是可学习重构参数；n为batch
‑
size的大小，即mini
‑
batch中的样本量；y
l
为归一化后的数据，作为下一层的输入数据；
[0092]
接下来，池化层选用最大池化操作，池化窗口大小为(1,3)，池化步长为(1,2)，最大池化操作保留特征窗口内的最大值。在池化层后添加dropout层，dropout(随机失活)是对一层的神经元设置一个随机置零的概率，以一定的概率将部分隐藏神经元的权重置零，而输入输出维度保持不变，dropout是另一个超参数；因此，各节点的“重要性”(即权重)会被平衡。对于保留下来的神经元，可以得到一个节点较少，规模较小的网络进行训练，从而达到避免过拟合的目的。本模型中采取的dropout参数为0.25(即任意结点的失活概率为25％)。输入层神经元相当于第0层，即a
(0)
，从计算公式上来看，对a
(z)
层加入了dropout，则有：
[0093]
a
(z)
＝r
(z)
*a
(z)
[0094]
其中，r
(z)
表示第z层神经元随机失活的概率分布，其服从于p＝1
‑
dropout的伯努利分布，a
(z)
表示加入了dropout的层。
[0095]
在编码器中，dropout层后继续添加卷积层与batch normalization，卷积核维度为(1,10)，数量为16，步长与上一层卷积层保持一致为(1,1)，目的是为了保留更多的信息；
[0096]
本发明所构建的编码器函数如下：
[0097]
e＝encode(x)＝maxpooling(σ(conv(maxpooling(σ(conv(x))))))
[0098]
其中，σ为激活函数，conv(x)为卷积函数，x为输入变量，maxpooling()为最大池化函数。
[0099]
本实施例卷积自编码器中解码器包含两层反卷积函数和两层上采样函数，序列和
结构解码器模型结构分别如图3(b)与图4(b)所示，第一层上采样层大小为(1,2)，反卷积层的卷积核维度为(1,12)，数量为16，步长为(1,1)，第二层上采样层大小为(1,2)，接下来反卷积层的卷积核维度为(4,10)，数量为1，步长为(1,1)，通过维度对比可以发现，重构完成后维度与输入维度相等；
[0100]
本发明所构建的解码器函数如下：
[0101]
decode(e)＝σ(deconv(upsampling(σ(deconv(upsampling(e))))))
[0102]
其中，σ为激活函数，deconv(e)为反卷积函数，e为编码器的输出，upsampling()为上采样函数；
[0103]
步骤3.2)使用随机梯度下降(sgd)优化器，sgd常用于分类任务，损失函数采用均方误差(mse)；除此之外，其余超参数batch_size设置为128，学习率为0.0001，迭代次数设置为20。
[0104]
步骤3.3)在微调阶段，提取步骤3.2)中编码器模型与权重参数，将序列和结构的编码器并行，增加两层lstm层，模型结构如图5所示；使用步骤2.3)中得到的数据集和标签，训练集与验证集的比例是8:2，使用独立测试集；最终实现结合位点预测；最终实现结合位点预测。
[0105]
步骤4)分别提取卷积自编码器中序列和结构编码器中卷积核的参数，分别与rna序列和结构的编码矩阵进行卷积运算，针对每个rbp，将所有滑动窗口中超过规定阈值的短序列拼接起来，使用工具weblogo 3对每组短序列进行序列比对，统计序列中每个位置上的一致性，作为motif。
[0106]
motif的提取与展示过程如下：
[0107]
步骤4.1)在卷积自编码器的微调阶段，模型在多次迭代之后，权重都被记录在卷积核中，如图3(a)和3(b)所示，提取第一层卷积核的参数，其中包含16个卷积核；
[0108]
步骤4.2)对于起点为i，将每个卷积核与步骤1)得到的序列矩阵在每个位点进行卷积运算，以长度为10的滑动窗口前进，对于每个碱基位点，可以得到中间值a
ij
，其中j表示是由第j个卷积核完成卷积操作的。对每个卷积核都设定了阈值threshold
j
，即该卷积核是所有短序列计算得到的中间值a
ij
的最大值的0.5倍，计算如下：
[0109][0110]
若a
ij
超过该阈值，则认为该起点为i，长度为10的短序列是有效的。由于motif长度一般为7，因此截取每条短序列的前7位；通过对每个卷积核滑动运算，可以得到16组对应的短序列集；
[0111]
步骤4.3)将16组短序列集拼接起来，得到一组含有共性的短序列集；
[0112]
步骤4.4)使用工具weblogo 3每个rbp的一组短序列集进行序列比对，统计出每个rbp的序列中每个位置上的一致性；
[0113]
步骤4.5)通过步骤4.4)所计算出的一致性，绘制序列标识图；序列标识图中每个碱基图形字符的大小与残基在该位置上出现的频率成正比，但纵轴数值一般不使用频率值，而是进一步计算后的结果。对于每个位置，每种碱基图形字符的高度height具体计算方式如下：
[0114]
height＝f
q
*r
[0115]
其中，q＝{a，u，g，c}，f
q
表示该种碱基的出现频率；r表示特定位点存在的信息量，
定义为最大可能熵和碱基分布熵之间的差值，计算如下：
[0116]
r＝log2(m)
‑
(entropy e
k
)
[0117]
entropy表示该位置的整体熵值，m表示碱基种类个数，不失一般性，对于蛋白质而言，m＝20，而对于本文的研究对象rna而言，m＝4。同时，计算中加入了校正因子e
k
，具体为小样本检验修正值，其中k表示用于统计的序列个数，计算如下：
[0118][0119]
e
n
＝(m
‑
1)/(2*ln2*n)
[0120]
对于结构motif，采用类似的处理方式也可以获取到。
[0121]
步骤5)对于每种rbp数据都训练一个对应的预测模型；对于一个待预测的rna序列，首先编码序列，然后将编码信息输入到已训练好的各个预测模型中，最终的结果即为特异rna与各个rbp对应的结合概率。

技术特征：
1.一种基于深度学习的rbp结合位点预测算法，其特征在于，包括如下步骤：步骤1)对于给定的rbp数据集，首先提取序列信息和独立测试集；再根据序列信息，预测出二级结构信息；根据序列和结构信息，分别构建序列和结构的编码矩阵；并构建独立测试集的编码矩阵，对应于每个rbp创建正样本和负样本；步骤2)构建包括生成网络和判别网络的生成对抗网络；在训练过程中，采用强化学习对判别网络优化，使优化过后的判别网络促使生成网络生成更逼真的合成数据，并在步骤1)中所得的序列和结构的编码矩阵上分别进行训练；步骤3)构建卷积自编码器预测模型，将步骤1)中所得到的编码矩阵和步骤2)中所得到的合成数据合并，分别得到数据规模更大的序列和结构矩阵，将序列和结构矩阵分别训练；步骤4)分别提取卷积自编码器中序列和结构编码器中卷积核的参数，分别与rna序列和结构的编码矩阵进行卷积运算，针对每个rbp，将所有滑动窗口中超过规定阈值的短序列拼接起来，对每组短序列进行序列比对，统计序列中每个位置上的一致性，作为motif；步骤5)对每种rbp数据均训练一个对应的预测模型；对于一个待预测的rna序列，首先编码序列，然后将编码信息输入到已训练好的各个预测模型中，最终的结果即为特异rna与各个rbp对应的结合概率。2.根据权利要求1所述的基于深度学习的rbp结合位点预测算法，其特征在于，所述步骤1)中构建序列和结构的编码矩阵过程如下：步骤1.1)设定rna序列长度为l，对rna序列和标签进行one
‑
hot编码；对于长度不足l的序列，采用“n”在末尾填充，“n”为等概率填充，“n”的向量表示为[0.25,0.25,0.25,0.25]；对于超过l的序列，将截取长度l作为训练序列；对于碱基a、c、g、u序列编码后编码矩阵维度为l
×
4；步骤1.2)rna二级结构为由focusfold算法预测得到的“点
‑
括号”格式，一对左右“括号”表示存在一对碱基配对，而“点”则表示未发生配对，每个符号分别按顺序对应于每个碱基位点，然后根据碱基配对原则，得到每个碱基位点的二级结构注释，所述二级结构注释包括堆叠s、自由端f、连接处j、发卡环h、内环与多环m，对rna结构进行one
‑
hot编码，结构编码后特编码阵维度为l
×
5。3.根据权利要求2所述的基于深度学习的rbp结合位点预测算法，其特征在于，所述步骤2)中生成合成数据的过程如下：步骤2.1)所述生成对抗网络包括生成网络和判别网络，其中，生成网络包括一层循环神经网络和一层全连接网络，全连接网络作为输出单元；判别网络包括一层循环神经网络、一层全连接网络和分类器；步骤2.2)将步骤1.1)和步骤1.2)中所得的序列和结构的编码矩阵分别输入到生成对抗网络训练，采用强化学习技术优化判别网络；被优化的判别网络的返回值反馈给生成网络，多次训练迭代，生成网络生成更逼真的合成数据；步骤2.3)提取步骤2.2)中已训练好的序列和结构生成网络，分别生成合成的序列和结构矩阵。4.根据权利要求3所述的基于深度学习的rbp结合位点预测算法，其特征在于，所述生成网络中，循环神经网络采用双向lstm，包含256个计算单元；全连接网络作为输出层，包含100个神经元；
所述判别网络中，循环神经网络采用双向lstm，包含256个计算单元；全连接网络作为中间层，包含100个神经元；分类器采用sigmoid激活函数。5.根据权利要求1所述的基于深度学习的rbp结合位点预测算法，其特征在于，所述步骤3)中卷积自编码器构建与训练过程如下：步骤3.1)所述卷积自编码器包括编码器与解码器，在预训练阶段，编码器与解码器共同工作，序列和结构分别训练，输入数据经过编码器的降维学习，得到中间抽象特征，之后解码器重构输入数据；步骤3.2)在微调阶段，提取编码器模型与权重参数，将序列和结构的编码器并行，增加两层lstm层，最终实现结合位点预测。6.根据权利要求5所述的基于深度学习的rbp结合位点预测算法，其特征在于，所述编码器包括两层卷积与两层池化；序列编码器的第一层卷积核大小为(4,10)，数量为16，步长为(1,1)，卷积核的长度为motif的1.5倍；在卷积层后面加入batch normalization，将每一层归一化到同一个正态分布，将其转化为均值为0、方差为1的标准正态分布；批归一化的函数运算如下所示:其中，x
l
是mini
‑
batch中的每个样本，将所有的x
l
归一化为标准正态分布，其中涉及到均值μ和方差δ两个参数；γ和β为可学习重构参数；n为batch
‑
size的大小，即mini
‑
batch中的样本量；y
l
为归一化后的数据，作为下一层的输入数据；所述池化层选用最大池化操作，池化窗口大小为(1,3)，池化步长为(1,2)，最大池化操作保留特征窗口内的最大值；在池化层后添加dropout层，对于保留下来的神经元，得到一个节点较少，规模较小的网络进行训练，以达到避免过拟合的目的；输入层神经元相当于第0层，即a
(0)
，从计算公式上来看，对a
(z)
层加入dropout，则有：a
(z)
＝r
(z)
*a
(z)
其中，r
(z)
表示神经元随机失活的概率分布，其服从于p＝1
‑
dropout的伯努利分布，a
(z)
表示加入了dropout的层；在编码器中，dropout层后继续添加卷积层与batch normalization，卷积核维度为(1,10)，数量为16，步长与上一层卷积层保持一致为(1,1)；构建的编码器函数如下：e＝encode(x)＝maxpooling(σ(conv(maxpooling(σ(conv(x))))))其中，σ为激活函数，conv(x)为卷积函数，x为输入变量，maxpooling()为最大池化函数；所述解码器包含两层反卷积函数和两层上采样函数，第一层上采样层大小为(1,2)，反
卷积层的卷积核维度为(1,12)，数量为16，步长为(1,1)，第二层上采样层大小为(1,2)，接下来反卷积层的卷积核维度为(4,10)，数量为1，步长为(1,1)；构建的解码器函数如下：decode(e)＝σ(deconv(upsampling(σ(deconv(upsampling(e))))))其中，σ为激活函数，deconv(e)为反卷积函数，e为编码器的输出，upsampling()为上采样函数。7.根据权利要求1所述的基于深度学习的rbp结合位点预测算法，其特征在于，所述步骤4)中motif的提取与展示过程如下：步骤4.1)在卷积自编码器的微调阶段，模型在多次迭代之后，权重被记录在卷积核中，提取第一层卷积核的参数；步骤4.2)对于起点为i，将每个卷积核与步骤1)得到的序列信息在每个位点进行卷积运算，以长度为10的滑动窗口前进，对于每个碱基位点，得到中间值a
ij
，其中j表示是由第j个卷积核完成的卷积操作；对每个卷积核设定阈值threshold
j
，即该卷积核是所有短序列计算得到的中间值a
ij
的最大值的0.5倍，计算如下：若a
ij
超过该阈值，则认为该起点为i，长度为10的短序列是有效的，将该序列加入motif序列集；通过对每个卷积核滑动运算，得到16组对应的短序列集；步骤4.3)将16组短序列集拼接起来，得到一组含有共性的短序列集；步骤4.4)使用工具weblogo 3将每个rbp的一组短序列集进行序列比对，统计出每个rbp的序列中每个位置上的一致性；步骤4.5)通过步骤4.4)所计算出的一致性，绘制序列标识图。8.根据权利要求7所述的基于深度学习的rbp结合位点预测算法，其特征在于，所述序列标识图中每个碱基图形字符的大小与残基在该位置上出现的频率成正比，但纵轴数值不使用频率值，而是计算后的结果；对于每个位置，每种碱基图形字符的高度height具体计算方式如下：height＝f
q
*r其中，q＝{a，u，g，c}，f
q
表示该种碱基的出现频率；r表示特定位点存在的信息量，定义为最大可能熵和碱基分布熵之间的差异，计算如下：r＝log2(m)
‑
(entropy e
k
)其中，entropy表示该位置的整体熵值；m表示碱基种类个数；e
k
表示校正因子，具体为小样本检验修正值，计算如下；e
k
＝(m
‑
1)/(2*ln2*k)其中，k表示用于统计的序列个数。
技术总结
本发明公开了一种基于深度学习的RBP结合位点预测算法，首先从RNA序列信息出发，增加二级结构特征作为输入；通过基于强化学习的生成对抗网络，对RBP进行数据增强，目的是通过增强训练集规模来提高预测模型的预测能力；提出卷积自编码器作为预测模型，预训练过程包括编码器和解码器的训练，提取已训练好的编码器微调；能够预测未知结合位点的RNA序列；通过卷积自编码器的卷积核，能够提取与RNA发生结合的序列motif与结构motif，探究不同RBP的序列和结构特异性；对于每种RBP数据都训练一个对应的预测模型；对于一个待预测的RNA序列，首先编码序列，然后将编码信息输入到已训练好的各个预测模型中，最终的结果即为特异RNA与各个RBP对应的结合概率。对应的结合概率。对应的结合概率。


技术研发人员：朱敏 闫建荣 明章强 王心翌
受保护的技术使用者：四川大学
技术研发日：2021.03.02
技术公布日：2021/6/24