信号肽预测之Signal-LMS法

引 言

蛋白质是生理功能的执行者及生命现象的体现者。随着后基因时代的来临，对于蛋白质功能的研究成为了生物信息学研究的核心内容[1]。从蛋白质序列出发可以预测与蛋白质功能相关的特征，信号肽则是其中重要的特征之一[2]。由于蛋白质数据库中的序列数量呈指数增长，蛋白质信号肽预测的任务越来越迫切，开发有效的方法势在必行[3～5]。信号肽位于分泌蛋白的 N 端，一般由 15～30 个氨基酸组成。它包括三个区域：一个带正电的 N 末端、一个中间疏水序列和一个较长的带负电荷的 C 末端[6]。蛋白质信号肽预测问题主要可以分为两个步骤。第一步：判断蛋白质是否含有信号肽部分，即是否为分泌蛋白；第二步：若该蛋白质为分泌蛋白，则预测其信号肽的剪切位点[7]。对于一个未知的蛋白质序列，第一步分类的结果是第二步预测的基础，它在一定程度上影响第二步预测的准确性。因此，本文主要针对第一步 (蛋白质信号肽分类预测) 提出新的方法。目前，为了避免信号肽序列长度差异所带来的问题，大多数算法[8～18]采取用滑动窗口将整个序列切割成长度相等的氨基酸序列的方式进行处理，如经典的权矩阵方法、子位点耦合模型、基于支持向量机和隐马尔科夫模型的方法，以及被广泛使用的 SignalP 方法等。而后，序列比对的方法[7,19,20]被应用到信号肽预测中，如经典的Needleman-Wunsch 全序列比对算法[21]，避免了借助滑动窗口带来的问题。全局比对是对给定的序列全长进行比较，这种算法适用于全局水平上相似性程度较高的两个序列[22]。然而，信号肽是蛋白质的功能位点之一，突出的是蛋白质序列局部片段所体现的生物特性。因此，局部比对往往比全局比对具有更高的灵敏度，其结果更具有生物意义[23～25]。在生物信息学领域，DNA 或蛋白质等分子序列存在局部的保守功能片段，这些功能片段作为整个序列的关键特征，可以用来度量序列间的相似性水平[26,27]。所以，在蛋白质信号肽预测问题上使用局部比对是较为合理的，它可能会揭示一些匹配的功能性序列片段，而全序列比对则很可能会被一些完全不相关的残基所淹没。因此，本文提出一种基于局部序列匹配相似度度量的方法来预测信号肽，首先采用了氨基酸疏水特性来编码蛋白质，然后搜索序列间的局部匹配子序列，并根据计分矩阵BLOSUM62来度量两个序列间的相似性，最后采用k最近邻算法来判断是否为信号肽。

材料与方法

在生物信息学研究中，我们的任务是从大量的生物数据中提取有价值的信息，而序列比对是最重要、最常用的研究手段。所谓序列比对，是指两个或多个序列按字母比较，反映它们之间的相似或相异性，从而发现与结构相关联的保守序列片段，或者从已知序列预测新序列的结构和功能。本文提出Signal-LMS 方法，其核心思想是利用氨基酸疏水特性来编码蛋白质，从而有效地搜索具有相似特性的局部匹配子序列；然后根据计分矩阵BLOSUM62 计算所有匹配子序列中氨基酸的相似性来反映序列的整体相似性；最后，在k最近邻算法的思想下预测信号肽。一般地，若两条序列具有越多、越长的匹配子序列，并且匹配子序列中的氨基酸相似性越高，则说明它们具有越高的整体序列相似性，具有相似的功能位点。

编码方案

