β-环糊精抑制α-半乳糖苷酶的动力学模拟

邹 伟，王 敏，陶佳妮，姚 娣，李 赟，王雪松，沈汪洋,2*

(1.武汉轻工大学 食品科学与工程学院，湖北 武汉 430023； 2.大宗粮油精深加工省部共建教育部重点实验室，湖北 武汉 430023)

β-环糊精抑制α-半乳糖苷酶的动力学模拟

邹 伟1，王 敏1，陶佳妮1，姚 娣1，李 赟1，王雪松1，沈汪洋1,2*

(1.武汉轻工大学 食品科学与工程学院，湖北 武汉 430023； 2.大宗粮油精深加工省部共建教育部重点实验室，湖北 武汉 430023)

为进一步了解β-环糊精与α-半乳糖苷酶的作用机制，在前、中期采用同源建模、分子对接的基础上，采用对复合物体系进行分子动力学模拟，并计算出均方根偏差、温度、密度以及总能量四个衡量指标，对体系动力平衡的稳定程度进行评价，得到了动力学模拟轨迹是稳定的结果，再采用MMPBSA能量计算及分解的方法，对β-环糊精与α-半乳糖苷酶之间的相互作用进行能量计算，确定了两者的结合是一个自发的过程，且能够达到一个稳定的状态，并得到了β-环糊精对α-半乳糖苷酶的抑制作用主要是由范德华相互作用能和静电作用能产生的结论。

β-环糊精；α-半乳糖苷酶；分子动力学模拟；MMPBSA能量计算

1 引言

β-环糊精(β-cyclodextrin)是由7个α-D-吡喃葡萄糖单元通过α-1,4糖苷键连接而成的环状低聚糖，具有一个环外亲水、环内疏水且具有一定尺寸的立体手性空腔[1]。α-半乳糖苷酶(α-galactosidase, EC 3.2.1.22)是一种外切糖苷酶，主要作用对象是具有α-半乳糖苷结构的碳水化合物[2-3]。

咖啡豆α-半乳糖苷酶是α-半乳糖苷酶的一种,从咖啡豆中提取得到[4-7]，咖啡豆α-半乳糖苷酶能把α-半乳糖基直接转移到β-环糊精的环上，形成新的分支环糊精：α-半乳糖基-β-环糊精[8-9](α-galactosyl-β-cyclodextrin)。1992年首次报道酶法合成了这种分支环糊精，目前咖啡豆α-半乳糖苷酶酶法合成是制备α-半乳糖基-β-环糊精主要的手段。α-半乳糖基-β-环糊精具有较高的水溶性，其分支上的半乳糖基在动物组织中具有特殊的识别受体，具有广阔的应用前景。

目前有关α-半乳糖基-β-环糊精的合成及应用报道并不多，主要原因是酶法合成α-半乳糖基-β-环糊精的产率不高(22 %左右)[10-12]，对分子机制了解的缺乏严重制约了该领域的深入研究。有研究发现，环糊精可以改变酶促反应的进程，可以与酶分子的特殊基团形成复合物，激活或抑制酶的活性，从而影响反应速率或反应平衡[13-14]。而酶法合成α-半乳糖基-β-环糊精的方法中，β-环糊精为合成反应的底物，α-半乳糖苷酶(咖啡豆来源)为反应酶。在反应过程中，环糊精和咖啡豆α-半乳糖苷酶共存于同一个反应体系中，可能是环糊精自身抑制了咖啡豆α-半乳糖苷酶的活性，从而导致酶法合成产率降低。

本实验从这一现象入手，前期采用同源建模的方法[15-17]，合理构建了咖啡豆α-半乳糖苷酶的结构模型，中期采用分子对接的方法，将α-环糊精、β-环糊精及γ-环糊精与咖啡豆α-半乳糖苷酶分别进行对接，结果表明环糊精抑制了酶的活性，β-环糊精的抑制作用较为明显。本文使用MMPBSA方法计算β-环糊精与α-半乳糖苷酶的结合自由能，研究二者的结合模式，探索β-环糊精抑制α-半乳糖苷酶的作用机制，可以为降低抑制效应，提高合成效率提供理论依据。

