時術(shù)華，張少龍，張慶剛

（1.山東建筑大學理學院，濟南250101；2.山東師范大學物理與電子科學學院，濟南250014）

1 引 言

近期的醫(yī)學臨床研究表明全世界眾多癌癥患者的體內(nèi)基本上都發(fā)現(xiàn)抑癌蛋白p53野生型功能的丟失.p53抑癌功能的丟失有兩種公認的機制：（1）p53 蛋白內(nèi)部的定點突變；（2）腫瘤蛋白MDM2通過副反饋的方式與p53相互作用，抑制了p53的功能.研究發(fā)現(xiàn)50%的癌癥患者體內(nèi)存在MDM2 的 過 度 表 達［1，2］.因 此 干 擾p53 與MDM2的結(jié)合成為癌癥治療的新途徑.

MDM2與p53復合物的結(jié)晶結(jié)構(gòu)表明p53的三個疏水性殘基Phe19，Trp23 和Leu26 插入MDM2疏水性裂縫，從而與MDM2 形成疏水性相互作用.用于治療癌癥的藥物基本上都模仿了這個作用模式.幾個課題組設(shè)計的肽類抑制劑能與MDM2產(chǎn)生較強的相互作用，擁有納摩爾數(shù)量級的抑制效果［3，4］.Phan等人設(shè)計合成的肽類抑制劑PDIQ 展現(xiàn)了不錯的抑制效果，其IC50值達到8nM［5］.該結(jié)果表明抑制劑PDIQ 對MDM2有很強的抑制效果.故闡明該抑制劑與MDM2的作用機制有助于新型高效抑制劑的研發(fā).

分子動力學模擬和結(jié)合自由能計算已經(jīng)成為研究抑制劑與蛋白質(zhì)相互作用的重要工具［6-8］.盡管Phan等人從結(jié)構(gòu)上解釋了抑制劑PDIQ 的抑制能力，但未從原子層次上闡明抑制劑PDIQ 與MDM2的作用機制.因此本工作將蛋白質(zhì)庫的晶格結(jié)構(gòu) （3JZS）用作動力學模擬的初始結(jié)構(gòu)，且采用已成功用于研究生物體系的MM-PBSA 方法計算抑制劑PDIQ 與MDM2的結(jié)合自由能，解釋二者的結(jié)合模式［9-12］，從原子層次上研究了PDIQ與MDM2 的作用機制和MDM2 內(nèi)部殘基運動.同時采用基于殘基的自由能分解方法計算抑制劑-殘基相互作用以便闡明PDIQ-MDM2復合物的結(jié)構(gòu)-親和能關(guān)系.我們期望該研究能為以p53-MDM2相互作用為靶標的抗癌藥物設(shè)計提供理論上的指導.

2 模型與計算

2.1 分子動力學模擬

用于動力學模擬的初始構(gòu)象取自蛋白質(zhì)庫的晶體結(jié)構(gòu) （3JZS）.圖1 給出PDIQ-MDM2 復合體的結(jié)構(gòu).復合物X-射線結(jié)構(gòu)中的水分子保留在初始化模型中，采用Amber 12 中的leap模塊添加晶體結(jié)構(gòu)中缺失的氫原子［13］；復合物的力場參數(shù)由Amber 12 中 的ff03 力 場 產(chǎn) 生［14］；PDIQMDM2復合體溶解在由5535個水分子組成的顯性水合子中，水分子的模為TIP3P；兩個氯離子添加到由水和復合體組成的系統(tǒng)中以確保整個系統(tǒng)維持電中性.

