流感病毒神经氨酸酶抑制剂的研究

本节将要介绍两个实例，分别是针对流感病毒神经氨酸酶NA的已知药物的改造以及新药的筛选与设计，系统地阐述CADD在新药开发中的应用。

10.3.1 神经氨酸酶简介

引起禽流感（AI）的病原为禽流感病毒（avian influenza virus），是负链RNA病毒，该病毒属正粘病毒科流感病毒属，对呼吸道系统都有致病性等，特别是这两种病毒对黏多糖和糖蛋白具有特殊的亲和力，尤其是对细胞表面的含唾液酸的受体具有更强的亲和力。神经氨酸酶NA是AIV的一个重要的表面抗原，由RNA片段6编码。NA由4个亚单位构成，即2个相同的糖蛋白的聚体以二硫键相连成聚合物，其具有唾液酸酶活性，并在病毒出膜时起到释放位于被感染区细胞表面的子代病毒作用。该功能不仅可以避免释放过程中病毒的聚集，而且还有可能有助于病毒穿透上皮组织的黏液层以完成对细胞的吸附。神经氨酸酶的这个特性也使其成为抗病毒药物的靶标，病毒的HA和NA之间有拮抗效应，两者相互协调地发挥作用是病毒有效吸附和释放的关键。

HA和NA均是易发生变异的2种免疫抗原。目前，根据它们的变异，将HA分为H1～H15共15个亚型，NA分为N1～N9共9个不同的亚型，当前，人类流感疾病在历史上主要由3种HA亚型（H1、H2、H3）和2种NA亚型（N1、N2）引起。但是最近人类感染的流感病毒又增添了新的HA亚型，包括H5、H7和H9。由H5N1禽流感病毒人感染病例，致死率超过50%，针对H5N1病毒的研究也最为广泛。

10.3.2 基于神经氨酸酶N1晶体结构对奥司米韦的修饰

10.3.2.1 奥司米韦的发展史

奥司米韦（Oseltamivir）是一种作用于神经氨酸酶的特异性抑制剂，通过抑制神经氨酸酶的作用，可以抑制成熟的流感病毒脱离宿主细胞，从而抑制流感病毒在人体内的传播以起到治疗流行性感冒的目的。Oseltamivir是基于结构的合理药物设计的成功案例，在这种药物的研发过程中大量应用了CADD的手段，根据受体靶酶的三维结构有针对性地设计了高效低毒专一性强的神经氨酸酶抑制剂。

奥司米韦是利用天然底物唾液酸类似物与流感病毒神经氨酸酶活性位点相结合的结晶X衍射结构，进行合理化药物设计而发现的。奥司米韦是通过对唾液酸类似物骨架进行修饰（包括加了亲脂侧链）发展而来的，这种修饰主要是为了使药物可用于口服。随着抗Oseltamivir病毒株的出现，而Zanamivir依然对这些抗药性病毒株保持着敏感性，所以对Oseltamivir的研究并没有结束。

图10-26为唾液酸、DANA、Oseltamivir、Zanamivr以及对Oseltamivir进行修饰用到的骨架分子的结构，从中可以看出这些分子的联系与不同：Oseltamivir与Zanamivir的最大不同就是4位和6位的基团，两个分子6位基团一个是疏水性强的戊烷基，一个是亲水性强的甘油基，而对Oseltamivir的修饰策略也正是希望找到一个疏水和亲水性比较适中的基团去适应N1口袋疏水环境的变化。

图10-26 唾液酸、DANA、Oseltamivir、Zanamivr以及对Oseltamivir进行修饰用到的分子骨架结构

10.3.2.2 H5N1-NA晶体结构分析

神经氨酸酶NA作为抗禽流感病毒的靶点一直是药物设计关注的焦点， H5N1-NA神经氨酸酶晶体：2HTY（PDB bank code）的解析成功成为抗禽流感道路上的又一重大事件，其为合理的新药物设计提供了N1靶点的三维结构信息以及成功的机会，因为现有的NA药物都是基于N9、N2和B型流感病毒的NA设计的。

