一种配体小分子关键子结构筛选方法及装置与流程

本发明涉及一种药物筛选方法及装置，尤其是一种基于分子指纹和机器学习的配体小分子关键子结构筛选方法及装置，属于药物虚拟筛选技术领域。

背景技术：

在新药发现过程中，通过其他各种途径首先发现的药物不一定是这类化合物中最理想的药物，但却提出了一个新的化学结构，我们通常称它为先导化合物。在先导化合物的分子结构上进行各种化学修饰，研究化学结构与生物效应的关系，通过药理实验与数学计算，确定最佳结构从而选定药物，这已经成为近代寻找新药的一种基本手段，而且取得了显著成绩。如何快速有效的发现能够实现特定功能的关键子结构已成为加速新药研发的关键环节。

药物设计过程中小分子化合物数量庞大，由于实验条件、人力以及资金等方面的限制，无法对所有小分子化合物进行实验测试，虚拟筛选(virtualscreening,vs)作为计算机辅助药物设计中常用的一种方法应运而生。目的是在小分子数据库中筛选出适合目标分子上的小分子结构。计算机的迅猛发展也使得虚拟筛选技术得到更加广泛的应用。虚拟筛选包括以下两类：

(1)基于受体的虚拟筛选是以受体分子的三维构象为基础，借助分子对接等方法将数据库中小分子与受体分子结合，预测小分子与受体分子之间的结合情况，借助打分函数进一步分析小分子与受体分子的活性位点氨基酸的匹配程度和亲和力强弱，然后对结果进行排名挑选潜在小分子药物。

(2)基于配体的虚拟筛选主要是对现有活性小分子的构象、大小、药效团、理化性质和活性关系进行分析，常用方法主要有相似性搜索(similaritysearching)，药效团匹配(pharmacophorematching)、子结构搜索(substructuresearching)构建定量关系和计算新化合物的亲和力。

基于配体的筛选方法最关键的是使用各种方法或分子指纹对化合物进行充分的描述，以研究化合物是否具有类似活性或致病机理，或归纳得到对化合物活性起关键作用的一些基团的信息。在基于配体的药物虚拟筛选中，配体分子的分子指纹表征是关键。它通过检测分子结构中一些特定的子结构是否存在，从而将分子结构转化成一系列二进制的指纹序列，这启发我们小分子的结构相似性可以用分子指纹的相似性来表示。另一方面，药物小分子与药物靶标结合时，只有局部位置相结合，这就意味着只有少数关键子结构对于其生物活性起决定性作用。若能发现这些与特定靶标结合起到特定功能的关键子结构对于新药研发有着极大的推动作用。

在最简单的分子指纹形式中，指纹中的每一位都表示对应分子是否具有(如果有，则该位设为1)某个特定的特征。maccs指纹就是这种指纹的典型代表，它共由166个二进制位组成，其中每一位都对应一个特定的结构属性信息。显而易见，这种类型的指纹都不可避免的具有自身的局限性，即：只有在字典中已定义的子结构才有可能被识别。而面对海量的化合物分子，小规模的子结构识别远远不能满足我们的需求，目前最流行的是使用ecfp(extended-connectivityfingerprints)分子指纹方法来进行表征。不同的ecfp指纹方法，如ecfp4、ecfp8、ecfp12等，后面的数字代表配体包含子结构的半径大小，一般来说，半径越大，生成的特征越多。由于在药物虚拟筛选中，化合物分子的数量巨大，基于ecfp方法产生的特征维度也是巨大的。使用最新的ecfp12分子指纹对配体分子进行表征，每个配体分子的特征维数可以高达上千万维。

配体分子主要是通过药效基团与靶标分子进行作用，而药效基团通常与配体分子的少数子结构有关，即大多数子结构与配体分子的活性无关。如果把配体分子数据集转化为矩阵，其每行表示一个样本，对应于一个配体分子，而每列表示一个特征，对应于一个特定的子结构。对配体海量特征进行筛选时，需要考虑特征具有的“稀疏性”，即矩阵中的很多列与配体分子的活性是没有关系的，通过特征选择去除这些列，则在实际的学习任务仅需要在较小的矩阵上进行，学习任务的难度可能有所降低，涉及的计算和存储开销会减少，学得模型的可解析性也会提高。在机器学习领域，生成的配体分子海量特征在实际任务中会陷入“维数灾难”问题，运算一次时间开销过大甚至无法拟合数据。

