HTLV—1蛋白酶與印地那韋識別的分子動力學(xué)模擬

常珊等

摘要：對HTLV-1蛋白酶和抑制劑印地那韋的復(fù)合物進行了11 ns的分子動力學(xué)模擬，從分子整體、氫鍵、結(jié)合自由能及運動性等方面進行了分析。結(jié)果表明，動力學(xué)模擬使得體系更為松弛，但大部分氫鍵仍得到保持；疏水相互作用有利于HTLV-1蛋白酶和印地那韋的結(jié)合，關(guān)鍵殘基主要分布在R10、L30～V39、V56～F67和W98～I1004個區(qū)域；結(jié)合自由能預(yù)測值與試驗數(shù)據(jù)吻合較好。印地那韋結(jié)合僅能部分抑制HTLV-1蛋白酶的Flap、top-site以及兩個exo-site位點柔性，解釋了印地那韋的HIV-1 PR抑制活性為何略大于HTLV-1 PR。該研究為基于受體結(jié)構(gòu)的抗HTLV-1藥物研發(fā)打下了一定的結(jié)構(gòu)和理論基礎(chǔ)。

關(guān)鍵詞：HTLV-1蛋白酶；印地那韋；分子動力學(xué)模擬；分子識別；結(jié)合自由能

中圖分類號：O641.3；R969.1 文獻標(biāo)識碼：A 文章編號：0439-8114（2015）12-3034-07

DOI：10.14088/j.cnki.issn0439-8114.2015.12.062

Molecular Dynamics Simulation of the Recognition between HTLV-1

Protease and Indinavir

CHANG Shan1，LIANG Li2，LIU Wei2，HU Jian-ping2，3，GOU Xiao-jun2

（1. Institute of Bioinformatics and Medical Engineering， Jiangsu University of Technology， Changzhou 213001， Jiangsu China；

2. Key Laboratory of Medicinal and Edible Plants Resources Development， Chengdu University， Chengdu 610106， China；

3. College of Chemistry， Leshan Normal University， Leshan 614004， Sichuan China）

Abstract： One 11 ns molecular dynamics simulation （MD） for the compocind system of Human T-cell leukemia virus type 1（HTLV-1） protease （PR） and its inhibitor indinavir was performed. The mutual recognition of the two molecules was analyzed from exploring the whole molecule，hydrogen bond，binding free energy and mobility，etc.After dynamics simulation，the system keeps more relaxed compared with the crystal structure， and most of the hydrogen bonds were still maintained. Both energy decomposition and free energy calculation showed that，hydrophobic interaction facilitates the association of HTLV-1 PR and indinavir，and key residues were mainly distributed in the four regions including R10，L30～V39，V56～F67 and W98～I100.Additionally，binding free energy prediction was in good agreement with the experimental data. Indinavir can only partially attenuate the flexibility of the flap，top-site and two exo-site in the PR system，which investigated why indinavir possess more HIV-1 PR inhibition ability than that of HTLV-1 PR. The study will provide some structural and theoretical basis for anti-HTLV-1 drug research based on the receptor structure.

Key words： HTLV-1 protease； Indinavir； molecular dynamics simulation； molecular recognition； binding free energy

I型人類T細胞白血病病毒（Human T-cell leukemia virus type 1，HTLV-1）參與了眾多種類癌癥，如成人T細胞白血病、熱帶痙攣性脊髓病以及HTLV-1相關(guān)的皮炎等[1，2]。目前，世界上約2 000萬人感染了該病毒，臨床治療策略主要是采用細胞毒性化療及干擾素、齊多呋定用藥，開展抗HTLV-1的新型藥物研發(fā)比較迫切[3]。