為消除系統(tǒng)中原子的不合理接觸，采用Amber12中的sander模塊對復合體體系執(zhí)行兩步的系統(tǒng)優(yōu)化：（1）采用150kcal／ （mol·?2）的力常數(shù)約束溶質(zhì)，以便優(yōu)化溶劑和中和離子；（2）無約束地優(yōu)化整個系統(tǒng)：每一步優(yōu)化均先執(zhí)行2000步的最陡下降優(yōu)化，接著執(zhí)行4000步的共軛梯度優(yōu)化.然后在200ps內(nèi)把系統(tǒng)從0K 加熱到300 K，隨后進行100ps的常溫300K，常壓1標準大氣壓的動力學平衡；最后是10ns的無約束分子動力學模擬.模擬期間采用SHAKE方法限制所有含氫原子化學鍵的伸縮［15］，模擬積分步長為2fs，PME方法用來計算長程的靜電相互作用，應(yīng)用周期性邊界條件以消除溶劑盒子的邊緣效應(yīng)，非成鍵相互作用的截斷值為10.0?.

圖1 PDIQ-MDM2 復合體的結(jié)構(gòu).MDM2 和PDIQ 分別用桔色和青色表示，殘基Phe3′，Trp7′和Leu26′用棍棒方式顯示Fig.1 Structure of PDIQ-MDM2 complex.MDM2and PDIQ are showed in orange and green，respectively，the residues Phe3′，Trp7′and Leu10′are displayed in stick mode

2.2 MM-PBSA方法計算結(jié)合自由能

采用MM-PBSA 方法計算PDIQ 與腫瘤蛋白MDM2的結(jié)合自由能.從動力學最后2ns的軌跡中每隔20ps抽取的100個構(gòu)象用于結(jié)合自由能計算，結(jié)合自由能計算的公式可表示為：

式中ΔEMM是氣相中的分子力學能，ΔGsol表示溶解自由能對分子結(jié)合的貢獻，TΔS表示熵變對結(jié)合自由能的貢獻.ΔEMM進一步分成兩部分：

方程中的ΔEele和ΔEvdw分別表示靜電相互作用和范德瓦爾斯相互作用.溶解自由能ΔGsol也分解為兩個組成部分：

上式中ΔGpb和ΔGsurf分別是極性溶解自由能和非極性溶解自由能.前者可使用Amber中的pbsa方法求解泊松-玻爾茲曼方程獲得，溶質(zhì)和溶劑的電介常數(shù)分別設(shè)為1.0和80.0.后者由如下經(jīng)驗方程計算：

式中γ和β值分別取為0.00542kcal／ （mol·?2）和0.92kcal／mol.TΔS是由于平動、轉(zhuǎn)動和振動自由度變化引起的熵變對結(jié)合自由能的貢獻，主要采用簡振模和傳統(tǒng)的熱力學計算TΔS.

表1 MM-PBSA 計算所得到的能量 （kcal／mol）aTable 1 Energies（in kcal／mol）calculated by MM-PBSA method a

2.3 相關(guān)矩陣計算

為了研究抑制劑結(jié)合誘導的蛋白質(zhì)內(nèi)部運動的相關(guān)聯(lián)程度，計算了相關(guān)矩陣 （Cij）.這一矩陣第i行j 列的元素Cij反映了第i，j殘基上Cα相對于其平均位置的漲落相關(guān)系數(shù)，通常由下列的方程確定：

式中尖括號表示對取樣的系綜平均，Δri表示第i個原子相對于其平均位置的偏差，Cij的變化范圍為-1到+1，正的Cij反映了殘基間的相關(guān)運動，而負的Cij表述了殘基間的反相關(guān)運動.

3 結(jié)果與討論

3.1 動力學平衡的穩(wěn)定性和結(jié)構(gòu)的柔性

為評估分子動力學平衡的穩(wěn)定性，采用Amber中的Ptraj程序計算了MDM2 主鏈原子相對于初始優(yōu)化結(jié)構(gòu)的均方根偏差 （RMSD）隨時間的函數(shù) （見圖2）.圖2表明動力學模擬4ns后，系統(tǒng)達到平衡.RMSD 的平均值為1.31?，漲落范圍低于0.50?.這表明動力學模擬平衡穩(wěn)定性是可靠的.

3.2 相關(guān)運動分析