根据免疫原神经氨酸酶NA表达不同，NA分为N1～N9共9个亚型，分成两组，组1包括N1、N4、N5和N8；组2包括N2、N3、N6、N7和N9。Russell等人的研究表明组1中的神经氨酸酶N1（2HTY）结构与组2中的N9（1F8B）有很大的不同。主要的区别在于活性口袋位置由氨基酸残基147～152构成的所谓“150-Loop”区域。为了突出这个位置特点，将N1：2HTY和N9：1F8B这两个不同组的同源蛋白的序列进行比对（见图10-27），三维结构进行了骨架叠合（见图10-28）。

图10-27 神经氨酸酶N9（Chain-1）与N1（Chain-2）序列比对（彩图见第410页）

图10-28 （a）神经氨酸酶N9（2HTY）与N1（1F8B）三维骨架叠合图；（b）奥司米韦（Os-eltamivir）在神经氨酸酶N9中靠近150-loop附近的活性口袋中的微观示意图；（c）奥司米韦在神经氨酸酶N1中的微观示意图（彩图见第410页）

从图10-27可看出，N1晶体的150-Loop由氨基酸残基片段“GTVKDR”组成， N9相应位置由“GTIHDR”组成。在图10-28中看到由于这一段氨基酸残基片段的变化，N1的晶体构象出现了没有预见到的情形，在临近保守性活性口袋的边缘这一段变化的片段向外打开，N1中的150-Loop中的氨基酸Val-149的Cα与N9中相应位置的Ile-149的Cα有7Å的距离，Asp151的羧酸基也指向口袋外，因此N1的口袋相对于N9的就成为一个敞开式构象（open-up conformation）[见图10-28（b）、图10-28（c）]。除此之外，N1中的保守性的功能性氨基酸残基Glu-119的侧链也与N9的取向不同。由于150-Loop的存在，N1中的酸性氨基酸残基Val-149和Asp151远离了Oseltamivir上4位的氨基（—NH2），在配体与受体其他活性位点相互作用相同的情况下，口袋150-Loop位置氨基酸残基构象的变化使得分子间的结合能力降低了。

Oseltamivir与N1晶体2HTY结合时一般可分两步，第一步是抑制剂与敞开式构象（open-up conformation）结合，然后口袋构象慢慢变化成为类似N9的关闭式构象（closed conformation）导致抑制剂与受体的结合力得到加强。总体来讲，在2HTY晶体结构中150-Loop的存在使得这种敞开式（open-upconformation）构象的能量低于关闭式（closed conformation）构象，所以一个低能量的敞开式构象更容易和抑制剂结合，其具有能量上的优势。研究同时也表明，Oseltamivir和Zanami-vir都有和2HTY敞开式构象结合的能力。所以对Oseltamivir的修饰希望能得到一个新的抑制剂让它更好地、选择性地与N1结合，尤其是在局部与150-Loop结构得到很好的匹配。

10.3.2.3 Oselatmivir的修饰

基于以上对2HTY晶体结构的分析，对Oseltamivir的修饰列出了其中53个典型的分子结构（见表10-12）。为了提高Oseltamivir C4—NH2的结合能力， C4—NH2分别变为NHC（ NH═ +2 ）和NHC（ N—CH═ +4 ）NH2；C8—戊烷基团改成部分疏水、部分亲水的基团，以便于这个基团能适应疏水环境的变化，修饰后的基团又与Zanamivir对应部分不同。考虑到Oseltamvir的羧酸基团与受体的三个精氨酸很好的结合能力，以及—NHCOCH3与由Trp178和Arg224形成的疏水位置的匹配，所以对Oseltamivir的这两个位置的基团未作改动。另外，由于C6—O—C8的柔性，C8的手性也未考虑改动，根据这个分子修饰的原则，表10-12中列出了修饰后的54个Oseltamivir的类似物分子，包括Oseltamivr本身在内。

表10-12 Ose Itamivir的相似化合物

（续表）

（续表）

（续表）

10.3.2.4 分子对接结果及讨论

为了研究分子中R和X取代基对分子间相互作用力的贡献问题，所有的分子与受体2HTY采用分子对接软件Dock6.1版本进行对接计算。新版本的Dock6.1对接计算对构象的评估打分采用两步进行。

第一步采用分子格点打分（Grid Score），即结合能由范德华作用能与静电作用能加和确定。即Ebinding=Evd W+Ees。