2 材料与方法

2.1 材料与仪器

分子动力学软件包Amber 12；分子显示软件PyMOL 1.7；计算机Intel Xeon E5-2650 v2

2.2 实验方法

2.2.1 分子动力学模拟

动力学模拟的蛋白质初始构象来自实验前期同源建模的咖啡豆α-半乳糖苷酶的结构模型，在AmberTools 12的Tleap模块中补全缺失的氢原子，由ff 10力场产生立场参数。小分子的初始构象来自蛋白质数据库的晶体结构(3CGT)中的配体β-环糊精, 在实验的前期工作中已对接完毕，用sqm模块中的AM1-BCC计算分配原子电荷。将复合物置于截角八面体水盒子中心，水分子采用TIP3P模型，水盒子壁到复合物最近原子的距离为10.0 Å。加入5个Na+保持模拟体系呈电中性。

采用Amber 12的PMEMD模块对模拟体系进行四个步骤的操作: (1)能量最小化:每一步优化先执行10000步的最陡下降优化，再执行20000步的共轭梯度优化；(2)升温:在50 ps内将体系温度从0K升到300 K；(3)平衡:执行50 ps的常温300 K的温度平衡；(4)动力学模拟:对体系进行6 ns的无约束分子动力学模拟，时间步长为2 fs，坐标文件每1 ps记录一次，用于后续分析。

2.2.2 计算结合自由能

自由能计算由Amber 12的MMPBSA模块完成，从分子动力学模拟的后2 ns每10 ps取一次构象，分子动力学轨迹共200个构象用来做自由能计算。

受体和配体的结合自由能由如下方程计算：

ΔGBind=ΔGComplex-ΔGReceptor-ΔGLigand

其中ΔGBind是结合自由能，ΔGComplex,ΔGReceptor和ΔGLigand分别是复合物，受体和配体的自由能。用MMPBSA方法把它们分解为三个部分：

ΔG=ΔEgas+ΔEsol

其中为ΔEgas是气相能量焓，ΔEsol是溶解自由能。

气相能量焓ΔEgas可分解为如下两项：

ΔEgas=ΔEvdw+ΔEele

其中为ΔEvdw是范德华相互作用能，ΔEele是静电相互作用能。

溶解自由能ΔEsol可分解为如下两项：

ΔEsol=ΔEcal+ΔEsurf

其中为ΔEcal是极性溶解自由能，ΔEsurf是非极性溶解自由能。

ΔEvdw和ΔEele即MM部分；ΔEcal可通过泊松-玻尔兹曼(Poisson-Bohzmann)方程求解，即PB部分；ΔEsurf由如下经验公式计算：

ΔEsurf=γSASA +β

其中SASA是溶剂可接触表面积(solvent-accessible surface)，即SA部分。

3 结果与讨论

3.1 动力学模拟体系的平衡

对复合物体系进行了6 ns的分子动力学模拟，为了评价体系动力学平衡的稳定程度，计算了均方根偏差(Root-mean-square deviation,RMSD)、温度、密度以及总能量这四个衡量指标,分别见图1、图2、图3及图4。

图1 MD模拟中体系的均方根偏差随时间变化图

图2 MD模拟中体系的温度随时间变化图

图3 MD模拟中体系的密度随时间变化图

图4 MD模拟中体系的总能量随时间变化图

从3.5 ns后，体系趋于稳定，在3.5～6 ns内均方根偏差RMSD未见明显波动，收敛于2.74 Å，涨落范围低于0.31 Å，体系达到了动态平衡。在整个模拟过程中体系的温度、密度及总能量随时间的变化曲线较为稳定，无异常变化。上述结果表明用于下一步计算分析的动力学模拟轨迹的稳定性是可靠的。