為了研究MDM2內(nèi)部殘基的相關(guān)運動，我們基于分子動力學模擬軌跡中原子位置，構(gòu)建描述原子位置漲落的相關(guān)矩陣 （見圖3）.依據(jù)圖3，整體而言抑制劑PDIQ 的結(jié)合并沒有誘導MDM2構(gòu)象的明顯變化.

圖2 MD 模擬中主鏈原子的均方根偏差Fig.2 Root-Mean-Square Deviation of the backbone atoms on PDIQ-MDM2 complex

從圖3 可以觀察到，子單元S1 和S2 內(nèi)部（對角方塊）發(fā)生強烈的相關(guān)運動，同時子單元S1內(nèi)部也有明顯的反相關(guān)運動.子單元S3 相對于S4以相關(guān)運動的趨勢運動.圖3表明子單元S2相對與S1和S3產(chǎn)生強烈的反相關(guān)運動，說明此處可能發(fā)生較大的構(gòu)象變化.此外，子單元S3也與S1發(fā)生較弱的反相關(guān)運動.

圖3 系統(tǒng)平衡后主鏈Cα原子相對于其平均位置漲落的相關(guān)矩陣Fig.3 Cross-correlation matrices of the fluctuations of the coordination for Cαatoms around their mean positions after the equilibrium of system

3.3 結(jié)合自由能計算

采用MM-PBSA 方法計算了抑制劑PDIQ 與MDM2的結(jié)合自由能，其各項作用成分的貢獻列在表1中.依據(jù)表1，PDIQ 與MDM2的結(jié)合自由能為-19.8kcal／mol，這表明PDIQ 與MDM2 產(chǎn)生較強的相互作用.

表1表明范德瓦爾斯作用能 （ΔEvdw）和非極性的溶解自由能 （ΔGsurf）分別是-67.76 和-7.80 kcal／mol.盡管氣相中的靜電相互作用 （ΔEele）也提供了-135.2kcal／mol的有利貢獻，但是完全被不利于抑制劑結(jié)合的153.20kcal／mol的極性溶解自由能所抵消，產(chǎn)生了不利于抑制劑結(jié)合的極性作用 （17.99kcal／mol）.由表1看到熵變對結(jié)合自由能的貢獻 （-TΔS=36.78）也不利于抑制劑結(jié)合，由熵變的計算得到有利的信息，抑制劑結(jié)合前后熵明顯減少，這說明抑制劑結(jié)合后結(jié)合位點原子的運動受到了很大約束.總之有利于抑制劑結(jié)合的成分中范德瓦爾斯作用占了89.7%，因此范德瓦爾斯作用驅(qū)動了PDIQ 與MDM2的結(jié)合.這個結(jié)果基本上吻合了幾個課題組先前的研究［16-19］.

3.4 結(jié)構(gòu)-親和能關(guān)系

為了定量地從原子層次上研究抑制劑PDIQ與MDM2的結(jié)合模式，采用基于殘基的自由能分解方法計算了抑制劑-殘基相互作用，相互作用譜如圖4所示.依據(jù)圖4，PDIQ 能與MDM2的8個殘 基 Leu54，Ile61， Met62， Tyr67， Gln72，His73，Val93和His96產(chǎn)生強于1.4kcal／mol的相互作用，其結(jié)構(gòu)親和能關(guān)系分析如下.

圖4 抑制劑PDIQ 中關(guān)鍵殘基與MDM2各個分離殘基的相互作用譜Fig.4 Residue-residue interaction spectrum between the inhibitorp DIQ and MDM2

