廖灵志 杨智欢 洪玲慧 李 瑾
(西南医科大学 医学信息与工程学院,四川 泸州646000)
分子对接[1]算法中的采样方法,即构象搜索方法对于分子对接的精度、尤其是效率有着至关重要的影响。早期构象搜索中应用的算法包括:遗传算法[2]、模拟退火[3]、禁忌搜索[4]等各种方法。随着群体智能算法的发展,一种新的群体智能算法——烟花算法[5,6]因其具有很强的优化问题求解能力而倍受关注。本文对分子对接中的构象搜索问题展开研究,提出一种基于烟花算法的构象搜索方法,我们的工作证明了烟花算法在解决蛋白质- 配体对接的构象搜索问题方面表现出的强大性能,具有更快的收敛速度,更加稳定的表现。
2.2.1 烟花初始化
2.2.2 爆炸算子
将烟花按适应度值进行排序,通过烟花的适应度值排序进行火花数目的分配,而不再依赖适应度值数值本身。采用公式(3)的方式进行火花数目的计算,公式如下所示:
其中,M为每个烟花最多爆炸的火花数目;ri表示烟花wi的适应度值排序。
我们将该构象搜索算法的解空间搜索范围分为三个层次。一是构象中心位置在蛋白质结合口袋中的变化;二是构象中心位置固定后,取向的变化;三是中心位置和取向固定后,可扭转键的键角变化。因此,对应解空间的这三个层次,将烟花爆炸的范围划分为如下三个区域:
Ω1:变化区域为可扭转键变化,范围为[-π, π];
Ω2:变化区域为取向变动,范围为[-π, π];
Ω3:变化区域为中心位置变动。
从Ω1 到Ω3 代表爆炸半径越来越大。
若爆炸后的火花为不可行解,则在解空间里随机生成一个火花。
2.2.3 变异算子
FWAVina 采用随机变异的方式来产生变异烟花,以增加种群多样性,避免陷入局部最优。从烟花种群(包含N 个初始烟花及S 个爆炸出的火花)中随机选择N 个烟花,每个烟花随机选择几个维度,将这些维度上随机产生变量得到变异个体。
2.2.4 烟花选择策略
当代种群中的N 个初始烟花、S 个火花、N 个变异烟花构成了候选集合K,从候选集合K 中选择N 个个体作为下一代初始烟花。
FWAVina 构象搜索方法的具体步骤可以分成如下几步。
2.3.1 烟花初始化。初始化N 个烟花(每个烟花代表一个配体的构象),将配体构象表示为解向量,设置AutodockVina 打分函数为适应度函数。
2.3.2 烟花爆炸。根据烟花爆炸范围的计算公式和产生火花数量的计算方法,将初始烟花进行爆炸,产生爆炸火花,其数量用S 表示。
2.3.3 烟花变异。从N 个烟花和S 个火花中随机选择N 个个体,每个个体按照烟花变异策略发生变异,形成变异烟花。
2.3.4 烟花选择。将候选集合(包括当代种群中的N 个烟花、S 个爆炸火花和N 个变异火花) 中的所有烟花按适应度函数值从小到大排序,按选择策略选择N 个烟花组成下一代种群。
2.3.5 判断适应度函数值是否满足收敛准则,如果是,转步骤2.3.7,否则,转步骤2.3.6。
2.3.6 判断是否达到最大迭代次数,如果是,转步骤2.3.7,否则,转步骤2.3.2。
2.3.7 输出近似最优配体构象以及分值。
为了测试FWAVina 在对接和虚拟筛选中的性能,我们使用了广泛应用的标准数据PDBbind 的核心集。对核心集中的每个复合物分别进行了30 次对接,然后每个复合物的预测结合能、RMSD 值及运行时间均取其平均值。201 个复合物的平均结合能、平均RMSD 值、平均运行时间见表1。我们分别对Vina 和FWAVina 预测的201 个复合物的结合能、RMSD 值以及运行时间,进行配对样本t 检验。p 值表明FWAVina 预测构象的结合能、RMSD 值以及运行时间均显著低于Vina。另外,表1 中展示了Vina 和FWAVina 30 次对接的平均准确率,即201 个复合物对接一次后计算准确率,一共对接30 次,准确率取其平均值。
在运行效率方面,与Vina 相比,FWAVina 的执行时间减少了52.5%,见表1。显然,FWAVina 的对接速度有了很大的提高。在基于分子对接的虚拟筛选应用方面,分子对接程序的速度对虚拟筛选的效率起决定性作用,因此一个快速的分子对接程序更适合应用于药物的虚拟筛选。
表1 在Coreset 数据集上进行30 次对接后Vina 和FWAVina 的对接性能对比
本文提出了一种基于烟花算法的分子构象搜索方法,并且在AutodockVina 的框架上予以实现,编写了分子对接程序FWAVina。其次,本文在PDBbind 数据集上进行了分子对接模拟。结果表明,FWAVina 与Vina 相比,分子对接的准确性略有提升,而对接效率提升较大,并且FWAVina 对于不同柔性的配体对接来说迭代次数不会呈现大幅上升。