3.2 MMPBSA能量计算及分解

采用MMPBSA方法，使用动力学模拟轨迹后2 ns的200帧构象，对β-环糊精与α-半乳糖苷酶之间的相互作用进行能量计算，结合自由能各成分的贡献均在表1中列出。

表1 MMPBSA计算结合自由能

组成成分复合物受体配体差值均值均值均值均值范德华作用能-1681．2489-1683．4017-9．8159-58．0313静电作用能-15302．5695-15063．6189-183．8766-55．0741气象能量焓-12312．0024-12144．1514-53．2335-144．6175非极性溶剂化能60．779861．95104．2188-5．3899极性溶剂化能-3581．3175-3622．9516-44．415286．0492溶解自由能-3520．5377-3561．0006-40．196480．6593结合自由能(PB法)-15832．5401-15705．1520-93．4299-33．9582结合自由能-33．9582

VDW为范德华作用能；ELE为静电作用能；GAS为气象能量焓，GAS=VDW+ELE；PBSUR为非极性溶剂化能；PBCAL为极性溶剂化能；PBSOL为溶解自由能，PBSOL=PBSUR+PBCAL；PBTOT为用PB方法计算的结合自由能，PBTOT= GAS+ PBSOL；ΔGBind为受体与配体的结合自由能，ΔGBind= PBTOT-TSTOT。因为该体系是一个能够达到动态平衡的状态，反应的熵变为零，所以最后得到的结合自由能ΔGBind即为PBTOT的值-33.9582。结合自由能是一个负值，表明两者结合的过程为一自发过程，且能够达到一个稳定的状态。从表1中可以看到体系中的范德华相互作用能有利于体系的稳定，并且相对于其它的能量项数值最小，因此为包合的主要驱动力。体系的静电作用能也有利于二者结合的形成，但相对于范德华相互作用要弱一些，可能是极性的溶解自由能抵消了静电能，有利于包合的这种趋势。非极性溶解自由能为负值，说明疏水效应也有利于两者结合的稳定性，但其绝对值较小，因此也不是形成结合的主要驱动力。

4 结论

在本研究中对β-环糊精与α-半乳糖苷酶复合物体系进行了分子动力学模拟，得到了动力学模拟轨迹是稳定的结果，再采用MMPBSA能量计算及分解的方法，对β-环糊精与α-半乳糖苷酶之间的相互作用进行能量计算，得到PBTOT的值-33.9582，确定了两者的结合是一个自发的过程，且能够达到一个稳定的状态，并得到了β-环糊精对α-半乳糖苷酶的抑制作用主要是由范德华相互作用能和静电作用能产生的。这为β-环糊精与α-半乳糖苷酶的相互作用机制的研究和开发该酶的一些特殊作用提供了相应的依据。

[1] 宫斌, 黄据芹, 王怀友. β-环糊精对中性红及荧光素包合性能的研究 [J].化学分析计算, 2010,02:19-23.

[2] 郝桂娟, 张凯, 王学智. α-半乳糖苷酶的研究进展 [J].中国畜牧兽医, 2013, 03:149-154.

[3] 陶佳妮, 李赟, 姚娣. 咖啡豆α-半乳糖苷酶的同源建模及分析 [J].食品工业科技, 2016,01:326-330.

[4] Marraccini P,Rogers W J. Biochemical and molecular characterization of from coffee beans[J].Plant Physiology and Biochemistry, 2005,43: 909-920.

[5] 沈汪洋, 金征宇. 发芽咖啡豆α-半乳糖苷酶的性质研究 [J]. 安徽农业科学, 2011,19:11958-11960.

[6] 侯廷军, 徐筱杰. 基于分子动力学模拟和连续介质模型的自由能计算方法 [J]. 化学进展, 2004,02:153-158.