圖4表明8個殘基中Leu54與MDM2產(chǎn)生了最強的相互作用，其相互作用能為-3.81kcal／mol.從結(jié)構(gòu)圖5來看，這個較強的相互作用主要源自兩個因素的貢獻：（1）Trp7′疏水性環(huán)上的氮原子與殘基Leu54 的羰基氧原子形成氫鍵 （表2）；（2）Trp7′的疏水性環(huán)與Leu54 的烷基形成CH-π 相互作用.抑制劑PDIQ 與Val93的相互作用也非常強，作用能是-3.05kcal／mol.從結(jié)構(gòu)上看，這個作用主要由Trp7′吲哚環(huán)與Val93的烷基間CH-π 相互作用和Leu10′的烷基與Val93烷基間的CH-CH 相互作用提供.依據(jù)圖5，Met62和Ile61的烷基與殘基Phe3′的苯環(huán)臨近形成的CHπ相互作用分別為PDIQ 與MDM2的作用分別提供了-2.57和-2.25kcal／mol的能量貢獻 （圖4）.基于軌跡的氫鍵動力學分析表明，殘基Gln72與殘基Phe3′形成一個氫鍵 （表2），這個氫鍵為PDIQ 與MDM2 的 相 互 作 用 提 供 了-1.99kcal／mol的貢獻.殘基His73與抑制劑PDIQ 的作用是-1.71kcal／mol，這與His73的5元環(huán)和Phe3′的苯環(huán)間的π-π相互作用吻合.殘基Tyr67的苯環(huán)在結(jié)構(gòu)上與Phe3′的苯環(huán)臨近且平行，故能形成平行的π-π相互作用 （圖5），該作用提供了大約-1.47kcal／mol的能量貢獻.殘基His96與PDIQ的作用是-1.5kcal／mol，該作用主要來自His96的疏水性5元環(huán)與Leu10′烷基間的CH-π 相互作用.上述分析結(jié)果基本上吻合幾個課題組先前的工作［17-21］.

表2 主要殘基的氫鍵分析Table 2 The hydrogen bonds of the key residues

圖5 pDIQ-MDM2復合體中主要殘基的相對位置以及形成的氫鍵Fig.5 Relative geometries of the key residues in PDIQ-MDM2complex and hydrogen bonds

綜合以上的結(jié)構(gòu)親和能關(guān)系可以發(fā)現(xiàn)，殘基Ile61，Met62，Tyr67，Gln72 和His73 形成抑制劑殘基Phe3′的疏水性結(jié)合口袋；殘基Leu54，His73和Val93形成抑制劑殘基Tyr7′的疏水性結(jié)合口袋；PDIQ 殘基Leu10′的疏水性結(jié)合裂縫由殘基Leu54，Val93和His96構(gòu)成.在這些結(jié)合口袋中主要的結(jié)合方式為CH-CH，CH-π和π-π相互作用.因此這些作用驅(qū)動了抑制劑與MDM2的結(jié)合.

4 結(jié) 語

基于動力學軌跡的相關(guān)矩陣分析表明抑制劑結(jié)合誘導了反相關(guān)運動的主體形式.結(jié)合自由能的計算結(jié)果表明范德瓦爾斯作用主導了抑制劑與MDM2的結(jié)合.基于殘基的自由能分解揭示了MDM2疏水性結(jié)合口袋內(nèi)的主要作用形式，CHCH，CH-π 和π-π 相 互 作 用 驅(qū) 動 了 抑 制 劑 與MDM2的結(jié)合.這個研究結(jié)果能為抗癌藥物的設(shè)計提供一定的理論啟示.

［1］ Chen H F，Luo R.Binding induced folding in p53-MDM2complex［J］.J.Am.Chem.Soc.，2007，129：2930.

［2］ Ding K，Lu Y，Nikolovska-Coleska Z，et al.Structure-based design of spiro-oxindoles as potent，specific small-molecule inhibitors of the MDM2-p53interaction［J］.J.Med.Chem.，2006，49：3432.

［3］ Czarna A，Popowicz G M，Pecak A，et al.High affinity interaction of the p53 peptide-analogue with human Mdm2and Mdmx［J］.Cell Cycle，2009，8：1176.

［4］ Pazgier M，Liu M，Zou G，et al.Structural basis for high-affinity peptide inhibition of p53interactions with MDM2and MDMX ［J］.Proc.Natl.Acad.Sci.USA，2009，106：4665.