蛋白质由 20 种氨基酸组成，每种氨基酸可以用字符表中的一个字母来表示。所以，我们一般将蛋白质序列看作有限字符集∑上的字符串。∑={A,C,D,E,F,G,H,I,K,L,M,N,P,Q,R,S,T,V,W,Y}组成蛋白质的 20 种氨基酸具有不同的特性，如酸碱性、亲水性、疏水性等。疏水性是氨基酸基本的理化性质之一，20 种氨基酸的相对疏水值见表 1。氨基酸的疏水性，是影响氨基酸溶解行为的重要因素，也是影响蛋白质和肽物理化学性质的重要因素[28]。因此，本文按照氨基酸的相对疏水值进行编码，疏水值大于0的为疏水性氨基酸(hydrophobic)，疏水值在－3.2～0的为中性氨基酸(neutral)，疏水值小于－3.2的为极性氨基酸(polar)。根据上述规则，可以将20 种氨基酸分为三类 (见表 1)：polar、neutral 和 hydrophobic，简写为：P、N 和 H。根据每个氨基酸的疏水性，可以将蛋白质的编码进行如下转换，例如一条蛋白质 p 的一级结构(由 20 种氨基酸组成) 表示为：MELSGATMARGLAVLLVLFLHIKNLPAQAA经疏水值编码后为：HPHNNHNHHPNHHHHHHHHHNHPPHNHPHH

相似度度量

考虑到蛋白质功能位点由较短的序列片段组成，尽管在序列的其它部位可能有插入、删除或突变，但是这些序列片段却具有相当大的保守性。因此，在信号肽预测问题上，我们采用局部序列匹配的方法来度量两个蛋白质序列间的相似度。相关的定义如下：定义1：如果 S 是一个序列，那么 |S| 表示 S 的字符长度，S[i]表示序列的第 i 个字符。如果序列S和T匹配，必须满足以下条件：①|S|=|T|；②S[i]=T[i]，(0<i≤|S|)定义2：如果 x 和 y 是两个字符，那么 b(x,y)表示 x 和 y 字符在进行比较时所得的分值，b为一个记分矩阵。定义 3：如果 S 和 T 是两个序列，那么 Ms,t表示序列S和T的一个匹配子序列，其中s奂S, t奂T。SlS,T表示所有Ms,t的集合，且满足以下条件：①|Ms,t|≥l②(坌Ms',t'∈SlS,T)∧(Ms',t'≠Ms,t)∧(|Ms',t'|＜|Ms,t|)圯(s'埭s)∨(t'埭t)条件①表示匹配子序列的长度应不小于最小匹配长度 l，条件②表示匹配子序列集合中不存在相互重叠的序列，即每个匹配子序列都应尽量长，其子序列不能作为新的匹配子序列。假设有两个蛋白质序列X=(x1,x2,…,xn)与Y=(y1,y2,…,ym)，X和Y分别由 20 种氨基酸字符组成，那么这两个序列的相似性得分计算步骤如下：1)将两个蛋白质序列X和Y根据氨基酸的相对疏水值编码后得到X'=(x'1,x'2,…,x'n)和Y'=(y'1,y'2,…,y'm)，X' 和 Y' 由氨基酸的 3 种疏水性 H、P 和 N 组成。2)根据下面给出的SMS (search matching subsequences)算法获得序列X'和Y'的所有匹配子序列集合SlX',Y'={M(1)x',y',M(2)x',y',…,M(h)x',y'}，h 为匹配子序列的个数。Algorithm 1. SMSInput: Two sequences X, Y and the minimum length of matching subsequence LOutput: Matching subsequences set SBegin1. For(i=1 to |X|－L＋1) do2. k=13. While(k<=|Y|－L＋1)4. If(X[i…i＋L－1]==Y[k…k＋L－1] and X[i－1]~=Y[k－1])5. M=X[i…i＋L－1]6. m=i＋l, n=k＋l7. While(m<=|X| and n<=|Y| and X[m]==Y[n])8. Add X[m] to M9. m=m＋1, n=n＋110. EndWhile11. Add M to S12. Clear M13. EndIf14. k=k＋115. EndWhile16. EndFor17. Return SEnd3)匹配子序列M(i)x',y'对应原序列X和Y中的两个序列，表示为M(i)x和M(i)y。4)计算两个序列的相似度Score。首先根据替换矩阵 BLOSUM62，计算 M(i)x和M(i)y两个序列相对应的氨基酸相似度之和，然后将所有匹配子序列的相似度累加，得到两个序列整体的相似度。其中，b表示替换矩阵，|M(i)x',y'| 表示匹配子序列的长度，h 表示匹配子序列的个数。

