基于深度学习的蛋白质卷曲螺旋结构特征预测实现方法与流程

1.本发明涉及的是一种基因工程领域的技术，具体是一种基于深度学习的蛋白质卷曲螺旋结构特征预测(cocopred)实现方法。


背景技术：

2.卷曲螺旋是蛋白质中的一种结构模式，有两个至多个α螺旋相互缠绕构成。普通的α螺旋的螺旋周期是7.2，即螺旋一周包含3.6个氨基酸。但是卷曲螺旋中螺旋的周期是7，这会导致螺旋具有一种超螺旋结构，使得多个螺旋缠绕在一起。三级结构中的卷曲螺旋对应到一维序列上的区域被称为卷曲螺旋区域。卷曲螺旋区域由重复的七元组构成，通常表示为[abcdefg]
n
，a
‑
g被称为七元标记位置。在缠绕的过程中，螺旋会以一种“把手进入洞”的方式结合在一起，即一个螺旋的几条侧链会插入到对向螺旋侧链形成的空洞之中。这种缠绕方式会大大提高卷曲螺旋的稳定性。构成卷曲螺旋的螺旋数目被称为寡聚体状态，据此可以将卷曲螺旋分为顺式二聚体，反式二聚体，三聚体，以及四聚体。
[0003]
根据领域内的报道，卷曲螺旋在很多生物过程中发挥重要的作用。准确预测卷曲螺旋的结构特征对于理解其功能具有重要的意义。尽管领域内已经出现大量的研究工作，卷曲螺旋特征预测仍然有改进的空间。首先，大多数工作只聚焦于某一种特征，很少有模型对卷曲螺旋的三种特征(卷曲螺旋区域、七元标记位置、寡聚体状态)同时建模。当要想对蛋白质序列做全面的预测，需要组合使用多种模型才能实现。其次，在寡聚体状态预测问题上，大多数模型只能区分有限的状态，只有logicoil可以区分全部四种寡聚体状态。最后，随着近年深度学习的广泛运用，以及蛋白质结构的快速累积，开发准确高效的卷曲螺旋特征预测模型成为可能。


技术实现要素：