第二步采用Amber Score来给第一步Grid Score筛选的最好构象再次进行打分。在AMBERScore中加入了GB/SA（generalizedborn/surfacearea）模型考虑了溶剂化效应，而且在这一步，受体不再是刚性的，活性口袋的氨基酸和配体都是柔性的，允许配体周围的氨基酸残基柔性，在对接计算过程中就产生了所谓“诱导契合”（inducedfit）的效果。

Amber Score由Ebinding=Ecomplex-（Ereceptor+Eligand）定义，其中Ecomplex代表复合物的能量，Ereceptor、Eligand分别代表受体和配体的能量。配体与受体的结合自由能E由下式得到：

E=EMM+（Ep-sol+Enp-sol）

EMM=Evd W+Ees+Eint——依据AMBER力场得到的势能函数

Ep-sol——溶剂化自由能的极性部分

Enp-sol——溶剂化自由能的非极性部分

其中EMM为结合自由能的势能部分。计算中可以增加或减小能量最小化（minimi-zation）和分子动力学（MD）模拟步长（steps）的大小，由于计算时间的问题并不是步阶越大越好。对不同的蛋白测试表明，能量最小化和分子动力学的步阶分别取100和3000计算得到的结果最好。

对于每一个分子的对接结果，在前100个构象里挑选出几何构型最佳匹配的那一个，然后将这54个分子按照取代基的不同分成a、b、c三个系列，每个系列的18个分子构象叠合在一起，如表10-13所示。a系列是0a～17a，X取代基是—NH+3的叠合；b系列是X取代基，是—NHC（ NH+═ 2 ）的叠合；c系列是X取代基，是—NHC（ N—CH+═ 4）NH2的叠合。修饰后的R基团，由于部分疏水和部分亲水的特点，与受体2HTY活性口袋内的疏水环境得到很好的匹配。修饰后的配体能否具有与受体更好的结合能力，从表10-13的具体数据便能说明问题，在结果中可以发现以下几个特征：

（1）X基团相同时，每个系列的配体从0a～17a的分子格点打分有数值减小的趋势，即结合能力趋强，这就说明分子的修饰达到了预期的效果。

（2）R基团相同时，分子格点打分有a系列＜b系列＜c系列的规律，c系列化合物X基团上的甲基（—CH3）由表可以看出，甲基深入到一个疏水孔洞里，所以延长C4-NH2的空间可能会加强分子间的亲和力。但是也发现，无论是-NHC （ NH+═ 2）还是-NHC（ N-CH═ +4）NH2，都无法越过由Asp151氨基酸残基形成的“障碍”（见图10-28）。因此，更好的办法是设计一个新的分子片段匹配150-Loop口袋，然后将这个片段与C4-NH2键和成为一个分子，这样极有可能会产生一个神经氨酸酶抑制剂新系列。

表10-13 化合物0c～17c、0b～17b和0c～17c与神经氨酸酶N9（2HTY）之间的结合自由能

（3）对于一个分子格点得分很低的配体来讲，它的Amber力场打分却不一定也低，如c系列分子的Amber力场打分就多数没有b系列的低，除了1c、2c和8c三个例外，这三个分子有一个共同点就是它们都比其他的分子少一个C原子，所以造成这个结果的原因可能是分子的体积问题。体积太大时可能去溶剂化比较困难，当体积合适时效果才好。因此，从两种打分情况看，b系列的化合物去溶剂化的效果要好。

按照分子格点打分+Amber力场打分值最大原则，以Oseltamivir得分为参考，在每个系列化合物中挑选出几个最有希望成为药物候选物的分子，在表10-13中以粗体显示，也许这种选择有不妥之处，因为任何一种药物候选物分子都要用随后的生物活性实验、生物利用度、药物溶解度、毒性测试等一系列实验来验证理论设计的分子是否能成功。

10.3.2.5 结论

通过对Oseltamivr分子的局部修饰，从结果上看出，修饰得到了预期的效果，绝大多数修饰后的分子与受体对接的结果都要好于Oseltamivir，所以基于受体的结构信息，有针对性地改造配体的部分片段，使得分子的药效团与受体口袋相应位点达到很好的匹配，这种直接药物设计的方法现在看来是非常有效也很流行的。