技术实现要素：

本发明的目的在于：针对现有技术存在的缺陷，提出一种配体小分子关键子结构筛选方法及装置，去除掉大部分的不相关特征，即子结构，一方面可以保留下与配体活性相关的关键子结构，另一方面在后续的学习过程中仅需要在小部分特征上构建模型，大幅提高模型的训练效率。

为了达到以上目的，本发明提供了一种配体小分子关键子结构筛选方法，采用基于配体的药物虚拟筛选方法，对海量的配体分子使用ecfp分子指纹对配体分子进行表征，每一位二进制编码对应了该配体分子是否具有特定子结构；采用基于加强的对偶多面体投影准则的lasso的特征筛选方法，保留下与配体活性相关的关键子结构，用于后续的学习过程构建模型，

其主要包括如下步骤：

步骤1：构建所需要的初始数据集d^s，所述的初始数据集d^s包括生成ecfp特征所需要的配体分子的smiles表示，以及配体分子与药物靶标之间的产生活性作用的浓度值standardvalues，并求其十进对数的相反数作为生物活性response值；

步骤2：给定初始数据集d^s，对初始数据集d^s进行处理，得到配体分子的ecfp特征即数据集d^t；

步骤3：基于数据集d^t，基于edpplasso方法的配体分子特征筛选；

步骤4：基于鲁棒性选择(stabilityselection)方法的配体分子特征选择；

步骤5：特征回溯，关键子结构可视化。

本发明的进一步限定技术方案为：步骤1中所述的生物活性response＝-log10v,其中v是standardvalues值，response值反映配体分子与gpcr作用的生物活性大小，值越小表示作用的活性越低。

进一步的，所述给定初始数据集其中是每个分子的smiles原子连接图，yi是每个样本的response值，所述配体分子的ecfp特征即数据集d^t＝{(xi,yi)|xi∈r^1*m,1≤i≤n}。

进一步的，所述的步骤3中，所述数据集d^t，对于满足条件(λ∈(0,λ0])的λ＝{λi|0≤i<k,λi>λi+1}，得到属于每个λ值的特征筛选结果τ＝{τi|τi∈r^1*m,0≤i<k}，其中τi值为1代表该特征保留，为0代表是无关特征，可删除。

且步骤3中，基于edpplasso方法的配体分子特征筛选，假定数据x∈r^n*m，n为样本个数，m为特征维数，则标准lasso问题为：

其中y为标记，β∈r^m，(1)式中通过对损失函数加上l1范数的惩罚函数，将变量β的系数进行压缩并使某些回归系数变为0，在优化损失函数的过程中同时达到特征选择的目的；

对于(1)式可转换求其对偶问题：

其中θ为对偶变量；

式(1)和式(2)的最优解是相同的，所以只要求其对偶问题即可；

为了方便起见，把优化问题(2)的解记为θ^*(λ)(同理优化问题(1)的解记为β^*(λ))；

根据kkt条件有：

y＝xβ^*(λ)+λθ^*(λ)(3)

其中[*]i代表第i个特征，由式(4)所示的kkt条件，得到

就是一个无关的特征(5)

先估计一个区域θ，其中包含θ^*(λ),(5)可以被写成如下形式：

就是一个无关的特征(6)

进一步的，所述的步骤4中，对于步骤3的特征筛选结果，把k个ti叠加累和，这就得到每一维特征被选择的频数，特征被选择的次数越多，代表其越有可能是相关的特征，选取被选择次数最多的前p个特征，得到鲁棒的特征选择结果,所述的p为自然数0到p。

进一步的，所述的步骤5中，回溯在步骤4中得到的特征选择结果，保留下最终特征,可视化出其具体的子结构图,对于新药研发具有重要意义。

一种配体小分子关键子结构筛选装置，包括如下模块：

数据预处理模块，用于整理数据库中下载的相关文件，具体包括相关文件的下载，数据清洗，对重复的配体分子求取response平均值，构建初始数据集d^s；

ecfp特征生成模块，依据所述数据预处理模块构建的初始数据集d^s构数据集d^t，用于生成ecfp分子指纹，一般情况下，指纹深度越深，生成的特征维数越多，子结构的分析粒度越精细。