HTLV-1蛋白酶（Protease，PR）是同源蛋白二聚體，每個鏈由125個氨基酸殘基，屬于天冬氨酸蛋白酶家族[4，5]。HIV-1蛋白酶被認定為抗AIDS藥物研發(fā)的一線靶點已有30年歷史，目前FDA已經(jīng)批準(zhǔn)10多種HIV-1蛋白酶抑制劑藥物上市[6]。而HTLV-1 PR和HIV-1 PR的序列同源性為28%，其中抑制劑結(jié)合區(qū)同源性更高，達到了45%。兩類PR盡管具有相似的整體折疊，但是有不同的抑制劑結(jié)合特異性[7]。大部分HIV-1 PR抑制劑沒有體現(xiàn)出抑制HTLV-1 PR的潛力，目前僅有HIV-1蛋白酶抑制劑印地那韋能夠在較低的微摩爾級層次上抑制HTLV-1 PR活性，抑制常數(shù)為3.5 μm[8]。

圖1給出了印地那韋的分子結(jié)構(gòu)，分子式為C36H47N5O4，由主鏈（C10-C11-C12-C13-C21-N4）和側(cè)鏈部分（比如吡啶基、哌嗪基、叔丁基，苯基，茚基）組成。綜合考慮HTLV-1 PR和HIV-1 PR的結(jié)構(gòu)類似性和印地那韋的雙重抑制活性，本研究利用分子動力模擬方法研究了印地那韋和HTLV-1 PR的相互作用，將有助于后續(xù)的印地那韋類抗HTLV-1 PR抑制劑的設(shè)計。

1 模擬方法

1.1 分子動力學(xué)模擬

基于HTLV-1蛋白酶結(jié)合抑制劑印地那韋的復(fù)合物晶體[9]結(jié)構(gòu)，僅保留A鏈和B鏈接的氨基酸部分以及底物印地那韋，以此作為分子動力學(xué)（Molecular dynamics，MD）模擬的初始構(gòu)象。MD模擬應(yīng)用AMBER程序[10]，采用的力場是AMBER力場，蛋白質(zhì)的力場參數(shù)均基于試驗值擬合[11]，模擬溫度為300 K，溶劑采用TIP3P水模型[12]。在溶質(zhì)印地那韋外圍加上1.0 nm的去頭八面體水盒子，總共加入9 913個水分子，含水后體系的總原子數(shù)為33 453。在MD模擬之前，對體系進行兩次能量優(yōu)化：①約束溶質(zhì)[約束力常數(shù)為2.09×103 kJ/（mol·nm2）]，用最陡下降法優(yōu)化5 000步，再用共軛梯度法優(yōu)化5 000步；②去約束后，再進行5 000 步最陡下降法優(yōu)化和5 000步共軛梯度法優(yōu)化，收斂條件為能量梯度小于4.18×10-4 kJ/（mol·nm）。

MD 模擬分為兩步：首先進行1 ns約束溶質(zhì)的MD 模擬[（約束力常數(shù)為41.8 kJ/（mol·nm2）]，溫度從0 K逐步升高到300 K；再進行10 ns的無約束恒溫MD模擬，用VMD 圖像顯示軟件實時跟蹤體系的構(gòu)象。采用SHAKE 算法[13]約束鍵長，MD模擬的積分步長設(shè)為2 fs，非鍵相互作用截斷半徑設(shè)為1.2 nm，每隔1 ps采集1次構(gòu)象，共采集11 000個構(gòu)象。

1.2 MM-PBSA自由能計算

用MM-PBSA（Molecular mechanics/Poisson-Boltzmann surface area）方法[14，15]計算了HTLV-1 PR和底物印地那韋的絕對結(jié)合自由能，該方法通過從復(fù)合物的MD模擬軌跡中每隔一定的時間間隔提取出體系結(jié)構(gòu)，計算出平均結(jié)合自由能，計算公式為：

ΔGbind=ΔH-TΔS≈ΔEMM+ΔGsol-TΔS