10.3.3 3D-QSAR在H5N1抑制剂模型中的应用

本节将根据上市效果最好的两种神经氨酸酶抑制剂扎那米韦（Zanamivir）和奥司米韦（Oseltamivir），来进行神经氨酸酶抑制剂的改造、筛选和设计，以获取高效、毒副作用小的新型流感病毒抑制剂。

定量构效关系（quantitative structure-activity relationship，QSAR）分析是搜寻新的潜在药物分子的一种通用的手段。在设计和改进神经氨酸酶抑制剂方面起着重要的作用。本小节分别运用偏最小二乘法（PLS），比较分子力场分析（Co MFA）法和分子相似性指数分析法（Co MSIA）对环己烯类衍生物类似神经氨酸酶抑制剂的活性进行了分组和整体两种三维定量构效关系（3D-QSAR）研究，分别建立了各自合理的预测模型，两个模型均可以有效地指导先导化合物结构改造和新的神经氨酸酶抑制剂设计，并通过3D-QSAR给出的等势面图加以分析，方便、直观地指导新型药物的改造工作。

10.3.3.1 数据来源

通过查找文献从中挑选的126个类似神经氨酸酶抑制剂化合物（neuramini-daseinhibitors，NIs），并经过Auto Dock软件对接后获得结合自由能最低的构象作为研究对象。这些类似NIs化合物的结构差异较大，为提高模型的预测性，在具体QSAR模型建构及数据分析时按共同骨架的异同将126种化合物分成了6类，各类中按比例挑选一部分分子作为本组训练集来构建COMFA和COMSIA模型，本组其余部分小分子做预测集用以对本组所建模型进行验证。它们的抑制活性指标用半数抑制浓度IC50的负对数-lg IC50即p IC50进行量度，数据依然是从文献中获得。

10.3.3.2 参数的选择及模型的拟合方法

在QSAR模型建构中，参数选择过多可能会造成共线性或过拟合的不良影响，模型的预测能力也会减弱。选择数量适当的参数，有助于建立预测能力强、稳定的模型。一般来说，若样本容量为N，变量数为M，变量的选择应遵循N/M≥5经验规则。本节实验样本数一共为126个，分为六组，各组变量依据经验原则数量不等。QSAR模型评价的统计学指标有：交叉验证系数Q2，相关系数R2，Fisher检验值（统计方差）F和均方根偏差S。R2和F值越大，S值越小，模型的拟合能力就越强。其中R2和S由系统自动产生，F值由式（10-1）计算。

式中：n为样本个数；r2为相关系数；k为描述值个数。

10.3.3.3 模型的预测能力评价

10.3.3.4 Co MFA和Co MSIA的原理

比较分子力场分析（Co MFA）是3D-QSAR常用的研究方法。它通过计算，定量地将分子的生物活性与其周围的力场（如静电力、范德华作用力等）联系起来，应用数学方法建立化合物的生物活性与周围各力场之间的空间分布关系的模型，为分子结构的修饰改造提供信息。其基本假设是：①在分子水平上，影响生物活性的因素主要是非共价性相互作用；②分子力场则侧重于处理立体场和静电场的非共价相互作用，可以很好地解释分子的各种理化参数性质。因此计算分子中原子的电荷就显得极为重要，不可或缺了。Co MFA方法是通过将一系列具有相同性质的分子按照其相同的几何作用点在三维空间进行叠加，计算这一组分子叠加的静电场和立体场，用某种探针原子对这些场进行作用，测试空间各点作用的能量，最后用偏最小二乘法（PLS）及交叉验证法（cross validated）获取预测模型，其中PLS是通用的主流方法。

比较分子相似性指数分析（Co MSIA）方法是在Co MFA方法基础上延伸发展起来的一种新颖的构建三维定量构效关系的方法。为了有效避免分子表面附近网格点上能量的显著变化，分子场采用了与距离相关的高斯函数形式，与Co MFA法相比，不需要定义能量阈值。Co MSIA方法考虑分子静电场、立体场、疏水场、氢键受体场、氢键给体场的贡献，用高斯函数表示，是对Co MFA的改进。

10.3.3.5 参数的设定和模型构建