k 最近邻算法预测信号肽

k最近邻算法是一种基于类比的机器学习方法，即通过将给定的未知元组与训练元组进行比较来学习，在模式分类中的应用非常广泛[29～31]。当给定一个未知元组时，应用 k 最近邻分类法搜索该模式空间，找出最接近未知元组的k个训练元组，通过这k个最近邻样本中的多数样本的类别来确定未知元组的类别。本文采用 k 最近邻算法思想来预测信号肽。假设一个蛋白质数据集S 包含 N 个蛋白质样本 p1,p2,…,pN，可以把S分成两个数据子集S＋和S－，S＋表示分泌蛋白数据集，S－表示非分泌蛋白数据集。对于一个未知类别的蛋白质样本p，按下面方法判断它属于哪个类。将待测样本p和数据集S中的N个样本组成N个序列对：用d1表示序列对 <p,p1>的相似度，d2表示序列对 <p,p2>的相似度，以此类推，我们可以得到p与N个蛋白质的相似度：根据k最近邻的思想，取出与p相似度最高的k个序列对，得到相应的k个类别：

结 果

数据集

为便于对比实验结果，本文选用文献[20]中对 Neilsen 公布的 SWISS-PROT (Version 50.7)进行筛选后的数据集来实验。该数据集包括分泌蛋白和非分泌蛋白的 N 终端部分，它的构造过程经过了严格的筛选，具有代表性。每个分泌蛋白包括信号肽部分和成熟蛋白质的前30个氨基酸残基，而每个非分泌蛋白包括了前 70 个氨基酸残基。本文采用了 6 组不同生物体的数据集。它们的组成如表 2。实验评价标准实验使用Matlab 编写仿真程序，采用的数据集为 Human、Plant、Animal、Eukaryotic、Gram-positive 和 Gram-negative。为了更客观地评价预测效果，本文采用马氏相关系数作为评价标准[19]：式中，Ntp和Ntn分别为正确预测的分泌蛋白和非分泌蛋白个数，Nfp为分泌蛋白预测为非分泌蛋白的个数，Nfn为非分泌蛋白预测为分泌蛋白的个数。实验运用 5 折交叉验证法验证预测结果。将每个数据集平均分成 5 份，取 1 份作为测试集、1 份作为检验集，其余 3 份作为训练集。训练集是已知样本的类别；检验集用于优化训练集，将训练集中影响预测效果的噪声样本删除；测试集用于计算预测精度。

实验参数确定

在序列间寻找匹配子序列时，最小的序列匹配长度L的设置是很重要的，它决定了匹配子序列的个数，从而影响到两个序列间的相似性度量值。若L的取值太小，则大大增加了序列间匹配子串的个数，同时增大了系统开销；反之，若L的取值太大，那么就会出现较多的序列对不存在匹配子串的情况，从而使得序列相似性度量效果差。因此，我们通过实验的方法对最小的序列匹配长度L进行确定。从数据集中随机选取分泌蛋白和非分泌蛋白样本各50个，在 100 个样本中进行两两比对，可以得到约 5000 个序列对。实验中，最小的序列匹配长度L分别取 4、5、6 和 7，并对所有的序列对进行匹配子串个数统计，实验结果见图 1。图 1 中，实线表示不同类序列对的统计结果，虚线表示同类序列对的统计结果。从图 1 可以看出，在总体情况下，不同类序列对的匹配子串个数比同类序列对要少。当L取值为4 时，匹配子串个数大部分集中在 25～35 个，在计算相似度时需要较大的系统开销；当L取 7 时，大部分的序列对不存在匹配的子串，使得大多数序列对的相似度为 0，不利于分类效果。当 L 取 5 和 6 时，匹配子串个数分别集中在 10～15 和 3～5，考虑到 L=6时部分不同类序列对匹配子串的个数为 0，即它们的相似度为 0，使得类别之间的区分性更强，而L=5时没有得到匹配子串数为 0 的序列对。因此，本文实验中最小的序列匹配长度L取值为 6。