式中，ΔGbind是結(jié)合自由能，ΔH是焓變，TΔS是絕對溫度和熵變的乘積；ΔEMM是真空下分子內(nèi)能差，包括分子非極性能差（VDWIN）和分子靜電能（ELEIN）差；ΔGsol是溶劑化自由能，ΔEMM可通過分子力學(xué)方法計算得到，ΔGsol由極性溶劑化自由能差（ELEPB）和非極性溶劑化自由能（VDWPB）差組成，極性部分通過求解有限差分Poisson-Blotzmann（PB）方程求得[16-18]，非極性部分通過估算溶劑可接近表面積（Surface area，SA）[19]擬合得到。TΔS的計算采用了正則模方法[20]。

1.3 能量分解

能量分解采用MM/GBSA（Molecular mechanics/generalized Born surface area）方法[14，15]，將每個殘基的能量貢獻近似分為用分子力學(xué)（MM）方法計算的真空分子內(nèi)能、用廣義波恩（GB）模型[21，22]計算的極性溶劑化能，以及用LCPO模型[23，24]計算的非極性溶劑化能，并且把能量分解到殘基的主鏈原子和側(cè)鏈原子上。其中，真空分子內(nèi)能是由極性部分（靜電相互作用）和非極性部分（范德華相互作用）組成的。LCPO模型的核心是非極性溶劑化能與分子溶劑可接近表面積正相關(guān)，通過MM/GBSA能量分解，可以考察蛋白質(zhì)中的主要殘基對于結(jié)合底物分子所做的貢獻。本研究基于HTLV-1蛋白酶與抑制劑印地那韋的復(fù)合物體系，MD 模擬平衡后6.5～11.0 ns 內(nèi)的軌跡（共4.5 ns），每隔100 ps采樣一次，共取45個構(gòu)象，通過計算獲得所有構(gòu)象的平均能量分解結(jié)果。

1.4 構(gòu)象成簇

對體系MD 模擬所獲得的11 000 個分子構(gòu)象進行了成簇分析，分析采用MMTSB工具[25]，計算過程依據(jù)Daura等[26]報道的策略。逐個計算各個構(gòu)象間的RMSD值，隨后建立RMSD矩陣（N×N，N為構(gòu)象數(shù)）。設(shè)定一個RMSD域值，如體系兩個構(gòu)象之間的RMSD小于該域值，則將它們歸為一簇，最后將每簇中能量最低構(gòu)象作為該簇的代表性結(jié)構(gòu)。本研究的成簇分析則是建立HTLV-1 PR-印地那韋復(fù)合物體系整體RMSD矩陣，即計算各個體系間的整體RMSD值，以分析基于整體構(gòu)象的成簇分布，進而依據(jù)構(gòu)象所占概率獲得分子整體所體現(xiàn)出的優(yōu)勢構(gòu)象。

2 結(jié)果與分析

2.1 體系的整體動力學(xué)

2.1.1 勢能及回轉(zhuǎn)半徑 圖2給出了HTLV-1蛋白酶與印地那韋復(fù)合物體系的勢能和回轉(zhuǎn)半徑（Radius of gyration）隨時間的變化。從圖2A可以看出，在約束溶質(zhì)的MD模擬初級階段（即0～1 ns之間），勢能平均值和標(biāo)準(zhǔn)偏差分別是-405 811×105 kJ/mol和6 477×105 kJ/mol，偏差率為1.6%。在去約束的MD模擬階段（即1～11 ns之間），勢能維持在-（414 201±681）×105 kJ/mol，偏差率為0.2%。在長時間全含水MD模擬中，體系勢能維持穩(wěn)定，模擬過程及軌跡可靠。從圖2B可以看出，HTLV-1蛋白酶結(jié)構(gòu)最初保持約束狀態(tài)（0～1 ns之間），回轉(zhuǎn)半徑數(shù)值穩(wěn)定，均值為1.895 nm；在去約束后的1.0～5.5 ns之間，回轉(zhuǎn)半徑漲落幅度變大，均值為1.899 nm；經(jīng)過5.5～6.5.0ns的回轉(zhuǎn)半徑（均值1.913 nm）快速增加階段，蛋白質(zhì)已經(jīng)去掉了晶體中固有的大部分約束和碰撞，結(jié)構(gòu)更松弛。在6.5～11 ns之間，回轉(zhuǎn)半徑維持在1.930 nm附近漲落。