［5］ Phan J，Li Z，Kasprzak A，et al.Structure-based design of high affinity peptides inhibiting the interaction of p53with MDM2and MDMX［J］.J.Biol.Chem.，2010，285：2174.

［6］ Chen J，Liang Z，Zhang Q，et al.Insight into Interaction Mechanism of the inhibitor pDI6W with MDM2Based on Molecular Dynamics［J］.J.At.

Mol.Phys.，2012，29：953（in Chinese）［陳建中，梁志強，張慶剛等.分子動力學研究抑制劑pDI6W與MDM2的相互作用機制［J］.原子與分子物理學報，2012，29：953］

［7］ Chen J，Zhang D，Zhang Y，et al.Computational studies of difference in binding modes of peptide and non-peptide inhibitors to MDM2／MDMX based on molecular dynamics simulations［J］.Int.J.Mol.Sci.，2012，13：2176.

［8］ Chen J，Zhang S，Liu X，et al.Insights into drug resistance of mutations D30N and I50V to HIV-1 protease inhibitor TMC-114：Free energy calculation and molecular dynamic simulation ［J］.J.Mol.Model.，2010，16：459.

［9］ Wu E L，Han K L，Zhang J Z H.Selectivity of neutral／weakly basic P1group inhibitors of thrombin and trypsin by a molecular dynamics study［J］.Chemistry-A European Journal，2008，14：8704.

［10］ Case D A，Cheatham III T E，Darden T，et al.The Amber biomolecular simulation programs［J］.J.Comput.Chem.，2005，26：1668.

［11］ Case D A，Cheatham T E，Darden T，et al.The Amber biomolecular simulation programs［J］.J.Comp.Chem.，2005，26：1668.

［12］ Chen J，Yang M，Hu G，et al.Insights into the functional role of protonation states in the HIV-1 protease-BEA369complex：molecular dynamics sim-ulations and free energy calculations［J］.J.Mol.Model.，2009，15：1245.

［13］ Case D A，Darden T A，Cheatham III T E，et al.AMBER 12.San Francisco：University of California，2012.

［14］ Duan Y，Wu C，Chowdhury S，et al.A pointcharge force field for molecular mechanics simulations of proteins based on condensed-phase quantum mechanical calculations ［J］.J.Comp.Chem.，2003，24：1999.

［15］ Coleman T G，Mesick H C，Darby R L.Numerical integration［J］.Ann.Biomed.Eng.，1977，5：322.

［16］ Cheng W，Chen J，Liang Z，et al.A computational analysis of interaction mechanisms of peptide and non-peptide inhibitors with MDMX based on molecular dynamics simulation ［J］.Computational and Theoretical Chemistry，2012，984：43.

［17］ Cheng W Y，Liang Z Q，Zhang Q G，et al.Insight into p53-MDM2interaction based on molecular dynamics simulation and molecular mechanics［J］.J.At.Mol.Phys.，2012，29：393（in Chinese）［程偉淵，梁志強，張慶剛，等.p53-MDM2 相互作用的分子力學和動力學研究［J］.原子分子物理學報，2012，29：393］

［18］ Chen J，Wang J，Xu B，et al.Insight into mechanism of small molecule inhibitors of the MDM2-p53 interaction：molecular dynamics simulation and free energy analysis ［J］.J.Mol.Graph.Model.，2011，30：46.

［19］ Ding Y，Mei Y，Zhang J Z H.Quantum mechanical studies of residue-specific hydrophobic interactions in p53-MDM2 binding ［J］.J.Phys.Chem.B，2008，112：11396.

［20］ Hu G，Wang D，Liu X，et al.A computational analysis of the binding model of MDM2 with inhibitors［J］.J.Comp.Aid.Mol.Des.，2010，24：687.

［21］ Lu S Y，Jiang Y J，Zou J W，et al.Molecular modeling and molecular dynamics simulation studies on pyrrolopyrimidine-based ［alpha］-helix mimetic as dual inhibitors of MDM2and MDMX［J］.J.Mol.Graph.Model.，2011，30：167.