[7] 沈汪洋, 金征宇. 植物多酚对α-半乳糖苷酶提取的影响 [J]. 食品工业科技, 2008,04:90-92.

[8] 沈汪洋, 金征宇, 高虹. α-半乳糖苷酶法合成α-半乳糖基-β-环糊精的研究 [J]. 湖北农业科学, 2011,19:4041-4043.

[9] 沈汪洋, 金征宇. 半乳糖基-β-环糊精的合成与结构研究 [J]. 食品科学, 2011,12:57-59.

[10] 崔波, 金征宇. 麦芽糖基(α-1-6)β-环糊精的酶法合成和结构鉴定 [J]. 高等学校化学学报, 2007,02:283-285.

[11] Harak , Fujitak , Kuwaharad, et al. Galaetosylation of cyclodextrins and branched cyclodextrins by α-galactosidase [J]. Bioscienee, Bioteelmology and Biochem istry, 1994, 4:652—659.

[12] Koizumi K, Tanimoto T, Fujita K. et al. Preparation, iIsolation, and characterization of novel heterogeneous branched cyclomalto-oligosaccharides having β-galactosyl residue(s) on the side chain [J].Carbohydrate Research, 1993,238:75一91.

[13] Seosj, Kimsh, Sasagawat, et al. Delivery of all trans-retinoic acid(RA) to hepatocyte cell line from RA/galactosyl α-cyclodextrin inclusion complex [J]. European Journal of Pharm aceutics and Biopharm aceutics, 2004,58: 681-687.

[14] 于博, 张焕新, 金征宇等. 环糊精及其衍生物对普鲁兰酶抑制机理的研究 [J]. 食品工业科技, 2011,03:109-111+115.

[15] 洪玮. 糖苷水解酶Fncel5A和GB1多肽的分子动力学模拟 [D].合肥：中国科学技术大学生物信息学，2012.

[16] 王加磊, 孟现美, 张庆刚. HIV-1蛋白酶异位抑制剂体系的长时问分子动力学模拟 [J]. 山东科学, 2013,02:1-6.

[17] 柯国涛. 用分子模拟方法研究HIV-1整合酶与病毒DNA的相互作用 [D].北京：北京工业大学生物物理学.2008.

Dynamics simulations of inhibition of β-cyclodextrin on α-galactosidase

ZOUWei1,WANGMin1,TAOJia-ni1,YAODi1,LIYun1,WANGXue-song1,SHENWang-yang1,2*

(1.School of Food Science and Engineering, Wuhan Polytechnic University, Wuhan 430023,China; 2. Key Laboratory of the Deep Processing of Bulk Grain and Oil Authorized by Ministry of Education, Wuhan 430023,China)

In order to further understand the mechanism of action of β-cyclodextrin and α-galactosidase, this test used homology modeling and molecular docking based on the molecular dynamics simulation of complex system in the early and middle stage, and calculated four indicators include rmsd deviation, temperature, density and total energy, and evaluated the degree of stability of the system dynamic equilibrium, aquired a stable result of dynamics simulation trajectory , then by using the MMPBSA energy calculation and decomposition method for calculating the energy of interaction between β-cyclodextrin and α-galactosidase, determined the combination of β-cyclodextrin and α-galactosidase was a spontaneous process, and could achieve a steady state, and aquired the result that the inhibitory effect of β-cyclodextrin on α-galactosidase was mainly composed of Fan Dehua interaction energy and the electrostatic interaction.

β-cyclodextrin；α-galactosidase；molecular dynamics simulation；MMPBSA energy calculation

2016-12-23

邹伟(1992-)，男，硕士研究生。E-mail：530294379@qq.com.

*通信作者：沈汪洋(1978-)，男，副教授，博士。E-mail: whwangyangshen@126.com

国家自然科学基金委员会；国家自然科学基金青年科学基金项目(31301415)

2095-7386(2017)01-0022-04

10.3969/j.issn.2095-7386.2017.01.004

TS 201.2

A