2.1.2 RMSD分析 圖3給出了復(fù)合物RMSD隨時間的變化。從圖3可以看出，HTLV-1 PR在0～1 ns的約束MD模擬過程中，RMSD非常穩(wěn)定，均值維持在0.01 nm；在1.0～5.5 ns區(qū)間對應(yīng)的數(shù)值為（0.13±0.01） nm，5.5～6.5 ns為RMSD快速增加區(qū)間，在6.5～11.0 ns區(qū)間，對應(yīng)的數(shù)值為（0.19±0.01） nm。A鏈和B鏈的RMSD走勢稍有差別，A鏈在1.0～5.5 ns區(qū)間的RMSD要低于B鏈，5.5 ns之后，二者走勢趨同。就印地那韋而言，從模擬開始到3.0 ns之間，RMSD迅速提高到0.3 nm，隨后逐漸保持平穩(wěn)，數(shù)值為（0.22±0.02） nm，略大于蛋白質(zhì)部分，可見小分子有較高的運動性。與圖2B所示的回轉(zhuǎn)半徑結(jié)果相比，蛋白質(zhì)部分RMSD和回轉(zhuǎn)半徑數(shù)值快速增加的時間節(jié)點相同，均為5.5～6.5 ns之間，蛋白質(zhì)RMSD走勢與回轉(zhuǎn)半徑有較強相關(guān)性。另外，蛋白質(zhì)A鏈與B鏈的運動性略有不同，相對而言A鏈略高于B鏈。

2.1.3 RMSF分析 方均根漲落（Root mean square deviation，RMSF）是體系原子相對于平均構(gòu)型的方均根漲落，反映了體系在整個模擬中原子的平均變化情況。RMSF數(shù)值大說明對應(yīng)位置的柔性大，反之相反。圖4給出了HTLV-1蛋白酶二聚體的A鏈和B鏈的Cα原子的RMSF及柔性分布。從圖4A可以看出，HTLV-1蛋白酶的兩條鏈RMSF分布類似，這里將RMSF值大于0.4 nm和小于0.2 nm的區(qū)域分別定義為柔性區(qū)和剛性區(qū)?？梢姡嵝詤^(qū)有S46～N53，G60～F67和F80～T83；剛性區(qū)有A29～V39，A99～A105。圖4B直觀給出了體系柔性的區(qū)域分布，柔性大的區(qū)域都是遠離結(jié)合部位的無規(guī)卷曲或β片；剛性區(qū)則是直接參與識別抑制劑IND，也包括附近剛性較大的α螺旋區(qū)。另外，A鏈的柔性略高于B鏈。

在X-ray試驗中，蛋白質(zhì)的B因子（B-factor）體現(xiàn)了晶體中原子電子密度的模糊度（Diffusion），模糊度實際反映了蛋白質(zhì)分子在晶體中的構(gòu)象狀態(tài)。B因子越高，模糊度越大，相應(yīng)部位的構(gòu)象柔性越大。將基于X-ray試驗獲得的B因子和MD模擬獲得的RMSF進行了相關(guān)性分析。HTLV-1蛋白酶A鏈和B鏈的Cα原子B因子及RMSF數(shù)據(jù)的相關(guān)性，相關(guān)性值分別為0.66和0.91（圖5A、圖5B）。在116個樣本數(shù)前提下，相關(guān)性非常高，可見蛋白質(zhì)A鏈和B鏈的柔性分布類似，也表明用X-ray和MD方法來表述蛋白質(zhì)柔性分布的有效性。蛋白質(zhì)A鏈及B鏈的B因子和RMSF數(shù)值的相關(guān)性，相關(guān)系數(shù)分別為0.56和0.71（圖5C、圖5D）。綜上，基于RMSF分析得到的蛋白質(zhì)柔性分布結(jié)果可信，模擬所得軌跡可以用于后續(xù)分子間識別的探討。