实验结果

本文进行比较的算法有 PrediSi[32]、Signal-3L[20]和全序列比对方法。其中，PrediSi 和Signal-3L 方法是近年来预测蛋白质信号肽及其剪切位点较为流行的方法，我们选用 PrediSi和 Signal-3L 中预测信号肽部分的算法与本文的方法进行比较，且它们采用了与本文相同的数据集，具有一定的可比性。全序列比对方法则采用了经典的 Needleman-Wunsch 算法作为序列比对方法，其余方法同本文的算法。比较结果见表 3。实验中，全序列比对方法和本文的方法均采用了k最近邻的方法进行分类，其中，k值的选择很重要，它直接影响分类的结果。若 k 值选择过小，得到的近邻数过少，会放大噪声数据的干扰；而k值选择过大，若待测样本属于训练集中数据较少的类，那么近邻中包含了更多不相似的样本，就会造成噪声增加而导致分类效果降低。因此，本文将 k 从 10 到100、以步长为 10 的方式取值，即 k 分别取 10、20、30、…、100，然后得到平均准确率。

讨 论

本文提出的 Signal-LMS 方法是基于局部序列匹配相似度的方法来预测信号肽，与其进行比较的 PrediSi 方法是通过统计蛋白质序列中氨基酸出现频率来构造权矩阵，而 Signal-3L方法是通过建立 OET-KNN 多分类器来预测信号肽。这两种方法均没有考虑到蛋白质序列间的相似性，然而相似的蛋白质序列往往能体现相同的生物特性。从表 3 中可以看出，除了 Gram-positive 数据集，在其他数据集上采用序列比对的方法，均比 PrediSi 和 Signal-3L方法的准确率要高。这说明，将序列比对相似性的方法应用于信号肽预测问题具有较好的效果。在平均准确率上，本文提出的 Signal-LMS 方法要略高于其他方法，而在 Plant、Animal和 Gram-negative 数据集上，效果较为突出。但是，在 Human 和 Eukaryotic 数据上，全序列比对的方法比本文的方法准确率高。为了更好地比较本文方法和全序列比对方法，我们做进一步讨论。本文的 Signal-LMS方法和全序列比对方法均采用了序列间相似性比对的算法，但相似性度量方法不同。Signal-LMS 方法强调的是序列间局部相似性，反映序列中的局部特性，而全序列比对方法则是从序列的全局特性考虑。因此，为了比较这两种相似性度量方法，我们在采用 k 最近邻思想进行分类时，对取不同k值而得到的准确率进行比较，结果见图2。从图2 看出，随着 k 值的增大，Signal-LMS 方法预测的准确率相对稳定，而全序列比对方法随着k值的增大，预测的准确率逐渐降低。这说明，随着近邻数的增大，全序列比对方法待测样本的近邻中包含更多与其不同类的样本，反映出该相似性度量方法对不同类样本区分不明显。而本文的局部序列匹配相似性的度量方法能够更准确地找到具有信号肽特征的同类蛋白质，体现了蛋白质序列的局部保守性。这说明全序列比对的结果不能很好地体现序列的相似性，会出现不同类蛋白质比同类蛋白质相似性高的情况，造成错误预测，而采用局部序列比对比全局序列比对更能突显序列的相似性。因此，在信号肽预测问题上，Signal-LMS方法比全序列比对方法具有一定的优势。

总 结

在预测信号肽问题上，本文采用了一种基于氨基酸相对疏水值的编码方案，通过搜索序列间的局部匹配子序列，并根据计分矩阵 BLOSUM62 来度量两个序列间的相似性，最后采用k最近邻算法进行分类。本文方法有效地避免了以往采用滑动窗口所引起的数据不平衡等问题，并且反映了蛋白质序列的局部保守性。在 SwissProt 数据集上对 6 组生物体数据进行实验分析，结果表明本文的Signal-LMS 方法是一种高效的信号肽预测方法。