本实例首先使用SYBYL中自带的优化功能对对接后所得的126个小分子结构进行几何上的优化，具体参数设置如下：使用Amber4.0立场中的Amber-all-at-om来设定电荷和将Powell中的根均方差能量梯度值设定为0.05kcal/（mol· Å），整个系统进行最小化收敛。在计算能量的过程中，使用绝缘系数5来模拟在蛋白酶环境中溶剂对抑制剂的溶剂化影响。并按照分组，依据共同骨架把这些构象和位置叠合在一起，即将药效团或力场等进行分子叠合，使每一个分子采用相同最大力场相似的空间取向进行计算。叠合时以本组最高活性化合物的最优构象为模板分子，采用Align Database方法，将其余化合物与模板分子叠合。为进行更深入的Co MFA和Co MSIA研究分析做准备。

本实例考虑到不同空间网格大小和不同的探针对3D-QSAR模型的预测能力都有一定的影响，所以对网格的大小及探针的类型对模型的影响都做了模拟测试，通过调试获取最佳匹配。模拟结果表明：空间的网格从0.3nm减小到0.2nm，对Co MFA的预测结果有所改善，其交叉验证的值有所提高。通过使用了各种不同的探针，包括带SP3杂化的碳原子、氮原子、氧原子及卤素离子和氧负离子，还有钠、钾离子等，发现不同的探针对预测结果均有一定的影响。其中当探针原子选择带1个正电荷的SP3杂化碳原子时交叉验证值最大，模拟效果最好，成为后续Co MFA分析中的空间场探针。模型构建的具体过程为：

（1）Co MFA：用SP3杂化的碳原子作为空间场探针，C+离子作为静电场探针来测定空间场和静电场的相互作用。在x，y，z轴方向上格子的长度设定为2.0Å， Co MFA将会自动生成打格区域，这些打格区域将包围所有的分子，并在每个方向上有4.0Å的延伸。静电场和空间场的最大能量阈值采用空白设置，为30kcal/mol。采用标准选项进行变量变换，留一法、SAMPLE方法作为交叉验证来确定OPC的数目。最后COMFA模型使用这个优化的OPC数在非交叉验证的情况下建立。

（2）Co MSIA：在Co MSIA模型构建过程中，五种物理化学特性，即静电场、立体场、疏水场、氢键的受体和供体，将被逐个计算。静电场特性可以表现为Gasteiger Hückel电荷，立体场的贡献可以通过原子半径的第三种力来描述。评估原子疏水性时按照Ghose及其合作者建立的参数进行操作。格子的空间要像在Co MFA中一样，大到足够包围整个配体分子。这里为避免Co MSIA中的单独格点位置依赖于化合物的物理化学特性，没有随意设置阈值。用半径为1Å的正一价疏水离子探针，衰减因子0.3来计算与前面相似的指数。统计分析采用方法与Co MFA相同。

10.3.3.6 结果与讨论

表10-14为126个小分子的活性实验值，还有各自通过本组训练集叠合，构建模型的模型预测值，为方便比较，突出分组构建模型的优越性，将它们与不分组的模型预测值集中于一个表格中。

表10-14 抑制剂的活性值及分组与不分组时Co MFA和Co MSIA预测值

（续表）

（续表）

（续表）

（续表）

注：有*的分子为测试集，无*的分子为训练集。

1）分组模型与整体模型的比较

（1）第一组分子：包括1～57号分子，以活性最高的是35号分子为模板，其中带有星号的分子为本组的预测集，不带星号的为训练集，其共同骨架为A（见图10-29）。

图10-29 第一组分子，包括1～57号分子

将其抑制剂的生物活性（p IC50）以及Co MFA方法和Co MSIA方法的预测值列于表10-14，通过表10-14，给出了各个组化合物的模型预测活性和实验活性的对比以及与整体模型的预测值。将表10-15中的第一组化合物统计表和已知的相关性一起考虑，可以推断Co MFA和Co MSIA模型对于骨架相同的类似化合物有较高的预测能力。为了进一步验证构建模型的可靠性，用这两种模型来预测测试集中5种抑制剂的抑制活性。结果也在表10-15中列出。对于测试集的Co MFA和Co MSIA模型的相关系数r2的预测值分别是0.895和0.920。这表明两种预测模型是可靠的，对于含有相同骨架的化合物的预测性较高，可以用来指导改造同类新型抑制剂。