2.2 分子識別

2.2.1 氫鍵分析 根據(jù)圖3的RMSD分布，分析了7.0～9.0 ns之間的氫鍵分布。氫鍵采用幾何判據(jù)[27]，供體和受體之間距離小于0.35 nm，供體-氫-受體的夾角大于135°，氫鍵占有率表示在整個模擬過程中氫鍵構(gòu)象出現(xiàn)的幾率。從表1可以看出，參與抑制劑印地那韋和蛋白酶識別的有6個氫鍵，氫鍵給體涉及到抑制劑的O2、O4、N1和N4原子，貫穿小分子的大部分結(jié)構(gòu)。氫鍵是促進印地那韋與HTLV-1蛋白酶特異性識別的重要作用力。X-ray試驗結(jié)果[28]表明，R10、D32、G34、D36、A59、W98能與印地那韋形成氫鍵，在全含水MD優(yōu)化后，大部分氫鍵得到了保持。

圖6給出了HTLV-1蛋白酶B鏈D32-OD2與印地那韋-O2-H18間最穩(wěn)定氫鍵的距離及角度隨時間的變化。在7～9 ns間，受體-供體距離和受體-氫-供體角度維持穩(wěn)定，漲落幅度較小，在2 000個構(gòu)象中，共有1 889個構(gòu)象均保持有印地那韋-O2-H18-D32-OD2氫鍵（占有率為94.45%）。

2.2.2 結(jié)合自由能計算 基于HTLV-1蛋白酶-印地那韋復(fù)合物7～9 ns的MD模擬軌跡，刪除水分子和4個平衡氯離子，然后從該軌跡中每隔200 ps取一次構(gòu)象，共11個構(gòu)象用于自由能計算，最終結(jié)果通過取平均得到。表2列出了各能量項對結(jié)合自由能的貢獻，用MM-PBSA方法計算得到的結(jié)合自由能模擬值-32.87 kJ/mol。印地那韋抑制HTLV-1蛋白酶的抑制常數(shù)Ki為3.5 μmol/L，相應(yīng)的ΔGbind試驗值為-31.33 kJ/mol，自由能計算值和試驗值非常吻合。

在形成蛋白酶-印地那韋復(fù)合物時，體系的范德華相互作用（VDWIN）和非極性溶劑化能（VDWPB）都有利于底物的結(jié)合；體系的靜電相互作用則是不利于底物的結(jié)合，主要表現(xiàn)在盡管真空條件下分子靜電能（ELEIN）有利于底物的結(jié)合，但是極性溶劑化能（ELEPB）極其不利于二者識別，并直接反轉(zhuǎn)了ELEIN有利于結(jié)合的趨勢。體系熵效應(yīng)（TS）是不利于底物結(jié)合的，這與分子間識別基本都是熵減效應(yīng)一致。

2.2.3 能量分解及結(jié)合模式 為了研究HTLV-1蛋白酶與抑制劑印地那韋的相互作用，同時考慮到用GB模型計算靜電省時間的因素，用GBSA方法計算了蛋白酶結(jié)合印地那韋所涉及到的重要殘基能量貢獻。