基于edpplasso的特征筛选模块，依据数据集d^t，用以通过设置λ值改变筛选的精细程度，去除无关特征；一般来说，λ值越小，保留的特征越多。

鲁棒性选择特征模块，进一步优选关键特征；用以设置关键特征保留数量，筛选出最重要n个子结构，所述的n为自然数0到n。

特征回溯，关键子结构可视化模块；该模块能够可视化出影响药物靶标与配体分子结合活性的关键子结构，指导先导化合物的优化。

所述各模块之间通过数据线传送连接。

本发明的有益效果：

1、解决“维数灾难”问题，快速并大幅去除不相关特征，得到鲁棒的相关特征，使得后续学习过程仅需在小部分特征上构建模型，大幅提高模型的学习效率。

2、特征回溯，可视化关键子结构，可以指导先导化合物的优化工作，使其克服缺陷并最终得到优良的候选药物。这一结果将极大加速新药发现过程。

附图说明

下面结合附图对本发明作进一步的说明。

图1为本发明关键子结构筛选方法流程示意图。

图2为本发明关键子结构筛选装置结构示意图。

具体实施方式

以下结合图1对本发明的方法流程作进一步的详细说明。

步骤1、先从uniprot数据库中下载得到7tmrlist文件，该文件含有所有共3092个g蛋白偶联受体(gpcr)的uniprotid号，我们筛选出人类的gpcr蛋白质，共825个。紧接着从glass数据库中下载allinteractiondata文件，该文件包含了uniquegpcr-ligandentries519,051个。我们按照作用配体的数量对825个人类gpcr蛋白质进行排序，并选取了配体数量最多的前25个gpcr蛋白质。这25个gpcr蛋白质中无任何三维结构的共8个，这8个gpcr的配体个数均大于3000，侧面说明其对于人体功能的重要性，我们选取这8个gpcr作为本文的实验对象。

在本文中，为了更为精确的研究，对这8个gpcr蛋白质，我们分别从chembl数据库中得到它们所有的作用配体并保存这些作用配体分子(smiles格式)及response值。

具体的：表1提供了一个整理后的数据集概况。

表1实验数据集

以表1为例，初始数据集中包括：

canonicalsmiles：配体分子的smiles原子连接图。

response：各配体分子的生物活性值。

步骤2、配体分子ecfp特征生成：

给定初始数据集其中是每个分子的smiles原子连接图，yi是每个样本的response值。对初始数据集进行进一步处理，得到描述样本的ecfp特征，即数据集d^t＝{(xi,yi)|xi∈r^1*m,1≤i≤n}。

把从数据库中得到的每个分子的smiles原子连接图及所需指纹半径输入ecfp生成软件，会得到每个配体分子固定长度的ecfp特征。由于数据集的特征是由全体分子生成的，分子之间会有共同的特征，也会有自身独特的特征，因此把所有分子的特征组合起来，删除重复的共有特征，把剩下来的特征作为最后的特征描述。

步骤3、基于edpplasso方法的配体分子特征筛选：对数据集d^t，应用edpp特征筛选准则，对于满足条件(λ∈(0,λ0])的λ＝{λi|0≤i<k,λi>λi+1}，得到每个λ值的特征筛选结果τ＝{τi|τi∈r^1*m,0≤i<k}，其中τi值为1代表特征保留，为0代表是无关特征，可删除。

其中，edpp特征筛选具体过程如下：

假定数据x∈r^n*m，n为样本个数，m为特征维数，则标准lasso问题为：

其中y为标记，β∈r^m，(1)式中通过对损失函数加上l1范数的惩罚函数，将变量β的系数进行压缩并使某些回归系数变为0，在优化损失函数的过程中同时达到特征选择的目的。

对于(1)式可转换求其对偶问题：

其中θ为对偶变量。式(1)和式(2)的最优解是相同的，所以只要求其对偶问题即可。为了方便起见，把优化问题(2)的解记为θ^*(λ)(同理优化问题(1)的解记为β^*(λ))。根据kkt条件有：

y＝xβ^*(λ)+λθ^*(λ)(3)

其中[*]i代表第i个特征，由式(4)所示的kkt条件，得到

就是一个无关的特征(5)

换句话说，可以充分利用式(5)来为lasso问题找出无关的特征。