表10-15 第一组分子的Co MFA和Co MSIA模型统计分析

（2）第二组分子：包括58～69号分子，以活性最高的60号分子为模板，以B为共同骨架（见图10-30）。

图10-30 第二组分子，包括58～69号分子

通过叠合构建第二组的Co MFA和Co MSIA，其对应的活性值和预测值见表10-14所示。通过观察预测值与实验值的残差、统计表（见表10-16）以及相关活性可以看到第二组的模型预测性也非常可靠，比起不分组模型的预测有很大的改善。为进一步验证模型的可靠性，也对本组中测试集的小分子进行预测，其相关系数分别为0.971和0.843，本组模型通过验证再次说明了其对骨架相同的抑制剂小分子的预测性是可信的。

表10-16 第二组分子的Co MFA和Co MSIA模型统计分析表

（3）第三组分子：包括70～74号分子，以70号为模板，72号做预测，共同骨架为C（见图10-31）。

图10-31 第三组分子，包括70～74号分子

由于实验得到含有相同骨架的小分子的活性值较少，训练集也比较少，故模型的统计表（见表10-17）和相关系数虽然比较高，其通过对72号分子的验证也和实验值非常贴近，但其统计学意义并不大，仅作参考，需要对含有相同骨架的类似分子做预测时，不分组的模型的预测性也比较可靠，具有更好的统计意义。

（4）第四组分子：包括75～94号分子，以90号分子为模板，其骨架为D（见图10-32）。

表10-17 第三组分子的Co MFA，Co MSIA模型统计分析表

图10-32 第四组分子，包括75～94号分子

叠合方法与以上三组相同，通过表10-14实验值与模型预测值比较看出，预测值和实验值非常靠近，并且正负误差分布均匀，绝对值较小，模型预测性较高。模型统计表（见表10-18）以及相关性也显示了相同的结果，再看本组模型对预测集预测的结果，其Co MFA和Co MSIA模型的相关性分别为0.918和0.957，预测集相关性也非常高，最终得出第四组模型完全可靠的结论。

（5）第五组分子：包括95～109号分子，以99号分子为模板，其骨架为E（见图10-33）。

（6）第六组分子：包括110～126号分子，以126号抑制剂为模板，以F为骨架（见图10-34）。

表10-18 第四组抑制剂分子的Co MFA和Co MSIA模型统计表

对第五、第六组采用与以上四组完全相同的叠合操作，分别构建各自的预测模型。其活性及预测值列于表10-14中，其中第五组的统计表为表10-19，第六组统计表为表10-20。参照它们的活性值与预测值以及各自的统计表，还有对各自测试集的预测可以看出，第五、第六组的模型搭建也是非常成功的，完全具备预测骨架相同的抑制剂分子的能力。

图10-33 第五组分子，包括95～109号分子

图10-34 第六组分子，包括110～126号分子

将126个小分子不分组的模型构建时，首先将各组活性最高的小分子的最优构象（36、60、70、90、99、126号）在SYBYL软件包中进行手动叠合，然后再将组内其他分子分别与所在组活性最高的抑制剂分子进行构象叠合，使它们处在同一个三维空间内，来构建Co MFA和Co MSIA模型，它们对小分子的预测值列于表10-14中。统计表如表10-21所示。

表10-19 第五组抑制剂分子的Co MFA和Co MSIA模型统计表

表10-20 第六组抑制剂小分子活性统计表

通过表10-15～表10-21的7个模型统计表可以看出，分组模型构建比整体模型构建的相关系数有所提高，其对骨架相同的分子预测能力更强。对抑制剂的预测能力更为可靠，为了更直观地比较其模型的可靠性，分别把它们的实验值与预测值作线性拟合，观察散点分布图（见图10-35和图10-36，分组模型实验活性值与预测值的散点图；图10-37、图10-38整合模型实验活性值与预测值的散点图），可以看出，图中直线均近似为通过原点、斜率为1的直线，但分组模型构建得到的直线A值更小，几乎完全靠近原点，B值更大，较之整合的Co MFA，Co MSIA模型更近于1，其散点分布也更为集中、密集、均匀地分布在直线上以及附近，较之整体模型离群点更少。所以对于骨架相同的抑制剂进行模型构建及分析预测时，分组构建更为可靠，其预测的活性值和实验以及真实的活性值更为接近，对新型药物筛选和设计的指导有非常重要的意义。但未知抑制剂的分子与以上骨架均不相同时，整体模型就会体现出其适用性广泛的优点。