表3給出了HTLV-1蛋白酶中與抑制劑印地那韋識別密切相關(guān)殘基的能量分解結(jié)果（能量小于-0.5 kJ/mol），表中的右上角撇號表示B鏈。從表3可以看出，利于HTLV-1蛋白酶與印地那韋結(jié)合的關(guān)鍵殘基主要分布在4個區(qū)域：R10、L30～V39、 V56～F67、W98～I100。其中A鏈和B鏈有利于結(jié)合印地那韋的關(guān)鍵殘基分別是10和11個。其中在A鏈和B鏈都出現(xiàn)的關(guān)鍵殘基有7個，具體為R10、G34、A35、V56、G58、W98和I100。盡管蛋白質(zhì)A鏈和B鏈空間上中心對稱，小分子也結(jié)合在正中間，但是由于小分子結(jié)構(gòu)并不對稱，具有一定的空間異構(gòu)性，A鏈及B鏈對于結(jié)合印地那韋的貢獻類似，但不盡相同，該結(jié)果與之前RMSD（圖3）及RMSF（圖4）分析結(jié)果一致。B鏈的F67也是促進蛋白質(zhì)與小分子識別的關(guān)鍵殘基，這個殘基在晶體結(jié)構(gòu)[28]中并沒有被提到。

基于上面的能量分解結(jié)果，圖7給出了HTLV-1蛋白酶與印地那韋的識別模式。圖中藍色和綠色stick模型分別表示關(guān)鍵殘基，蛋白質(zhì)整體用flat-ribbon模型表示。

2.3 分子的運動性

用MMTSB軟件包[25]對MD模擬軌跡11 000個構(gòu)象進行了成簇分析。圖8給出了基于相互RMSD對比的復(fù)合物體系整體構(gòu)象成簇分析結(jié)果，其中RMSD閾值設(shè)為0.13 nm。MD模擬構(gòu)象共分為5簇，2.5 ns之前基本為第1簇，共2 460個構(gòu)象，所占概率為22.36%；2.5～6.3 ns之間為第2簇，共3 484個構(gòu)象，占比31.67%；第3簇則主要是出現(xiàn)在9.0 ns后，構(gòu)象數(shù)和占比分別為1 447和13.15%；第4簇和第5簇的時間區(qū)間基本為6.0～9.0 ns之間。值得一提的是，第3簇、第4簇和第5簇構(gòu)象有較為明顯的交叉躍遷現(xiàn)象，尤其第4簇和第5簇之間。

圖9給出了HITLV-1蛋白酶與抑制劑印地那韋復(fù)合物體系的每簇代表性構(gòu)象的疊落結(jié)果。代表性結(jié)構(gòu)指的是每個體系各簇中能量最低構(gòu)象，即第1～5簇代表性構(gòu)象代碼分別為1 754、4 220、10 998、8 410和7 745。圖9A給出了基于蛋白質(zhì)主鏈Cα原子重疊的5個代表性構(gòu)象疊落結(jié)果。蛋白質(zhì)A，B鏈部分的Flap、top-site以及兩個exo-site位點疊落較差，相應(yīng)區(qū)域的柔性非常大。藥物篩選研究試驗[29]表明，flap、top-site、exo-site等長loop區(qū)較高的柔性是抑制劑與蛋白酶結(jié)構(gòu)不夠穩(wěn)定的結(jié)果。印地那韋抑制劑與HTLV-1蛋白酶抑制劑的結(jié)合力要弱于與HIV-1蛋白酶的識別，這能夠部分解釋loop區(qū)柔性為何較高。

另外，在代表性構(gòu)象中首先保證小分子哌嗪環(huán)（N1-C1-C2-N3-C8-C9）重疊，然后進行疊落操作，結(jié)果如圖9B所示。蛋白質(zhì)相應(yīng)區(qū)域都是中心對稱，由于印地那韋的不對稱導(dǎo)致其與蛋白酶結(jié)合和運動具有一定的復(fù)雜度。印地那韋的吡啶、哌嗪環(huán)和叔丁基主要與B鏈結(jié)合，運動幅度較小，這與之前的RMSD和RMSF分析一致；印環(huán)基則主要靠近A鏈，疊落效果最差，運動幅度也最高；而苯環(huán)基主要靠近B鏈的R10殘基，柔性居中。

3 結(jié)論