[0004]
本发明针对现有技术存在的上述不足，提出一种基于深度学习的蛋白质卷曲螺旋结构特征预测实现方法，采用自注意力机制使得模型自动关注卷曲螺旋中对预测至关重要的位置，并采用多类别交叉熵作为损失函数。模型最终的损失函数是三个网络损失函数的加权，以及网络参数的正则项。
[0005]
本发明是通过以下技术方案实现的：
[0006]
本发明涉及一种基于深度学习的蛋白质卷曲螺旋结构特征预测实现方法，使用两层卷积网络和两层双向长短时记忆网络分别提取序列的局部特征以及全局特征，经拼接后通过残差网络预处理，再分别经卷曲螺旋区域预测网络、七元标记位置预测网络以及寡聚体状态预测网络预测得到卷曲螺旋的卷曲螺旋区域(coiled
‑
coil domain)、七元标记位置(register)以及寡聚体状态(oligomeric state)。
[0007]
所述的卷曲螺旋区域预测网络、七元标记位置预测网络以及寡聚体状态预测网络均包含卷积层和全连接层。技术效果
[0008]
本发明针对卷曲螺旋特征预测，使用卷积网络和双向长短时记忆网络，提取序列的局部特征和全局特征，使用注意力机制预测寡聚体状态，以及使用条件随机场预测七元标记位置。采用条件随机场刻画七元标记位置之间存在的先后顺序关系。同时考虑到卷曲螺旋区域的长度不一样，在输出至条件随机场之前会对张量进行切割，去除非卷曲螺旋区域的影响。对于寡聚体状态预测网络，本发明使用注意力机制，使得模型自动的关注对预测更加重要的区域。注意力机制也很好的解决卷曲螺旋区域长度不一致问题，使模型对任意长度的卷曲螺旋都能预测一个固定长度的向量。另外七元标记位置网络的输出会作为特征之一输出至寡聚体状态预测网络，提升模型的预测能力。
附图说明
[0009]
图1为实施例流程图；
[0010]
图2为实施例方法与logicoil的比较示意图；
[0011]
图3为实施例中七元标记位置混淆矩阵示意图；
[0012]
图4为实施例中七元标记位置转移矩阵示意图；
[0013]
图5为实施例卷曲螺旋预测的混淆矩阵示意图；
[0014]
图6和图7为实施例正确预测的卷曲螺旋和错误预测的卷曲螺旋对应的频率分布示意图。
具体实施方式
[0015]
如图1所示，为本实施例涉及一种基于深度学习的蛋白质卷曲螺旋结构特征预测实现方法，包括以下步骤：
[0016]
s1、从领域内工作提供的数据集、卷曲螺旋数据库cc 中整理出训练集、验证集和测试集，具体包括：
[0017]
s11、领域内已有科研人员整理测试数据，用来评价目前广泛使用的几种卷曲螺旋区域预测软件。测试数据包括718条蛋白质，其中279条不包含卷曲螺旋，439条包含卷曲螺旋；
[0018]
s12、从公开的卷曲螺旋数据库cc 中获取卷曲螺旋蛋白质，具体为：首先下载所有包含至少一个卷曲螺旋的蛋白质。然后剔除满足以下任一条件的蛋白质：(i)分辨率太差(ii)序列长度小于25或者大于700；(iii)序列中包含20种常见氨基酸之外的氨基酸；(iv)序列中的卷曲螺旋长度小于8；(v)序列中存在参与多个卷曲螺旋的氨基酸；(vi)序列中的卷曲螺旋的寡聚体状态不属于以下四种类别：顺式二聚体，反式二聚体，三聚体，以及四聚体。在剩余的序列中，去除与测试序列相似度大于30％的蛋白质。再对自身去冗余至50％，即得到本实施例的训练正样本，共2327条序列。
[0019]
s13、将deepcoil模型训练集中来自于pdb数据库的负样本补充到训练集中，共9358条序列。
[0020]
s2、使用hhblits工具对整理数据集中的序列提取hmm特征；对数据集中的序列做独热(one
‑
hot)编码；使用eisenberg scale作为序列的疏水指标特征。
[0021]
所述的独热编码作为表征蛋白质序列最直接的一种方式，具体是指：将序列中的氨基酸编码成长度为20的向量，其中只有一个元素值为1，其他元素值等于0。对于一个长度
为l的蛋白质序列来说，独热编码的维度为l
ⅹ
20。
[0022]
所述的hmm特征来源于隐马尔可夫概率模型，本实施例使用hhblits序列比对工具搜索uniclust30数据库得到hmm特征，执行命令及参数为：hhblits
–
i sequence.fasta
‑
n 3
‑
e 0.001
‑
d uniclust30_2017_10
‑
cpu 6
‑
ohhm sequence.hmm
‑
diffinf
‑
id 99
‑
cov 50。hmm特征有30维。
[0023]
所述的疏水指标特征，即氨基酸侧链的疏水程度，本实施例使用eisenberg scale疏水指数。对于一个氨基酸残基来说，疏水指数特征共1维。
[0024]
s3、使用两层卷积网络和两层双向长短时记忆网络分别提取序列的局部特征以及全局特征，在特征维度上拼接，之后输出至一层残差网络中得到更有效的特征，最终将该特征进一步输入卷曲螺旋区域预测网络、七元标记位置网络以及寡聚体状态预测网络，分别得到卷曲螺旋的三种结构特征。
[0025]
本实施例通过将训练集等分成五份，使用5折交叉验证框架使用同一套参数分别在五个模型各自的训练数据上训练，根据模型在验证集上的性能对模型进行参数调整。
[0026]
所述的卷曲螺旋区域预测网络包括：一层卷积层和一层全连接层，该网络使用交叉熵作为损失函数。
[0027]
所述的七元标记位置网络包括：一层卷积层和一层全连接层，其中：全连接层的输出根据卷曲螺旋区域进行切割；该网络使用条件随机场作为损失函数，即切割结果输出至条件随机场中得到损失函数；在预测序列的七元标记位置时，该网络使用维特比算法得到最优标记位置序列。
[0028]
由于七元标记位置之间存在某种关系，例如a后面往往紧跟着b，因此完全独立地预测每个标记是不合适的。本实施例使用条件随机场刻画标记之间的关系，使得一个卷曲螺旋上七元标记位置之间的预测被绑定在一起。
[0029]
所述的七元标记位置网络的输出采用梯度停止技术，这样寡聚体状态预测网络的训练就不会直接影响到七元标记位置网络。
[0030]
所述的寡聚体状态预测网络使用注意力机制，其输入为七元标记位置网络中全连接层的输出，经残差层拼接得到，具体包括：一层卷积层和一层全连接层，其中：卷积层的输出根据卷曲螺旋区域进行切割，然后在特征维度上拆分成四份，分别输出至四个注意力头上，全连接层将全部注意力头的输出融合，得到最终结果，该网络采用多分类的交叉熵作为损失函数。
[0031]
本实施例的最终损失函数为三个子网络损失函数的加权和，以及网络中可训练参数的正则项。
[0032]
s4、针对不包含任何标记的蛋白质序列进行预测，得到其是否包含卷曲螺旋、卷曲螺旋的范围、七元标记位置以及寡聚体状态，或者针对已知卷曲螺旋区域预测，得到其七元标记位置和寡聚体状态。
[0033]
所述的针对不包含任何标记的蛋白质序列进行预测是指：分别使用hhblits软件对序列数据库uniclust30进行搜索得到hmm特征、根据序列得到独热(one
‑
hot)编码以及eisenberg scale疏水指标，将这些特征输入模型预测出序列是否包含卷曲螺旋以及卷曲螺旋对应的区域，进一步预测每个卷曲螺旋区域的寡聚体状态以及七元标记位置。
[0034]
所述的针对已知卷曲螺旋区域预测是指：直接在标注的区域信息基础上预测每个
卷曲螺旋的寡聚体状态以及七元标记位置，当正确预测出卷曲螺旋区域时，很可能正确预测其寡聚体状态，因此是否正确预测卷曲螺旋区域一定程度上表明本实施例能否正确预测对应的寡聚体状态。
[0035]
训练过程中的batch size为64，即64条序列。其中32条包含卷曲螺旋，32条不包含卷曲螺旋。由于四种寡聚体状态的卷曲螺旋数量极度不平衡，本实施例还对部分蛋白质进行复制，缓解样本不平衡问题。
[0036]
所述的训练，模型的学习率设为5e
‑
4。预测卷曲螺旋的最短长度为6。模型使用指数线性单元elu作为非线性激活函数，使用层归一化调整神经元的输出。为防止过拟合，模型还在长短时记忆网络和卷积网络上使用dropout层。
[0037]
除此之外，为保证模型训练的稳定性，本实施例采用梯度裁剪技术防止出现梯度爆炸问题。模型使用tensorflow框架，在训练集上训练100轮，参数根据模型在5折交叉验证框架下验证集上的性能进行调整。
[0038]
本实施例对上述预测得到的结果进行的评价指标分为两类，一类在氨基酸层面上做评价，一类在片段层面上做评价。一个卷曲螺旋区域被认为是正确预测，当且仅当预测的片段和真正的片段重叠区域的长度大于两者最长长度的一半以上。除此之外，本实施例还采用准确率召回率曲线下面积，片段重叠度，以及混淆矩阵等评价指标。
[0039]
所述的评价指标包括：
[0040][0041][0042][0043][0044][0045][0046]
将本实施例提出的算法(cocopred)与领域中已有的卷曲螺旋区域预测算法在测试集上进行比较，对比结果如表1所示。由于deepcoil算法只能对长度在500以下的蛋白质做预测，因此只用长度小于500的测试蛋白评测本实施例与deepcoil，结果如表2所示。本实施例提出的算法在几乎所有指标上都明显优于领域内其他算法。
[0047]
表1.本实施例与领域内已有算法在测试集上的效果对比
[0048]
表2.本实施例与deepcoil算法在测试集上的效果对比
[0049]
如图2所示，为本实施例(cocopred)与领域内最好的寡聚体状态预测算法logicoil的比较结果。为在较大规模的数据集上比较，保证结果的显著性，本实施例使用5折交叉验证中的验证集评测结果。可以看出本实施例在顺式二聚体，反式二聚体，以及三聚体上的性能都显著好于logicoil。在四聚体的预测上，本实施例受限于有限的数据量，没有达到较好的效果。未来随着蛋白质结构数量的不断增加，本实施例的效果也会进一步提升。
[0050]
如图3所示，为本实施例(cocopred)在七元标记位置预测任务上的混淆矩阵，评测数据为5折交叉验证框架下的验证集。相较于卷曲螺旋区域和寡聚体状态，七元标记位置的重要性稍小，领域内基本没有预测七元标记位置的软件。因此如图3仅所示，为本实施例的预测结果。可以看出，混淆矩阵中对角线上的数字是所在行所有数字中最大的，表明本实施例可以很好地预测所有类型的标记位置。
[0051]
如图4所示，为本实施例(cocopred)中七元标记位置网络的转移矩阵。转移矩阵中每个元素与对应标记位置之间的转移概率呈正相关关系。可以看出，a到b,b到c,c到d,d到e,e到f,f到g对应的转移概率是最大的，表明转移概率矩阵很好地抓住不同标记位置之间的顺序关系。
[0052]
根据数据集中的卷曲螺旋是否被cocopred正确预测出，卷曲螺旋区域可以被分成两类：正确预测的卷曲螺旋区域和错误预测的卷曲螺旋区域。如图5使用cocopred分别预测这两类卷曲螺旋的寡聚体状态，并绘制对应的混淆矩阵。矩阵的每个位置上有上下两个数字。下面的数字表示对应卷曲螺旋的数目，上面数字表示占数据集中对应寡聚体状态样本数目的比例。可以看出矩阵中每行比例之和等于1。从如图5可以发现，左边矩阵的对角线上的比例比右边矩阵的大不少，特别是对于顺式二聚体，三聚体，和四聚体。反式二聚体的比例基本持平。这表明对于正确预测的卷曲螺旋区域，cocopred有更高的概率正确预测其寡聚体状态。
[0053]
针对正确预测的卷曲螺旋区域和错误预测的卷曲螺旋区域，记录cocopred预测寡
聚体状态时的注意力向量。只保留状态预测正确的卷曲螺旋对应的注意力向量，同时对于一个卷曲螺旋，只保留最大值和次大值差别最大的注意力向量。最后对四种寡聚体状态的卷曲螺旋分别统计权重最大的氨基酸对应的标记位置的频率。如图6和如图7分别展示正确预测的卷曲螺旋和错误预测的卷曲螺旋对应的频率分布。可以看出，标记位置a，b，和e在如图6中的比例最大，表明对于正确预测的卷曲螺旋区域，这三种标记位置对于预测寡聚体状态至关重要。同时发现对于四种寡聚体状态，频率分布非常类似。如图7没有如图6中的规律明显，考虑到这些卷曲螺旋区域没有被cocopred正确识别，如图7中的分布可能没有正确反映出每种标记位置的重要性。
[0054]
本实施例使用三个子网络对卷曲螺旋的三种结构特征建模，能够区分顺式二聚体、反式二聚体、三聚体以及四聚体四种卷曲螺旋寡聚体状态。领域内的绝大部分工具只能区分其中的两种，例如二聚体和三聚体、顺式二聚体和三聚体。目前已知的能够区分全部四种状态的工具是logicoil。
[0055]
与现有logicoil技术相比，本发明在卷曲螺旋区域上的预测情况，能够反映其进一步预测寡聚体状态的性能，对于正确预测的卷曲螺旋区域，本方法正确预测其寡聚体状态的可能性更高，同时由于采用注意力机制，在卷曲螺旋区域上的预测结果能够给出额外的信息。实验表明，a,b,e三种标记位置在寡聚体状态预测上发挥更加重要的作用。
[0056]
上述具体实施可由本领域技术人员在不背离本发明原理和宗旨的前提下以不同的方式对其进行局部调整，本发明的保护范围以权利要求书为准且不由上述具体实施所限，在其范围内的各个实现方案均受本发明之约束。