但是，因为θ^*(λ)是一个不知道的值，不能直接利用式(5)来找出无关的特征。所以可以先估计一个区域θ，其中包含θ^*(λ)。所以，式(5)可以被写成如下形式：

就是一个无关的特征(6)

综上，只要能找到一个区域包含θ^*(λ)，使得区域内每一个θ与乘积的绝对值小于1即可，这样式(6)就可以作为一条准则用来为lasso问题找出无关的特征。可以看出，找到的区域θ越小，θ^*(λ)的估计就越正确，那么有更多的无关特征就可以通过规则筛选出来。

步骤4、基于鲁棒性选择方法的配体分子特征选择：对于步骤3的特征筛选结果，把k个ti叠加累和，这就得到每一维特征被选择的频数，特征被选择的次数越多，代表其越有可能是关键特征，选取p个被选择次数最多的特征，得到鲁棒的特征选择结果。

步骤5、将步骤4得到的鲁棒特征选择结果进行特征回溯，可视化即为关键子结构，指导后续的新药发现实验。

而上述步骤2和步骤3的特征筛选过程是通过两步优化，首先根据选取对于满足条件的一定数量的λ＝{λi|0≤i<k,)i>)i+1}，本实施例中一定数量选取100。对应于相同数量的特征筛选结果τ＝{τi|τi∈r^1*m,0≤i<k}，其中τi值为1代表特征保留，为0代表是无关特征，可删除。然后综合考虑这些特征筛选结果，把k个ti叠加起来，这就得到每一位特征被选择的频数，特征被选择的次数越多，代表其越有可能是相关的特征，选取p个被选择次数最多的特征，得到鲁棒的特征选择结果，所述的p为自然数0到p,从而避免了使用单一参数值可能带来模型效果差的问题。

本实施例还公开了一种新的配体小分子关键子结构筛选装置，如图2所示，包括：

数据预处理模块，用于整理数据库中下载的相关文件；包括相关文件的下载，数据清洗，对重复的配体分子求取response平均值。

ecfp特征生成模块，用于生成ecfp分子指纹；ecfp特征生成模块可以设置指纹深度，一般情况下，指纹深度越深，生成的特征维数越多，子结构的分析粒度越精细。

基于edpp的特征筛选模块，去除无关特征；该模块可以通过设置λ值改变筛选的精细程度，一般来说，λ值越小，保留的特征越多。

鲁棒性选择特征模块，进一步优选关键特征；可以设置关键特征保留数量，并筛选出最重要n个子结构,所述的n为自然数0到n。

特征回溯，关键子结构可视化模块；该模块能够可视化出影响药物靶标与配体分子结合活性的关键子结构，对于指导先导化合物的优化有极其重要的作用。

同时，各模块之间通过数据线传送连接。

本发明的有益效果概括如下：

(1)本发明解决了配体分子海量特征的“维数灾难”问题，一方面，生成配体分子海量特征能够解决某些分子指纹特征维数局限的问题，另一方面，使得后续学习过程仅需在小部分特征上构建模型，大幅提高模型的学习效率。

在这个数据爆炸的年代，高维度的数据随处可见。例如在许多生物相关的问题中，数据的维度非常高，由于收集数据需要昂贵的实验，因此可用的训练数据相当少。这个时候就会出现特征维度远大于样本数的情况，如果不做其他假设或限制的话，模型难以构建，同时也会引起过拟合的问题。而lasso方法通过构建一个惩罚函数剔除不相关的特征，解决了“维数灾难”问题，使得后续学习过程仅需在小部分特征上构建模型，大幅提高模型的学习效率。

(2)本发明利用基于edpp准则的lasso方法，快速并大幅去除不相关特征，得到鲁棒的相关特征，利于理解与配体活性相关的子结构，通过特征回溯模块给出了具体的关键子结构，对后续的新药发现工作有重要意义。

前述对本发明的具体示例性实施方案的描述是为了说明和例证的目的。这些描述并非想将本发明限定为所公开的精确形式，并且很显然，根据上述教导，可以进行很多改变和变化。对示例性实施例进行选择和描述的目的在于解释本发明的特定原理及其实际应用，从而使得本领域的技术人员能够实现并利用本发明的各种不同的示例性实施方案以及各种不同的选择和改变。本发明的范围意在由权利要求书及其等同形式所限定。