基于HTLV-1蛋白酶和抑制劑印地那韋的復(fù)合物長時間顯含水的分子動力學(xué)模擬，分析了二者的相互作用機制。體系的回轉(zhuǎn)半徑和方均根偏差分析表明，與晶體結(jié)構(gòu)相比，HTLV-1蛋白酶-印地那韋復(fù)合物構(gòu)象在水溶液中逐漸松弛，且A鏈柔性略高于B鏈。蛋白質(zhì)A鏈及B鏈的B因子和RMSF數(shù)值有很高的相關(guān)性，相關(guān)系數(shù)分別達到0.56和0.71，表明柔性分析合理。柔性大的區(qū)域都是遠離結(jié)合部位的無規(guī)卷曲或β片；剛性區(qū)則是直接參與識別抑制劑印地那韋，也包括附近剛性較大的α螺旋區(qū)。能量分解表明，HTLV-1蛋白酶與印地那韋結(jié)合的關(guān)鍵殘基主要分布在R10、L30～V39、V56～F67和W98～I100等4個區(qū)域。用MM-PBSA方法獲得HTLV-1蛋白酶和抑制劑印地那韋的結(jié)合自由能計算值（-32.87 kJ/mol）與試驗值（-31.33 kJ/mol）吻合較好，其中，體系的范德華相互作用和非極性溶劑化能都有利于二者結(jié)合運動性分析表明，蛋白質(zhì)A，B鏈部分的Flap、top-site以及兩個exo-site位點疊落較差，相應(yīng)區(qū)域的柔性非常大；而印地那韋抑制劑的印環(huán)基則主要靠近柔性較大的A鏈，疊落效果最差，運動幅度也最高。模擬結(jié)果為基于HTLV-1蛋白酶受體結(jié)構(gòu)的印地那韋類抑制劑設(shè)計有一定的指導(dǎo)意義。

參考文獻：

[1] LI M， LACO G S， JASKOLSKI M， et al. Crystal structure of human T cell leukemia virus protease， a novel target for anticancer drug design[J]. Proc Natl Acad Sci USA， 2005， 102： 18332-18337.

[2] MCGILL N K，VYAS J， SHIMAUCHI T， et al. HTLV-1-associated infective dermatitis： Updates on the pathogenesis[J]. Exp Dermatol， 2012，21：815-821.

[3] ISHIKAWA T. Current status of therapeutic approaches to adult T-cell leukemia[J]. Int J Hematol， 2003， 78： 304-311.

[4] WLODAWER A. Crystal structures of inhibitor complexes of human T-cell leukemia virus （HTLV-1） protease[J]. J Mol Biol， 2010， 401： 626-641.

[5] SHUKER S B， MARIANI V L， HERGER B E， et al. Understanding HTLV-1 protease[J]. Chem Biol，2003，10：373-380.

[6] MEHELLOU Y， DE CLERCQ E. Twenty-six years of anti-HIV drug discovery： Where do we stand and where do we go [J]. J Med Chem， 2010， 53： 521-538.

[7] KADAS J， WEBER I T， BAGOSSI P， et al. Narrow substrate specificity and sensitivity toward ligand-binding site mutations of human T-cell leukemia virus type 1 protease[J]. J Biol Chem， 2004，279： 27148-27157.

[8] BAGOSSI P， KADAS J， MIKLOSSY G， et al. Development of a microtiter plate fluorescent assay for inhibition studies on the HTLV-1 and HIV-1 proteinases[J]. J Virol Methods， 2004， 119： 87-93.

[9] MAREN K， HOLGER S， WIBKE E， et al. Structural basis for HTLV-1 protease inhibition by the HIV-1 protease inhibitor indinavir[J]. J Med Chem， 2014， 57： 6266-6272.

[10] WANG J M， CIEPLAK P， KOLLMAN P A. How well does a restrained electrostatic potential （RESP） model perform in calculating conformational energies of organic and biological molecules[J]. J Comput Chem， 2000， 21（12）： 1049-1074.