表10-21 整体抑制剂的Co MFA和Co MSIA模型统计表

图10-35 分组Co MFA预测值与实验活性值散点图（A=0.34843，B=0.94554）

2）Co MFA和Co MSIA模型的比较

最终模型的取舍及评价可通过抽一法交叉验证的相关性系数q2作为依据。一般而言，当q2＞0.6时，模型的可靠性较好；若q2＜0.6，而相关性系数r2较大，则模型面临较高的过拟合风险。传统验证系数r2的值较高是由于采用了不同的描述变量（空间场、静电场、疏水场、氢键受体和供体），这些变量适合描述神经氨酸酶和抑制剂间的相互作用，同时也适用阐述抑制剂周围的力场特性。通过对分组以及整体模型的7个统计表的对比分析可以看出，除Co MFA模型中的吡咯烷衍生物和Co MSIA模型中唾液酸衍生物外，在这两种模型中其他各类化合物的交叉验证系数都大于0.6，因此建立的两个QSAR模型都有一定的预测能力。通过更为深入的比较分析，可以发现在Co MSIA模型中立体场和疏水场贡献的和与Co MFA模型中立体场的贡献大体相同；暗示了这两种模型是相互联系并且内在一致的。而Co MFA中静电场的贡献相当于Co MSIA中静电场和氢键场贡献和的60%，这表明了极性的相互作用是抑制剂与神经氨酸酶反应的必要条件。这与神经氨酸酶的键合位点有10个电荷和只有4个有疏水性特性的氨基酸表现出高度的极性完全吻合。

图10-36 分组Co MSIA预测值与实验活性值散点图（A=0.17762，B=0.97224）

图10-37 整合Co MFA预测值与实验活性值散点图（A=0.64443，B=0.89882）

图10-38 整合Co MSIA预测值与实验活性值散点图（A=0.6263，B=0.90167）

通过整合Co MFA和Co MSIA模型直观结果可以通过St Dev*Coeff图选项获得，可视化的等势图可按照贡献因子得到。支持水平设定为0.8，不支持水平设定为0.2。在图10-39中，给出了Co MFA的等势图。

图10-39 COMFA的势能图[嵌合1号分子（GS4071）]（彩图见第411页）

GS4071与抑制剂在4Å距离内的神经氨酸酶残基都可在辅助的可视化中看到。在这个模型中的静电场项，有两个高极性区域是值得注意的（强极性区域）。在抑制剂的C1取代位置（羧酸的官能团处），有一个高的静电偏爱区。在这里抑制剂的羧酸与神经氨酸酶的2个氨基酸残基Arg292、Arg371有强烈的静电相互作用。在这个模型C4胍基附近有显著的亲电性倾向，这可以满足有3个氨基酸残基Glu119、Asp151、Glu227产生的负电荷环境的要求。显然，静电和氢键的相互作用都与这两个强极性区域相关，同时这两个区域是有效的神经氨酸酶抑制剂与神经氨酸酶键合所必需的，是药物设计时必须考虑的两个最重要的部分，羧酸和氨基或胍基是保持抑制剂生物活性必要的官能团。

除了上面提到的两处高极性区域外，在Co MFA等势图上，在C2和C3处也观察到强极性区域，这意味着在这些位置有正电荷的取代基可以增加NI的活性。例如，70号帕拉米韦的活性就比它的类似物74的要高，这是因为后者C2羟基上缺少带正电的氢。然而，为了维持这些位置正电荷的需求，用亚氨基代替2-CH2，得到一系列的新型吡咯化合物95～109，也具有高效的神经氨酸酶抑制活性。但是在等势图中，靠近C′2乙基正丁基出现了静电作用区域与1号GS4071的显示图像不符。事实上，这个区域与另外一种高效药物达菲的极性甘油的侧链相关，因此，这个模型实际上反映了在这个区域通过静电相互作用而得到的强的键合作用。

