多肽抑制剂抑制淀粉质多肽42构象转换的分子动力学模拟和结合自由能计算
2012-12-11董晓燕都文婕刘夫锋
董晓燕 都文婕 刘夫锋
(天津大学化工学院生物工程系,天津300072;天津大学,系统生物工程教育部重点实验室,天津300072)
多肽抑制剂抑制淀粉质多肽42构象转换的分子动力学模拟和结合自由能计算
董晓燕 都文婕 刘夫锋*
(天津大学化工学院生物工程系,天津300072;天津大学,系统生物工程教育部重点实验室,天津300072)
应用分子动力学模拟和结合自由能计算方法研究了多肽抑制剂KLVFF、VVIA和LPFFD抑制淀粉质多肽42(Aβ42)构象转换的分子机理.结果表明,三种多肽抑制剂均能够有效抑制Aβ42的二级结构由α-螺旋向β-折叠的构象转换.另外,多肽抑制剂降低了Aβ42分子内的疏水相互作用,减少了多肽分子内远距离的接触,有效抑制了Aβ42的疏水塌缩,从而起到稳定其初始构象的作用.这些抑制剂与Aβ42之间的疏水和静电相互作用(包括氢键)均有利于它们抑制Aβ42的构象转换.此外,抑制剂中的带电氨基酸残基可以增强其和Aβ42之间的静电相互作用(包括氢键),并降低抑制剂之间的聚集,从而大大增强对Aβ42构象转换的抑制能力.但脯氨酸的引入会破坏多肽的线性结构,从而大大降低其与Aβ42之间的作用力.上述分子模拟的结果揭示了多肽抑制剂KLVFF、VVIA和LPFFD抑制Aβ42构象转换的分子机理,对于进一步合理设计Aβ的高效短肽抑制剂具有非常重要的理论指导意义.
分子动力学模拟;阿尔茨海默病;淀粉质多肽;多肽抑制剂;构象转换
1 引言
阿尔茨海默病(AD)是一种以进行性记忆和认知功能损伤为特征的退行性神经系统疾病,目前确切病因尚不清楚.临床表现为认知和记忆功能不断恶化,日常生活能力进行性减退,并伴有各种神经症状和行为障碍.1,2AD是老年人的多发病和常见病,尤其是随着世界人口老龄化的加剧,已成为人类的第四号杀手,将对社会和经济发展造成巨大的负面影响.
目前AD的病理变化存在多种假说,如Aβ假说、胆碱能损害假说和炎症反应假说等.3,4由于AD的主要病理特征为脑神经细胞外出现β-淀粉样多肽(Aβ)聚集形成的老年斑和脑神经细胞内Tau蛋白异常聚集形成的神经纤维缠结等,因此Aβ假说得到了深入的发展.越来越多的证据表明,Aβ是引起AD的原发性病理因子.5Aβ假说认为,Aβ的错误折叠和聚集形成以β-折叠为主的聚集体,会引发一系列复杂的反应,包括突触变化和炎症反应等,最终出现神经元功能失调、斑块形成、神经原纤维缠结等病理现象.6
Aβ是淀粉质前体蛋白的水解产物,通常为含有39-43个氨基酸残基的多肽,以Aβ40和Aβ42为主,其中又以Aβ42的聚集沉淀趋势最强.7研究表明,Aβ分子首先从初始的α-螺旋转变成以β-折叠为主的二级结构,然后相互聚集形成寡聚体和纤维.8-10早期研究认为,只有Aβ聚集形成的淀粉质斑或纤维才具有神经毒性,而近来越来越多的研究证明,可溶性Aβ寡聚体和原细纤维体的神经毒性更大.11但无论是Aβ寡聚体还是成熟纤维(或者是两者),只要阻止Aβ的构象转变或聚集就可从根本上消除Aβ所引发的神经毒性.因此,目前对于AD的治疗策略主要着眼于Aβ聚集抑制剂的开发.
现有的Aβ聚集抑制剂主要有:有机小分子、蛋白质和多肽.12-15由于多肽类抑制剂具有较好的生物相容性、毒副作用小、分子量小、易于透过血脑屏障等优点,因此近年来备受关注.16-18目前研究开发的多肽类抑制剂大多为Aβ的片段(KLVFF和VVIA)或其衍生物(LPFFD).19-21大量研究表明,这些多肽类抑制剂可以与Aβ发生结合从而阻止Aβ聚集和毒性的产生.22-26但由于Aβ构象转换迅速,用目前的实验技术无法检测,导致多肽类抑制剂抑制Aβ构象转换的作用机理无法用现有的实验方法进行研究.然而分子动力学模拟的出现弥补了现有实验研究的不足,已广泛应用于Aβ构象转换和聚集及其抑制的研究中.27例如,Xu等10利用全原子分子动力学模拟探明了Aβ40从初始的α-螺旋到β-折叠构象转换的分子机理,并发现四个甘氨酸对于上述的构象转换起到非常重要的作用,将其突变成丙氨酸后就能有效地抑制其构象转换.Wei等28利用复制交换分子动力学模拟解析了Aβ25-35聚集成二聚体和纤维的结构形态.
此外,分子动力学模拟也常用来研究抑制剂抑制Aβ构象转换和聚集的分子机理.作者前期利用分子动力学模拟解析了海藻糖抑制Aβ16-22聚集和茶多酚EGCG抑制Aβ42构象转换的分子机理.9,29Yang等19以分子对接获得的Aβ42-LPFFD复合物为初始结构,利用全原子分子动力学模拟解析了多肽抑制剂LPFFD抑制Aβ42构象转换的分子机理.研究结果表明LPFFD通过破坏Aβ42的盐桥D23-K28来抑制Aβ42的构象转换.Viet等30利用分子动力学模拟和结合自由能计算方法解析了短肽抑制剂KLVFF和LPFFD影响Aβ16-22聚集和Aβ40构象转换的分子机理,并解析了抑制剂与Aβ之间的亲和力与其抑制能力之间的关系.但该研究没有解析短肽抑制剂KLVFF和LPFFD抑制Aβ42构象转换的分子机理.现有大量实验研究表明,位于Aβ42 C端的多肽片段VVIA能够有效抑制Aβ的聚集沉淀.31但至今还没有利用分子动力学模拟解析位于Aβ42 C端的多肽VVIA抑制Aβ42构象转换的分子机理.此外,从这三种多肽抑制剂的氨基酸序列可以看出,KLVFF和LPFFD除了含有疏水性氨基酸残基外还含有带电氨基酸残基,而VVIA仅含有疏水氨基酸残基.此外,LPFFD中脯氨酸残基的引入对抑制Aβ42构象转换的影响尚不清楚.因此,利用分子动力学模拟和结合自由能计算方法研究比较这三种多肽抑制剂抑制Aβ42的构象转换以及它们和Aβ42之间的亲和作用,进一步解析多肽抑制剂中带电氨基酸和脯氨酸的引入对于短肽抑制剂抑制Aβ构象转换的影响有助于理解短肽抑制剂抑制Aβ构象转换和聚集之间的相互关系,对于多肽类AD抑制剂的合理设计具有重要的理论指导和实际意义.
本文以Aβ42为研究对象,采用全原子模型的分子动力学模拟和结合自由能计算方法,考察了常见的三种多肽类抑制剂KLVFF、VVIA和LPFFD对Aβ42构象变化的影响,通过分析在不同抑制剂存在下Aβ42的二级结构变化、Aβ42的侧链接触、抑制剂和Aβ42之间的接触距离、抑制剂和Aβ42之间的结合自由能和氢键分析,来进一步解析三种多肽抑制剂抑制Aβ42构象转变的作用机理.
2 模拟方法
2.1 模拟体系构建
Aβ42的初始结构从蛋白质数据库(PDB)获得, PDB号为1IYT.32Aβ42的初始结构含有两个α-helix,分别定义为helix-1(8-25位残基)和helix-2(28-39位残基)(图1A).其氨基酸序列为:DAEFRHDSGYEVHHQKLVFFAEDVGSNKGAIIGLMVGGVIA.
三种多肽类抑制剂KLVFF、VVIA和LPFFD的化学结构式如图1(B,C,D)所示,它们均利用SYBYL软件绘制得到.这些抑制剂的N末端用乙酰基修饰,C末端用亚氨基修饰.
首先将Aβ42放入一个立方体盒子(8 nm×8 nm× 8 nm)中,然后将多肽抑制剂放入该盒子,使其随机分布,并与Aβ42之间的最小距离为0.5 nm.随后将盒子中加满水,并去除与多肽或抑制剂重叠的水分子;最后,用不同数目的Na+或Cl-替代同等数目的水分子使模拟体系保持电中性.其中,水溶液中添加3个Na+,抑制剂KLVFF溶液中添加5个Cl-, VVIA溶液中添加3个Na+,LPFFD溶液中添加11个Na+.采用三次能量最小化优化该体系:首先,固定Aβ42和抑制剂分子的结构和位置不变,仅让水分子的结构和位置发生变化;其次,仅固定Aβ42的结构和位置,使抑制剂和水分子的结构和位置发生变化;最后,使体系内所有分子都可以自由运动,利用能量最小化来优化上述模拟体系.
2.2 分子动力学模拟
图1 Aβ42的初始结构(A)与三种多肽抑制剂KLVFF(B),VVIA(C)和LPFFD(D)的化学结构Fig.1 Initial structure ofAβ42(A)and chemical structures of the three peptide inhibitors KLVFF(B),VVIA(C),and LPFFD(D)The main chain ofAβ42 is shown by a NewCartoon model and its side chains are represented by a line.
将上述优化好的体系进行分子动力学模拟,每个体系模拟时间为80 ns.各模拟体系的详细信息,如表1所示.模拟过程中,为了使模拟结果更具有代表性,对相同抑制剂条件下的同一体系通过改变动力学模拟参数中的随机数,使体系中各原子获得不同的初始速度,从而得到三条不同的模拟轨迹.对于图2-7,除了图2,3和5之外,文中的所有分析都基于三次模拟结果的平均.图2和图3是利用GROMACS软件自带的程序基于一条模拟轨迹计算获得.图5中Aβ42和多肽抑制剂之间相互作用的典型构象利用VMD1.9.1软件制作.
表1 模拟体系的参数Table 1 Parameters of the simulation systems
本文采用GROMACS 4.0.5分子动力学模拟软件,33选择GROMOS 96力场,34水分子采用单点电荷(SPC)模型.35利用Lenard-Jones函数计算范德华作用力,非键截断距离设定为1.4 nm,非键作用原子列表每5个步长更新一次;利用particle-mesh Ewald方法计算静电相互作用.36采用Verlet leapfrog算法37对每一步的运动方程进行求解,经过积分得到新时刻各原子的坐标,积分步长为2 fs.同时,采用LINCS算法38对所有成键原子之间的相对距离进行固定,以减小计算量.所有模拟均在等温等压(NPT)系综中进行,温度300 K,通过V-rescale方法控制温度,时间常数设为0.1 ps.压力为1.013×105Pa,压力控制采用Berendsen方法,39压力耦合常数为0.5 ps.所有分子动力学模拟计算均在曙光TC2600刀片服务器上完成.
2.3 数据分析方法
2.3.1 Aβ42的二级结构分析
采用GROMACS软件中的附属程序do_dssp计算Aβ42的二级结构,考察多肽的二级结构随模拟时间的变化.40
2.3.2 Aβ42侧链接触图计算
采用GROMACS软件中的附属程序g_mdmat计算多肽侧链接触图.截断距离设为0.9 nm.
2.3.3 抑制剂和Aβ42之间接触数和接触距离的计算
抑制剂和Aβ42之间的接触数和接触距离均采用GROMACS附属程序g_mindist进行计算.其中,计算抑制剂和Aβ42之间的接触数时所用的截断距离为0.6 nm.
2.3.4 结合自由能计算
抑制剂和Aβ42之间的结合自由能(ΔGbind)采用分子力学-泊松-玻尔兹曼-溶剂可及表面积(MMPBSA)方法进行计算分析.41,42计算公式如下:
其中,ΔEMM为分子内能,通过分子力学方法计算获得,ΔGsol为溶剂化自由能,ΔS为熵.分子内能主要包括范德华作用能(ΔEvdw)和静电相互作用能(ΔEelec).而溶剂化自由能主要包含极性溶剂化自由能(ΔGPB)和非极性溶剂化自由能(ΔGnps).极性溶剂化自由能利用泊松-玻尔兹曼方程计算获得,溶质和溶剂的介电常数分别为1和80.非极性溶剂化自由能由溶剂空穴作用项和溶质-溶剂范德华作用项组成,如公式(2)所示.
其中,γ和b分别为2.27 kJ·mol-1·nm-2和3.85 kJ· mol-1.SASA为溶剂可及表面积.
为了详细分析抑制剂和Aβ42之间的结合自由能,作者将抑制剂和Aβ42之间的自由能分解为极性作用能ΔGPBelec和非极性作用能ΔGnonpolar.在接下来的分析中,ΔGPBelec主要包含静电相互作用和极性溶剂化自由能,因此可将ΔGPBelec认为是抑制剂和Aβ42之间的静电作用贡献.而ΔGnonpolar主要包含范德华作用(ΔEvdw)和非极性溶剂化自由能,因此可以将其作为抑制剂和Aβ42之间的疏水作用贡献.
2.3.5 氢键分析
本文中抑制剂和Aβ42之间的氢键作用利用GROMACS附属程序g_hbond进行分析.当供体原子D与受体原子A之间的距离小于0.35 nm,且D―H―A之间的夹角大于120°时,即认为供体原子D和受体原子A之间形成氢键.最理想的氢键的键角是D―H―A为180°.
3 结果与讨论
3.1 Aβ42的二级结构分析
图2为Aβ42在水和不同多肽抑制剂溶液中的二级结构随模拟时间的变化.在模拟过程中,利用不同的颜色来表示多肽的二级结构.其中图2A为Aβ42在水溶液条件下的二级结构变化.从该图可以看出,在水溶液中随着模拟的进行,Aβ42的N和C端的初始α-螺旋(α-helix)结构被破坏,而转换成转角(turn)、无规卷曲(coil)和β-折叠(β-sheet)等结构;但在抑制剂KLVFF存在情况下,虽然Aβ42 C端的25-42位残基的初始α-螺旋被破坏并转变成转角,弯曲(bend)和无规卷曲等结构,并在之后的模拟过程中一直稳定存在;但其N端6-24位残基的α-螺旋结构得到了很好的保持(图2B所示).其主要原因是抑制剂上的带正电荷的赖氨酸(K)能够和6-24位残基中的带负电荷氨基酸(D7,E11,E23和D24)之间形成较强的静电相互作用,从而维持了N端6-24位残基的初始α-螺旋构象.图2C为Aβ42在抑制剂VVIA存在情况下的二级结构变化.从该图可看出,Aβ42 N端8-12位氨基酸从45 ns开始逐步转化为富含5-螺旋(5-helix),无规卷曲和弯曲等构象,并在随后的模拟中始终稳定存在;C端的构象在5 ns左右被破坏,出现少量的β-折叠和转角结构,但稳定性较差, 40 ns开始逐步被转角和弯曲等结构所替代;13-25位残基在模拟过程中始终维持着初始的α-螺旋结构(图2C所示).图2D为抑制剂LPFFD抑制Aβ42二级结构变化的情况.从该图可以看出,N末端的构象维持较好,30 ns左右8-13位残基开始出现转角结构,并一直稳定维持到80 ns;22-28位残基的构象发生转变,在35 ns以后出现较多的5-螺旋结构.此外,C端32-37位残基在模拟的80 ns内均维持初始的α-螺旋结构,但C末端38-40位残基在模拟初始时就出现转角结构且在80 ns的模拟过程中一直稳定存在;且在模拟过程中没有出现β-折叠二级结构.从图2D可以看出,抑制剂LPFFD除了能够抑制Aβ42 N端2-7位残基的构象转换,还能够维持C端32-37位残基的初始α-螺旋构象.其主要原因是抑制剂上的带负电荷天冬氨酸(D)可以和5位精氨酸(R)和28位赖氨酸(K)之间形成较强的静电相互作用,从而维持其相关区域的初始构象.上述模拟结果表明,三种多肽抑制剂都能够抑制Aβ42 N端的构象转换,维持Aβ42的初始α-螺旋构象,并抑制β-折叠结构的形成.此外,多肽抑制剂上的带电氨基酸残基对于抑制Aβ42的构象转换和维持其初始构象具有非常重要的作用.
图2 Aβ42的二级结构随模拟时间的变化Fig.2 Secondary structures as a function of simulation time forAβ42(A)without inhibitor;(B)KLVFF;(C)VVIA;(D)LPFFD.Secondary structure is color-coded.
3.2 Aβ42侧链接触图分析
多肽侧链的接触图(contact map)常被用来表征多肽分子内疏水相互作用.43,44图3为Aβ42在水溶液和不同抑制剂溶液中的侧链接触图.Aβ42中残基侧链之间接触距离的大小用图中不同的颜色表示,其中白色和黑色分别表示Aβ42侧链之间的距离为0和0.9 nm.颜色越接近于白色,表示相关残基侧链之间的距离越近;由于对角线表示的是相同残基侧链之间的接触距离,其值为0,因此在图中显示白色.
图3 Aβ42的侧链接触图Fig.3 Contact maps of the side chains ofAβ42(A)no inhibitor,(B)KLVFF,(C)VVIA,(D)LPFFD;Each square provides the mean average distance between any atoms of two residues ofAβ42.
可以看出,在水溶液中,Aβ42的远距离侧链之间的接触数比较多,这主要是由于Aβ42在分子内疏水作用下发生塌缩使其N末端和C末端靠近的结果(图3A),这与作者早期的研究结果44类似.而在三种多肽抑制剂存在的情况下Aβ42的远距离侧链的接触数明显减少,如图3B,3C和3D所示.这证明了抑制剂的存在抑制了Aβ42的疏水塌缩,从而降低了N和C末端之间的接触数.此外,图3B中Aβ42的氨基酸侧链之间的远距离和近距离接触数都比图3C和3D少,说明KLVFF抑制剂对Aβ42疏水塌缩的抑制作用效果要好于另两个抑制剂VVIA和LPFFD.由于Aβ42是疏水性很强的多肽,其分子内的疏水相互作用在其构象转换和聚集过程中发挥重要的作用.因此,抑制Aβ42的疏水塌缩,即抑制Aβ42分子内的疏水相互作用,是抑制其构象转换的前提.因此,多肽抑制剂的存在减少了Aβ42分子内远距离的接触,有效抑制了Aβ42的疏水塌缩,降低了Aβ42分子内的疏水相互作用,从而起到稳定其初始构象的作用.
3.3 多肽抑制剂和Aβ42之间的接触数和接触距离分析
为了进一步分析三种抑制剂抑制Aβ42构象转换的能力,计算了三种抑制剂和Aβ42的接触数和接触距离.抑制剂和Aβ42之间的接触数和接触距离可以用于表征抑制剂和Aβ42之间的作用力大小.抑制剂与Aβ42之间的接触数多和距离短,表明抑制剂和Aβ42之间的作用力强;反之,则说明抑制剂和Aβ42之间的作用力弱.图4为Aβ42和抑制剂之间的接触数随模拟时间的变化.从图4可以看出,抑制剂KLVFF和Aβ42之间的接触数上升最快,在模拟初始的1 ns内就急剧上升到3300左右,并且随着模拟的进行而缓慢增加,在15 ns处达到最大(~6000),而后会缓慢降低至4000左右.而抑制剂VVIA和LPFFD与Aβ42之间的接触数在模拟的初始阶段增加比较平缓,且两种抑制剂和Aβ42之间的接触数量相当.40 ns后,LPFFD和Aβ42的接触数大于VVIA和Aβ42之间的接触数.在模拟的后期, LPFFD和Aβ42之间的接触数继续增加到3600以上.图5为三种多肽抑制剂在80 ns时结合Aβ42的典型构象.从该图可以看出,三种抑制剂除了和Aβ42发生直接相互作用外,三种抑制剂自身团聚发生相互作用.因此,如果能够大大降低多肽抑制剂自身之间的团聚可以有利于多肽抑制剂和Aβ42之间的相互作用,从而大大提高其抑制能力.
图4 Aβ42和抑制剂之间的接触数随模拟时间的变化Fig.4 Number of contacts betweenAβ42 and peptide inhibitors as a function of simulation time
图6 Aβ42和抑制剂之间的平均接触距离的概率分布Fig.6 Distribution of the average contact distances betweenAβ42 and peptide inhibitors
图6为Aβ42和抑制剂之间的平均接触距离的概率分布.从该图可以看出在抑制剂KLVFF存在下,Aβ42中氨基酸残基与抑制剂间的距离分布范围较小,主要位于0.13-0.18 nm,且分布的峰值位于0.13-0.17 nm之间.说明KLVFF与Aβ42所有残基之间的作用情况较稳定,作用区域分布无特异性,抑制剂分子分布在Aβ42的整个表面,如图5A所示. Aβ42与抑制剂LPFFD的接触距离的概率分布的范围较宽,峰值位于0.14-0.18 nm.与抑制剂KLVFF相比较,峰值向远距离方向移动,且概率值分布较平均,说明抑制剂LPFFD与Aβ42之间结合的紧密程度比抑制剂KLVFF差(图5C).主要原因是脯氨酸的引入破坏了该肽段的线性结构,即仅有三个残基FFD能够与Aβ42发生相互作用,从而大大降低了LPFFD和Aβ42之间的结合作用.从而使抑制剂LPFFD和Aβ42之间的接触距离比抑制剂KLVFF长.抑制剂VVIA和Aβ42的氨基酸之间的接触距离的最高峰虽然也位于0.16 nm,但其整体概率分布范围很宽,几乎遍布0.13-0.27 nm,个别残基与抑制剂之间的距离达到0.27 nm.这说明抑制剂与部分氨基酸之间的作用力较弱,或不易结合,其原因可能有以下两点:其一,抑制剂VVIA为四肽而另外两种抑制剂为五肽,序列较短使得抑制剂无法覆盖Aβ42分子的整个序列,如图5B所示;其二,抑制剂VVIA全部由疏水性残基组成,在模拟过程中出现了明显的抑制剂自聚现象(图5B),从而影响了抑制剂与Aβ 42之间的结合,进而影响了VVIA的抑制效果.相反,由于抑制剂KLVFF和LPFFD中含有带电氨基酸,从而大大降低了它们的自聚能力.综上所述,抑制剂KLVFF与Aβ42之间的作用力要强于另外两种多肽类抑制剂VVIA和LPFFD.
图5 80 ns时Aβ42和多肽抑制剂之间相互作用的典型构象Fig.5 Representative structures of interaction ofAβ42 and peptide inhibitors at 80 ns(A)KLVFF,(B)VVIA,(C)LPFFD.The three snapshots extracted from the corresponding MD trajectories.For clarity,water molecules are not shown.The main chain ofAβ42 is shown by a NewCartoon model,and the peptide inhibitors are represented by a Licorice model.
3.4 抑制剂和Aβ42之间的相对结合自由能分析
为了进一步解析三种抑制剂抑制Aβ42构象转换的能力,作者利用MM-PBSA方法计算了抑制剂和Aβ42之间的结合自由能.利用模拟过程平衡阶段的最后10 ns(70-80 ns)的模拟轨迹进行分析,并根据分子间作用能(范德华作用能和静电作用能)和溶剂化作用能(溶剂化非极性作用能和溶剂化静电作用能)对多肽抑制剂-Aβ42的结合自由能进行分解,结果如表2所示.本文的主要目的是分析比较三种多肽抑制剂抑制Aβ42构象转换的能力.而计算结合自由能时所用的抑制剂-Aβ42复合物的结构是从分子动力学模拟后期稳定阶段(70-80 ns)获得的.此时,抑制剂和Aβ42的构象基本没有什么变化.因此,抑制剂和Aβ42之间的分子间作用能和溶剂化作用能对于多肽抑制剂抑制Aβ42构象转换的能力密切相关,而熵的关系不大.例如,Viet等30计算所得的多肽抑制剂KLVFF和LPFFD与Aβ40之间的熵基本相等,均为138.2 kJ·mol-1·K-1.该现象也被很多研究发现.45因此多肽抑制剂-Aβ42之间熵的贡献没有计算.本文中多肽抑制剂和Aβ42之间的自由能为相对结合自由能(忽略熵),而不是绝对结合自由能.另外,已有大量文献表明,相对结合自由能完全能够用来表征分子之间的结合作用力强弱.46-49从表2可以看出三种抑制剂和Aβ42之间的非极性相互作用能和静电相互作用能均为负值,说明在模拟过程中,Aβ42和多肽抑制剂之间的疏水和静电相互作用均有利于抑制剂和Aβ42结合,从而抑制Aβ42的构象转换.
从表2可以看出Aβ42与抑制剂KLVFF形成复合物的ΔGbind值最低,为-494.0 kJ·mol-1,说明Aβ42与KLVFF的结合作用最强,结合体最稳定.而与VVIA和LPFFD形成复合物时产生的ΔGbind相似.说明这两种抑制剂抑制Aβ42的构象转换能力相差不大.上述结果与侧链接触(图3)以及抑制剂和Aβ42之间的接触数(图4)和接触距离(图6)分析的结果一致.
表2 多肽抑制剂-Aβ42复合物的结合自由能分解Table 2 Binding free energy components of peptide inhibitors-Aβ42 complex
对比表中各种作用能发现,ΔEvdW,ΔGnps和ΔEelec均为负值,说明多肽抑制剂和Aβ42之间的范德华作用、非极性溶剂化自由能和静电作用对多肽抑制剂和Aβ42之间的结合均起促进作用.Aβ42与抑制剂KLVFF形成复合物时,ΔGPBelec在ΔGbind中占70%,因此静电作用在抑制剂KLVFF结合Aβ42的贡献大.而抑制剂LPFFD与Aβ42形成复合物时,ΔGnonpolar在ΔGbind中占主要地位(>78%),即疏水作用对抑制剂LPFFD与Aβ42之间的结合能的贡献大.对比抑制剂KLVFF和LPFFD与Aβ42之间的静电相互作用能可以发现,KLVFF和Aβ42之间的静电相互作用远大于LPFFD.主要是因为抑制剂KLVFF中的赖氨酸(K)带正电荷,可以与Aβ42上的带负电荷的氨基酸残基(3个天冬氨酸和3个谷氨酸残基)产生较强的静电相互作用.此外,带负电荷的氨基酸残基较多且分布较广,从而有利于和抑制剂KLVFF中的赖氨酸发生作用.因此Aβ42和KLVFF之间的静电作用较大.与此相反,抑制剂LPFFD带负电荷,虽然可以与Aβ42上的带正电荷的氨基酸残基(1个精氨酸和2个赖氨酸残基)产生较强的静电相互作用.但由于Aβ42中的带正电荷的氨基酸分布集中,且Aβ42整体也带负电荷,使得LPFFD和Aβ42之间的静电作用能小于KLVFF.
综上所述,多肽抑制剂与Aβ42之间的静电和疏水相互作用均有利于它们之间的结合.抑制剂KLVFF与Aβ42结合过程中分子间作用以静电作用为主,而抑制剂VVIA和LPFFD与Aβ42之间的作用以疏水相互作用为主.
3.5 抑制剂和Aβ42之间的氢键分析
早期研究表明,Aβ42构象转换的推动力主要有疏水和静电(包括氢键)相互作用.29此外氢键网络还是维持寡聚体和纤维规则结构的重要作用形式.50图7显示了Aβ42和三种抑制剂之间的氢键数量随模拟时间的变化.比较三种不同的抑制剂发现,Aβ42与抑制剂KLVFF相互作用时生成的氢键个数明显多于与另两种抑制剂之间产生的氢键个数.在模拟的初始阶段抑制剂KLVFF和Aβ42之间的氢键个数迅速上升,在25-65 ns期间氢键达到13个左右,65 ns后氢键数量略有下降,最终氢键个数维持在12个左右.抑制剂VVIA和LPFFD与Aβ42相互作用形成的氢键个数在整个模拟过程中缓慢增长,到70 ns后VVIA与Aβ42之间的氢键个数趋于稳定(约7个),而LPFFD和Aβ42之间的氢键个数在5个以下.
图7 Aβ42和抑制剂之间的氢键数量随模拟时间的变化Fig.7 Number of H-bonds between peptide inhibitors and Aβ42 as a function of simulation timeData are fitted by a logarithm function.
KLVFF是Aβ42的疏水核心,大量的实验和分子模拟结果显示它与Aβ42具有非常强的结合能力.30,51KLVFF含有四个疏水性残基和一个带电残基.该抑制剂除了通过主链的羰基和亚氨基与Aβ42形成氢键外,带电残基K含有较多的氢键供体和受体,它与Aβ42之间也能形成大量的氢键.相反,抑制剂VVIA的所有残基均为疏水性氨基酸,它们只能通过主链的羰基和亚氨基与Aβ42之间形成氢键,所以它与Aβ42之间的氢键数量要小于KLVFF.抑制剂LPFFD是KLVFF的衍生物,其氨基酸组成和KLVFF类似,均含有四个疏水性残基和一个带电残基.其与Aβ42之间的氢键数量要远小于KLVFF的主要原因是脯氨酸(P)的引入,破坏了该肽段的线性结构(图1D).从而使其与Aβ42之间的氢键数远小于KLVFF.在抑制剂LPFFD中的脯氨酸,可以抑制其它Aβ分子继续聚集,但同时也降低了其与Aβ的结合能力.综上所述,多肽抑制剂和Aβ42形成氢键数量的能力与其抑制Aβ42构象转换的能力正相关.多肽抑制剂中含有带电残基有利于其与Aβ42形成氢键,从而更有利于抑制其构象转换和聚集.
4 结论
利用全原子分子动力学模拟方法考察了三种多肽抑制剂KLVFF、VVIA和LPFFD对Aβ42构象变化的抑制作用,并分析了三种抑制剂与Aβ42结合过程的差异,得到如下结论:三种多肽抑制剂均能够有效抑制Aβ42的构象变化,稳定其初始的α-螺旋结构,尤其是N端α-螺旋构象,从而阻止了β-折叠构象的生成.多肽抑制剂还使Aβ42分子内远距离的接触减少,有效抑制了Aβ42的疏水塌缩,从而降低了Aβ42分子内的疏水相互作用,并起到稳定其初始构象的作用.抑制剂与Aβ42之间的疏水和静电作用(包括氢键)均有利于它们抑制Aβ42的构象转换.由于抑制剂KLVFF和Aβ42之间的疏水和静电作用均强于其它两种多肽抑制剂(VVIA和LPFFD),因此抑制剂KLVFF抑制Aβ42构象转换的效果最好.另外,抑制剂中含有的带电氨基酸残基可以增强其和Aβ42之间的静电相互作用(包括氢键),并降低了抑制剂之间的聚集,从而对Aβ42构象转换的抑制效果较好.在多肽抑制剂中引入脯氨酸会破坏多肽的线性结构,从而大大降低多肽抑制剂与Aβ42之间的结合作用力.
(1) Yan,H.;Jiang,F.C.Acta Phys.-Chim.Sin.2006,22,359. [鄢 浩,姜凤超.物理化学学报,2006,22,359.]doi:10.3866/ PKU.WHXB20060321
(2) Jakob-Roetne,R.;Jacobsen,H.Angew.Chem.Int.Edit.2009, 48,3030.doi:10.1002/anie.200802808
(3) Hardy,J.A.;Higgins,G.A.Science 1992,256,184.doi: 10.1126/science.1566067
(4)Luo,Z.W.;Wang,D.D.;Lai,L.H.;Xu,X.J.;Li,C.X.Acta Phys.-Chim.Sin.1995,11,419.[骆兆文,王丹丹,来鲁华,徐筱杰,李崇熙.物理化学学报,1995,11,419.]doi:10.3866/ PKU.WHXB19950507
(5)Karran,E.;Mercken,M.;De Strooper,B.Nat.Rev.Drug Discov.2011,10,698.doi:10.1038/nrd3505
(6) Hardy,J.;Selkoe,D.J.Science 2002,297,353.doi:10.1126/ science.1072994
(7) Esler,W.P.;Wolfe,M.S.Science 2001,293,1449.doi:10.1126/ science.1064638
(8) DaSilva,K.A.;Shaw,J.E.;McLaurin,J.Exp.Neurol.2010, 223,311.doi:10.1016/j.expneurol.2009.08.032
(9) Liu,F.F.;Ji,L.;Dong,X.Y.;Sun,Y.J.Phys.Chem.B 2009, 113,11320.doi:10.1021/jp905580j
(10) Xu,Y.;Shen,J.;Luo,X.;Zhu,W.;Chen,K.;Ma,J.;Jiang,H. Proc.Natl.Acad.Sci.U.S.A.2005,102,5403.doi:10.1073/ pnas.0501218102
(11) Rochet,J.C.;Lansbury,P.T.,Jr.Curr.Opin.Struct.Biol.2000, 10,60.doi:10.1016/S0959-440X(99)00049-4
(12) Mason,J.M.;Kokkoni,N.;Stott,K.;Doig,A.J.Curr.Opin. Struct.Biol.2003,13,526.doi:10.1016/S0959-440X(03) 00100-3
(13) Bartolini,M.;Andrisano,V.ChemBioChem 2010,11,1018.doi: 10.1002/cbic.200900666
(14)Mason,J.M.Future Med.Chem.2010,2,1813.doi:10.4155/ fmc.10.259
(15) Salomone,S.;Caraci,F.;Leggio,G.M.;Fedotova,J.;Drago,F. Br.J.Clin.Pharmacol.2012,73,504.doi:10.1111/ j.1365-2125.2011.04134.x
(16) Sciarretta,K.L.;Gordon,D.J.;Meredith,S.C.Methods Enzymol.2006,413,273.doi:10.1016/S0076-6879(06)13015-3
(17) Funke,S.A.;Willbold,D.Curr.Pharm.Des.2012,18,755.doi: 10.2174/138161212799277752
(18) Li,H.;Zemel,R.;Lopes,D.H.;Monien,B.H.;Bitan,G. ChemMedChem 2012,7,515.doi:10.1002/cmdc.201100584
(19) Yang,C.;Zhu,X.;Li,J.;Shi,R.J.Mol.Model.2010,16,813. doi:10.1007/s00894-009-0594-y
(20) Veloso,A.J.;Kerman,K.Bioelectrochemistry 2012,84,49.doi: 10.1016/j.bioelechem.2011.08.007
(21) Jagota,S.;Rajadas,J.Int.J.Pept.Res.Ther.2012,18,53.doi: 10.1007/s10989-011-9278-4
(22) Tjernberg,L.O.;Naslund,J.;Lindqvist,F.;Johansson,J.; Karlstrom,A.R.;Thyberg,J.;Terenius,L.;Nordstedt,C. J.Biol.Chem.1996,271,8545.doi:10.1074/jbc.271.15.8545
(23) Soto,C.;Kindy,M.S.;Baumann,M.;Frangione,B.Biochem. Biophys.Res.Commun.1996,226,672.doi:10.1006/ bbrc.1996.1413
(24) Soto,C.;Sigurdsson,E.M.;Morelli,L.;Kumar,R.A.;Castano, E.M.;Frangione,B.Nat.Med.1998,4,822.doi:10.1038/ nm0798-822
(25) Adessi,C.;Frossard,M.J.;Boissard,C.;Fraga,S.;Bieler,S.; Ruckle,T.;Vilbois,F.;Robinson,S.M.;Mutter,M.;Banks,W. A.;Soto,C.J.Biol.Chem.2003,278,13905.doi:10.1074/jbc. M211976200
(26) Hetenyi,C.;Szabo,Z.;Klement,E.;Datki,Z.;Kortvelyesi,T.; Zarandi,M.;Penke,B.Biochem.Biophys.Res.Commun.2002, 292,931.doi:10.1006/bbrc.2002.6745
(27) Yang,C.;Li,J.Y.;Li,Y.;Zhu,X.L.J.Mol.Struct.-Theochem 2009,895,1.doi:10.1016/j.theochem.2008.10.003
(28)Wei,G.;Jewett,A.I.;Shea,J.E.Phys.Chem.Chem.Phys. 2010,12,3622.
(29) Liu,F.F.;Dong,X.Y.;He,L.;Middelberg,A.P.;Sun,Y. J.Phys.Chem.B 2011,115,11879.doi:10.1021/jp202640b
(30)Viet,M.H.;Ngo,S.T.;Lam,N.S.;Li,M.S.J.Phys.Chem.B 2011,115,7433.doi:10.1021/jp1116728
(31) Fradinger,E.A.;Monien,B.H.;Urbanc,B.;Lomakin,A.;Tan, M.;Li,H.;Spring,S.M.;Condron,M.M.;Cruz,L.;Xie,C.W.; Benedek,G.B.;Bitan,G.Proc.Natl.Acad.Sci.U.S.A.2008, 105,14175.doi:10.1073/pnas.0807163105
(32) Crescenzi,O.;Tomaselli,S.;Guerrini,R.;Salvadori,S.; DʹUrsi,A.M.;Temussi,P.A.;Picone,D.Eur.J.Biochem. 2002,269,5642.doi:10.1046/j.1432-1033.2002.03271.x
(33) van Der Spoel,D.;Lindahl,E.;Hess,B.;Groenhof,G.;Mark, A.E.;Berendsen,H.J.J.Comput.Chem.2005,26,1701.doi: 10.1002/jcc.20291
(34) van Gunsteren,W.F.;Billeter,S.R.;Eising,A.A.; Hünenberger,P.H.;Krüger,P.;Mark,A.E.;Scott,W.R.P.; Tironi,I.G.In Biomolecular Simulation:The GROMOS96 Manual and User Guide;Zürich,Switzerland,Groningen, Holland,1996.
(35) Berendsen,H.J.C.;Postma,J.P.M.;van Gunsteren,W.F.; Hermans,J.Intermolecular Forces;Pullmann,B.Ed.;Reidel: Dordecht,Holland,1981.
(36) Darden,T.;York,D.;Pedersen,L.J.Chem.Phys.1993,98, 10089.doi:10.1063/1.464397
(37) Verlet,L.Phys.Rev.1967,159,98.doi:10.1103/PhysRev.159.98
(38) Hess,B.;Berendsen,H.J.C.;Fraaije,J.G.E.M.J.Comput. Chem.1997,18,1463.doi:10.1002/(SICI)1096-987X(199709) 18:12<1463::AID-JCC4>3.0.CO;2-H
(39) Beredsen,H.J.C.;Postma,J.P.M.;van Gunsteren,W.F.;Di Nola,A.;Haak,J.R.J.Chem.Phys.1984,81,3684.doi: 10.1063/1.448118
(40) Kabsch,W.;Sander,C.Biopolymers 1983,22,2577.doi: 10.1002/bip.360221211
(41) Hu,J.P.;Sun,T.G.;Chen,W.Z.;Wang,C.X.Acta Chim.Sin. 2006,64,2079.[胡建平,孙庭广,陈慰祖,王存新.化学学报,2006,64,2079.]
(42) Kollman,P.A.;Massova,I.;Reyes,C.;Kuhn,B.;Huo,S.; Chong,L.;Lee,M.;Lee,T.;Duan,Y.;Wang,W.;Donini,O.; Cieplak,P.;Srinivasan,J.;Case,D.A.;Cheatham,T.E.,III. Accounts Chem.Res.2000,33,889.doi:10.1021/ar000033j
(43) Li,H.;Luo,Y.;Derreumaux,P.;Wei,G.J.Phys.Chem.B 2010, 114,1004.doi:10.1021/jp908889q
(44) Liu,F.F.;Dong,X.Y.;Sun,Y.Acta Phys.-Chim.Sin.2010,26, 1643.[刘夫锋,董晓燕,孙 彦.物理化学学报,2010,26, 1643.]
(45)Hu,J.P.;Gong,X.Q.;Su,J.G.;Chen,W.Z.;Wang,C.X. Biophys.Chem.2008,132,69.doi:10.1016/j.bpc.2007.09.008
(46)Huang,B.;Liu,F.F.;Dong,X.Y.;Sun,Y.J.Phys.Chem.B 2012,116,424.doi:10.1021/jp205770p
(47)Huang,B.;Liu,F.F.;Dong,X.Y.;Sun,Y.J.Phys.Chem.B 2011,115,4168.doi:10.1021/jp111216g
(48)Hou,T.;Wang,J.;Li,Y.;Wang,W.J.Comput.Chem.2011,32, 866.doi:10.1002/jcc.21666
(49) Liu,F.F.;Liu,Z.;Bai,S.;Dong,X.Y.;Sun,Y.J.Chem.Phys. 2012,136,145101.doi:10.1063/1.3702195
(50) Santini,S.;Wei,G.;Mousseau,N.;Derreumaux,P.Structure 2004,12,1245.doi:10.1016/j.str.2004.04.018
(51)Kokkoni,N.;Stott,K.;Amijee,H.;Mason,J.M.;Doig,A.J. Biochemistry 2006,45,9906.doi:10.1021/bi060837s
May 14,2012;Revised:July 16,2012;Published on Web:July 16,2012.
Molecular Dynamics Simulation and Binding Free Energy Calculation of the Conformational Transition of Amyloid Peptide 42 Inhibited by Peptide Inhibitors
DONG Xiao-Yan DU Wen-Jie LIU Fu-Feng*
(Department of Biochemical Engineering,School of Chemical Engineering and Technology,Tianjin University,Tianjin 300072,P.R. China; Key Laboratory of Systems Bioengineering,Ministry of Education,Tianjin University,Tianjin 300072,P.R.China)
The molecular mechanisms of the conformational transition of amyloid β-peptide(Aβ)42 inhibited by the peptide inhibitors KLVFF,VVIA,and LPFFD were studied by using molecular dynamics simulations and binding free energy calculations.These studies confirmed that the conformational transition of Aβ42 from its initial α-helix to β-sheet structure is prevented by these three peptide inhibitors. The calculations also demonstrated that the intra-peptide hydrophobic interactions of Aβ42 are weakened, and its quantity of long range contacts decreased by these inhibitors.Consequently,the hydrophobic collapse of Aβ42 is alleviated and its initial structure is maintained well.Both hydrophobic and electrostatic interactions,including hydrogen bonding,were found to favor the binding of these peptide inhibitors to Aβ42.Moreover,the charged residues of the inhibitors were shown to enhance the electrostatic interactions including hydrogen bonding,decreasing the capacity of the peptide for self-assembly,and increasing the inhibition effect.It was also determined that interactions between the inhibitors and Aβ42 are reduced when proline residue is introduced into the peptide inhibitor,since its linear structure is disrupted. In general,this work has allowed a better understanding of the molecular mechanisms of the effects of the peptide inhibitors KLVFF,VVIA,and LPFFD on the conformational transition of Aβ42 and will assist in the systematic design of high efficiency peptide inhibitors of Aβ aggregation.
Molecular dynamics simulation;Alzheimerʹs disease;Amyloid peptide;Peptide inhibitor;Conformational transition
10.3866/PKU.WHXB201207162
∗Corresponding author.Email:fufengliu@tju.edu.cn;Tel:+86-22-27406590.
The project was supported by the National Natural Science Foundation of China(20906068,21076149),National Key Basic Research Program of China(973)(2009CB724705),and Natural Science Foundation of Tianjin from Tianjin Municipal Science and Technology Commission,China (10JCYBJC04500).
国家自然科学基金(20906068,21076149),国家重点基础研究发展规划项目(973)(2009CB724705)和天津市科委自然科学基金(10JCYBJC04500)资助
O641