[11] WANG J M， WOLF R M， CALDWELL J W， et al. Development and testing of a general amber force field[J]. J Comput Chem， 2004， 25（9）： 1157-1174.

[12] JORGENSEN W L， CHANDRASEKHAR J， MADURA J D， et al. Comparison of simple potential functions for simulating liquid water [J]. J Chem Phys， 1983， 79 （2）： 926-935.

[13] RYCKAERT J P， CICCOTTI G， BERENDSEN H J C. Numerical integration of the cartesian equations of dynamics of n-alkanes[J]. J Comput Phys， 1977， 23（3）： 327-341.

[14] KOLLMAN P A， MASSOVA I， REYES C， et al. Calculating structures and free energies of complex molecules： Combining molecular mechanics and continuum models[J]. Acc Chem Res， 2000， 33： 889-897.

[15] WANG W， DONINI O， REYES C， et al. Biomolecular simulations： Recent developments in force fields， simulations of enzyme catalysis， protein-ligand， protein-protein， and protein-nucleic acid noncovalent interactions[J]. Annu Rev Biophys Biomol Struct， 2001， 30： 211-243.

[16] HONIG B， NICHOLLS A. Classical electrostatics in biology and chemistry[J]. Science，1995，268：1144-1149.

[17] HECHT J L， HONIG B， SHIN Y K， et al. Electrostatic potentials near the surface of DNA： Comparing theory and experiment[J]. J Phys Chem， 1995， 99： 7782-7786.

[18] SANNER M F， OLSON A J， SPEHNER J. Reduced surface： An efficient way to compute molecular surfaces[J]. Biopolymers， 1996， 38： 305-320.

[19] SITKOFF D， SHARP K A， HONIG B. Accurate calculation of hydration free-energies using macroscopic solvent models[J]. J Phys Chem， 1994， 98： 1978-1988.

[20] KOTTALAM J， CASE D A. Langevin modes of macromolecules-applications to crambin and DNA hexamers[J]. Biopolymers， 1990， 29： 1409-1421.

[21] SIMONSON T. Macromolecular electrostatics： Continuum models and their growing pains[J]. Curr Opin Struct Biol， 2001， 11（2）： 243-252.

[22] BASHFORD D， CASE D A. Generalized born models of macromolecular solvation effects[J]. Ann Rev Phys Chem， 2000， 51： 129-152.

[23] STILL W C， TEMPCZYK A， HAWLEY R C， et al. Semianalytical treatment of solvation for molecular mechanics and dynamics[J]. J Am Chem Soc， 1990， 112（16）： 6127-6129.

[24] WEISER J， SHENKIN P S， STILL W C. Approximate atomic surfaces from linear combinations of pairwise overlaps （LCPO）[J]. J Comput Chem， 1999， 20（2）： 217-230.

[25] MICHAEL F， John K， Charles L， et al. MMTSB tool set： Enhanced sampling and multiscale modeling methods for applications in structural biology[J]. J Mol Graph Model， 2004， 22： 377-395.

[26] DAURA X， GADEMANN K， JAUN B， et al. Peptide folding： When simulation meets experiment[J]. Angew Chem Int Ed， 1999， 38： 236-240.

[27] HU J P， GONG X Q， SU J G， et al. Study on the molecular mechanism of inhibiting HIV-1 integrase by EBR28 peptide via molecular modeling approach[J]. Biophy Chem， 2008， 132： 69-80.

[28] MAREN K， HOLGER S， WIBKE E D. Structural basis for HTLV-1 protease inhibition by the HIV-1 protease inhibitor indinavir[J]. J Med Chem， 2014， 57： 6266-6272.

[29] PERRYMAN A L， ZHANG Q， SOUTTER H H， et al. Fragment-based screen against HIV protease[J]. Chem Biol Drug Des， 2010， 75： 257-268.