至于空间相互作用，在Co MFA模型中显示，在C1-羧基旁边存在一个大的空间阻碍作用区域，这意味着，在这些位置上有体积大的取代基对活性是不利的。这也可以解释为什么43号GS4071的R2位置上的氢被甲基取代导致44号类似物活性急剧下降。在C1′-乙酰胺基下面的两个黄色区域（即第一组R1的位置）暗示包含这些基团是受空间限制的，位于由Arg152、Trp178、Ile222形成了疏水性的小口袋入口处，此处取代基尺寸要满足几何匹配（即把口袋口填满），这就是37～40号分子生物活性急速降低的原因，虽然42号分子的R1取代基也很小，但R2取代基的强极性作用很好地补充了其几何不匹配的弱势，这也正是产生38号分子与神经氨酸酶的结合能很低且活性不高的矛盾的原因，与分子对接结果的分析完全吻合，与之类似还有115号分子的生物活性比110号分子要低很多。

在Co MFA模型中，还有两个又长又窄的与空间等势面相关的区域，它们都靠近1号的GS4071C2′原子上的乙醛基。在Co MSIA模型中，这两个区域被预测为疏水区（见图10-34）。根据两亲性原理，增加疏水的烷基的长度将会提高NIs的抑制活性，这与60＞59＞58的实验值活性排序相符合。

对于Co MSIA中的氢键受体和供体图（见图10-40），正如预料，这些都是与对接得到抑制剂-酶键合模型相一致的。在NA—NI复合物中，Glu119和Asp151形成了作为氢键受体与NI结构中胍基或氨基部分形成氢键，而Arg292和Arg371作为氢键供体与抑制剂羧基部分形成氢键。在Co MSIA立场中的贡献因子和对接得到NIs键合位点的拓扑结构特性的相符性表明了Co MSIA结果的可靠性。Co MSIA中立体场和静电场等势图与Co MFA模型是总体相似的。

图10-40 COMSIA模型等势面图[嵌套1号分子（GS4071）]（彩图见第411页）

通过把Co MFA模型和Co MSIA模型的等势图与NI键合位点的3D拓扑结构结合起来很容易得到键合模型。这个模型不仅很好地符合键合口袋的环境特性，而且可以为进一步设计和改造NI结构提供帮助。

QSAR模型对训练集数据库有较好表现，但是对于新的数据库表现不是十分理想。但对于Co MFA和Co MSIA模型会更加完善，因为工作被分子的结构类型和叠合规则（共同骨架）所约束。接下来研究的目的之一就是得到更加普遍和有力的模型来设计新型的神经氨酸酶抑制剂，所以要进行更有说服力的验证。

除了关心模型的普遍性之外，两种QSAR模型有内部较好的验证性表现，说明这些模型是相互联系的。为了阐明这个问题，分组所建的模型可以应用Y随机化进行验证。这包含了多次随机转换Y向量之后的重复计算。经过多次重复验证，结果显示它们的q2和r2值都较低，这意味着用这些分组所得模型获得成功并不是偶然的或是依赖训练集中的结果。

10.3.3.7 结论

在本研究中，108个抑制剂类似物选择出来作为训练集分别进行分组和整体构建Co MFA和Co MSIA模型，18个带星号的抑制剂作为测试集进行外部模型验证。来自训练集和外部测试集的统计结果都表明分组所建模型较整体模型预测准确性有所提高，但受叠合规则和骨架结构的限制更多，其广泛性有所降低，即分组模型预测准确性高，整体模型更具普遍意义，它们各有优缺点，但都是成功的。当新的未知抑制剂骨架与分组内骨架结构相同，则采用分组的模型预测会获得更为准确的预测结果，当其骨架相似性不高，则采用整体进行预测，很多时候整体和分组模型预测相结合，会获得更准确的结果。Co MFA和Co MSIA力场图都为研究人员提供了有用的信息来帮助他们进行抑制剂类似物的结构改造。因此，这些3D-QSAR模型都是快速有效地预测及新药设计必不可少的通用工具，可以用来指导新型抗禽流感药物的筛选和设计，减少日渐蔓延的可能将在全球爆发的流感疫情造成的恶劣影响，也为应用支持向量积提高3D-QSAR模型预测性、更深入地了解和设计神经氨酸酶抑制剂提供重要的理论基础和支持。