β分泌酶是治療阿爾茨海默病(AD)的理想作用靶點。采用以R基團技術為核心的Topomer CoMFA研究HEA類β分泌酶抑制劑的三維定量構效關系(3D-QSAR),構建了擬合與預測能力良好的3D-QSAR模型,得到擬合、交叉與外部驗證的復相關系數分別為r2=0.928,qloo2=0.605和rpred2=0.626。通過3D-QSAR模型搜索ZINC化合物結構片段源,得到活性貢獻提高的R基團并結合公共骨架設計得到15個新穎化合物,其預測活性值均優于訓練集中的活性最高分子。用分子對接研究新設計化合物與β分泌酶的相互作用模式,結果表明,氫鍵和疏水性是影響親和力的重要因素。研究表明,基于R基團的Topomer CoMFA與Topomer Search可以有效地篩選和設計HEA類β分泌酶抑制劑,所設計的分子為AD藥物的研發提供了新的候選物。
引用本文： 史博智, 劉永瀾, 李月婷, 王貴學, 梁桂兆. R基團搜索策略用于HEA類β分泌酶抑制劑的三維定量構效關系研究與分子虛擬篩選. 生物醫學工程學雜志, 2014, 31(1): 196-204. doi: 10.7507/1001-5515.20140038 復制
引言
阿爾茨海默病 (Alzheimer’s disease，AD)，俗稱老年癡呆癥，是老齡人口中發病率最高的疾病之一[1]，可造成患者記憶、思維、情緒等方面的障礙。然而，目前尚無特效藥物治療AD。AD的特征性病理變化為大腦皮層萎縮，并伴有β淀粉樣蛋白(amyloid-β,Aβ)沉積、神經原纖維纏結，以及老年斑的形成等[2]。近年來，隨著對AD病理機制研究的逐漸深入，越來越多的證據表明Aβ蛋白在AD致病分子機制中扮演著重要角色[3-5]。Aβ的異常增多可致細胞損傷乃至死亡。因此，抑制Aβ的產生，是具有前景的AD治療藥物研發策略。
Aβ是由細胞膜上淀粉樣蛋白前體(amyloid precursor protein,APP)經γ分泌酶和β分泌酶水解產生的一系列長短不等的多肽。由于γ酶本身具促細胞生長等生理作用，抑制γ酶的活性是否會妨礙其它重要功能尚未可知[6]。相比而言，β分泌酶抑制劑發展相對成熟，且已有實驗表明抑制小鼠腦內β分泌酶活性不會產生嚴重的毒副作用[7]。BACEl和BACE2是β分泌酶的兩種形式[8]。BACE2在腦組織中表達很少，在Aβ的產生中不占重要地位[9]；而BACEl是腦內β分泌酶的主要存在形式，具有β分泌酶的所有特性，并控制催化Aβ產生的初始和控速步驟[8]。因此，BACE1成為近年來尋找AD治療藥物的重要靶標[10]。
BACE1屬于天門冬氨酸蛋白酶家族[8]，根據其催化裂解機制，目前發現的抑制劑有肽和小分子兩類[6]。肽類BACE1抑制劑具有毒副作用小、生物相容性好等特點，在BACE1抑制劑研究中占據重要地位[11]，按其作用的核心結構不同，可以分為Hydroxyethylene[12]、Hydroethylamine(HEA)[13]和Carbinamine[14]等幾類。這些抑制劑大部分是BACE1催化過渡態的類似物[8]，通過與活性中心兩個天冬氨酸殘基及其他基團相互作用而抑制其活性。
HEA[13]類BACE1抑制劑表現出高效和細胞膜滲透性良好等優勢。Freskos等[15]利用基于結構的藥物設計，改善了之前抑制劑代謝穩定性不足的缺點，并提高了對BACE1的選擇性。同屬HEA類的GSK188909是由Hussain等[16]研發的有效抑制劑，已成功使APP轉基因小鼠腦內Aβ減少，不僅具有更強的抑制效果，而且對同是天冬氨酸酶家族的BACE2以及CAT-D有較好的選擇性。隨后，根據電子等排原理和藥物的構效關系理論，Clarke等[17-20]也設計出類似的具有口服生物活性的高效抑制劑。近年來，越來越多的HEA類BACE1抑制劑被不同研究小組所發現，都具有相當不錯的抑制效果。因此，HEA類是最具前景的β分泌酶抑制劑之一。
HEA類抑制劑一般具有相似的結構骨架以及可變的R基團，因而可以根據基于片段的藥物設計方法對其進行研究。為了闡明HEA抑制劑中R基團的結構與活性的關系以及獲取在理論上更優的化合物，本文采用Topomer CoMFA對HEA類分子進行三維定量構效關系(three-dimensional quantitative structure-activity relationship，3D-QSAR)研究，在此基礎上利用Topomer Search進行虛擬篩選并設計新的分子。通過Surflex-dock對新設計的配體分子與BACE1的作用機制進行合理解釋。最終為HEA類β分泌酶抑制劑的分子設計及新型AD藥物的研發提供參考。
1 原理與方法
1.1 數據集
所選的HEA類抑制劑的結構及活性數據來自文獻[17-20]，該數據集集合了Clark等同一研究小組的4篇研究成果，在相同的實驗條件下測定活性。52個分子具有共同的結構骨架以及影響分子活性的3個R基團(圖 1)。該類抑制劑在體外試驗已達到納摩爾水平，采用其IC50的負對數值(pIC50)作為活性標度，其中IC50值為化合物對BACE1的體外半數最大抑制濃度 (單位：nmol)。將pIC50按照從小到大排序，按照隔二取一法，得到訓練集分子35個，測試集分子17個。分子的三維結構由Sybyl 8.1創建，分子的能量優化采用最陡梯度下降法，力場為Tripos力場，加載Gasteiger-Huckel電荷，迭代1 000次，能量收斂限定為0.001 kJ/mol，優化后得到穩定的構象。除特殊說明外，其它參數均采用Sybyl 8.1默認設置。
圖1
HEA 類抑制劑的結構骨架
Figure1.
Sketch of structure of HEA derivatives
1.2 Topomer CoMFA建模
采用Topomer CoMFA[21]對訓練集建模。為了保證核心基團的一次性疊合以得到較好的建模結果，把分子切割為一個公共骨架(結構骨架中除去靠近酰胺鍵的苯環)和Ra、Rb兩個基團。以訓練集中的5號樣本為模板，切割樣式如圖 2所示。切割完成后Topomer CoMFA自動計算靜電和立體場參數，經分子疊合及偏最小二乘構建3D-QSAR模型。采用留一法(leave one out,LOO)交互檢驗和外部驗證評價模型的內外部預測能力。
圖2
樣本5的Topomer CoMFA切割式樣
Figure2.
Cutting style of sample 5 the green,blue and red represent the core group,Ra and Rb in sample 5,respectively
1.3 Topomer Search虛擬篩選
基于Topomer CoMFA模型進行Topomer Search研究。分別以Ra、Rb為提問結構，設定評價結構相似性程度的Topomer距離值在185以內命中，搜索對象為ZINC(2012)[22]數據庫中的Leadlike standard類[23]，原數據庫共有4 554 059個分子。該類片段源分子理化特征滿足以下標準：2.5≤logP (脂水分配系數對數)＜3.5、250≤mol.wt. (分子量)＜350、5≤rotatable bond count (可旋轉鍵計數) ≤7。選擇原數據庫中編號0-9的文件共534 597個化合物，作為本文的目標搜索數據庫，進行虛擬篩選。
1.4 Surflex-dock分子對接
采用Sybyl8.1中的Surflex-dock方法研究新設計的配體分子與BACEl的作用模式和機制。對接研究采用的BACEl晶體結構下載于PDB數據庫(PDB code：2VJ6)，自帶配體VG5_1447。對蛋白質結構進行預處理：刪除所有結晶水，抽取出配體分子。修正蛋白質晶體結構的化學鍵：選擇末端質子化狀態為charge，并為受體加全氫，氫取向設為H-bonding，以便形成氫鍵。以自帶配體生成結合口袋。其他均為默認設置。小分子結構構建方法同訓練集。將新設計的分子與準備好的蛋白質受體對接，每個分子都能輸出結合后20個構象的經驗打分。打分包括Total score、G-score[24]、D-score[25]、Chem-score[26]、PMF-score[27]以及反映前五者打分一致性的C-score(consensus score)。通過打分函數高低結合MOLCAD作圖模塊，對分子與配體結合構象進行研究。
2 結果與討論
2.1 Topomer CoMFA建模結果
作為第二代的CoMFA[28]技術，Topomer CoMFA可以不需人工疊合配體，減少了3D-QSAR分析的準備需求，避免了CoMFA模塊中分子疊合的弊端[21]。Topomer[29]能夠自動根據三維構象生成碎片，由一定經驗快速完成片段的構象選擇和疊合。再結合CoMFA原理，以分子場參數為自變量，分子生物活性(pIC50)為因變量，采用偏最小二乘構建3D-QSAR模型。
以35個訓練集樣本建立Topomer CoMFA模型，得到其主成分數(components)為8，截距(intercept)為5.66，建模的復相關系數r2=0.928，LOO的復相關系數qloo2=0.605，估計標準誤差SEE為0.255，F值為41.943。其中，截距為一個參考值，其與每個基團活性貢獻值之和即為分子的活性預測值。模型對17個測試集樣本預測的復相關系數qpred2=0.626。可以看出該模型具有良好的擬合及預測能力。
2.2 Topomer CoMFA 模型等勢圖分析
圖 3為Topomer CoMFA模型的立體場和靜電場的三維等勢圖，為了清晰地顯示，Ra、Rb等勢圖只選取了部分活性較高的分子進行疊合。可以根據等勢圖提出相應的改造意見以提高活性：立體場中黃色區域表示適宜引入體積小或空間位阻小的取代基，綠色區域表示適宜引入體積大或空間位阻大的取代基。靜電場中藍色區域表示適宜引入帶正電性的取代基，紅色區域表示適宜引入帶負電的取代基。
圖3
Topomer CoMFA模型的等勢圖
(a)Ra立體場等勢圖；(b)Ra靜電場等勢圖；(c)Rb立體場等勢圖；(d)Rb靜電場等勢圖
Figure3. 3D contour plot of Topomer CoMFA model(a) steric field map of Ra; (b) electrostatic field map of Ra; (c) steric field map of Rb; (d) electrostatic field map of Rb
結合圖 2切割樣式分析，在Ra立體場等勢圖a中，R1苯環下方及R2遠離公共骨架一側各有一片較小和較大的黃色區域；R1苯環上方及R2靠近公共骨架一側各有一片較大和較小綠色區域。在Ra靜電場等勢圖b中，R1一側有大片的紅色區域；R2一側有兩片藍色區域。在Rb立體場等勢圖c中，接公共骨架的遠端有一片大面積的綠色區域，近端一側有幾片較小的黃色區域。在Rb立體場等勢圖c中，接公共骨架的近端有幾片黃色區域，遠端有一片大面積的綠色區域。Rb靜電場等勢圖d中，接公共骨架遠側各有幾片紅色區域，近側有一片藍色區域。由上述分析可知，R1適宜引入體積大、帶負電性的基團；R2適宜引入帶正電性的基團；R3近端適宜引入體積小、帶正電的基團，遠端適宜引入體積大、帶負電的基團。
疊合圖中不同性質基團對應的色塊位置符合等勢圖分析。R1為大體積的環狀基團，對應大塊的綠色區域，且電負性大的氧原子所對方向正好為較大的紅色區域。R2基團中大部分為脂肪胺或包含脂肪烴，為斥電子基團，該區域帶正電荷，對應藍色區域。R3遠端三氟化碳體積大，帶負電，分別被綠色與紅色色塊包裹。
通過對比也可驗證上述分子場分析與基團活性貢獻值的關系。對比樣本5與樣本45發現，R2結構相同，樣本45的R1取代基增加了一個側鏈，增大了體積。因此樣本45的 Ra貢獻值更大。對比樣本40與樣本51可知，樣本51的三氟甲基使遠端取代基體積增大，且帶負電，都有利于活性提升，因此樣本51的Rb貢獻值明顯大于樣本21。
2.3 Topomer Search
Topomer Search[30]是一種基于3D配體的虛擬篩選工具，能夠基于Topomer CoMFA模型，根據結構相似性從含有大量片段源的數據庫中預測識別出有利于生物活性提高的R基團。Topomer Search 的結果中主要包含兩項：Topomer距離和R基團的活性貢獻值。Topomer距離通過計算值評價化合物與提問結構基團的結構相似性程度；而R基團活性貢獻值則是基于Topomer CoMFA模型對R基團活性值的預測打分。一般情況下，在設定的結構相似性程度范圍內都會優先考慮貢獻值的大小，以用于基于配體的分子優化設計以及虛擬篩選。
基于建立的Topomer CoMFA模型對指定的分子數據庫進行R基團的搜索，分別得到2 039個Ra結構和5 000個Rb結構。以活性最高的樣本50的Ra、Rb貢獻值對命中分子進行過濾。分析Ra命中結果，發現命中結果中的貢獻值并未明顯優于樣本50的Ra貢獻值：所有結果中僅有5個命中分子貢獻值大于樣本50 的Ra貢獻值，且最高結果僅比樣本50提高12.7%。為了還原樣本的結構骨架，篩選結果還應包含苯環且苯環上含有兩個互為間位的取代基。在滿足條件的結果中，按貢獻值大小，篩選得到滿足條件的前三個基團a1、a2和a3。對于Rb命中的結果，發現5 000個命中結果貢獻值均比樣本50的Rb貢獻值高，最低結果提高29%，最高結果提高135%。根據搜索結果，Rb有很大的優化改造空間。按貢獻值大小，篩選得到排名前5的基團b1、b2、b3、b4和b5。篩選出的Ra、Rb的結構，Topomer距離，Ra及Rb貢獻值如表 1所示。
表1
新設計分子的結構、Topomer距離、R貢獻值和活性預測值
Table1.
Structures of newly-designed molecules,Topomer distances,contribution values of R groups,and predicted pIC50
2.4 分子設計
以切割時的公共骨架為設計模板，用篩選后的Ra、Rb基團分別對應模板分子的Ra和Rb基團。得到15個新分子的結構及預測活性標度pIC50(表 1)，可看出新分子pIC50均明顯大于文獻中活性最高的樣本50分子的活性。
經觀察可以發現,表 1中新設計分子的R基團與樣本分子的R基團有明顯變化。R1基團引入苯環或環狀化合物使體積明顯增大，引入了高電負性氧原子使取代基帶負電。R3都是大體積的取代基，且遠端都包含N、O、Cl等電負性高的原子形成帶負電取代基。這些與Topomer CoMFA模型等勢圖的分析相一致。
表 1中部分新化合物Rb的貢獻值，搜索完成時的結果與組合后分子活性預測時的結果不一致，體現在b2、b3、b4、b5的Rb均小于模板分子。這說明雖然篩選時著重考慮基團貢獻值大小，但改造后整體分子的活性預測值降低，表明修改R1這個基團對于整個分子活性的提高是有風險的，可能是它與核心基團不匹配的原因，或者是組合后的分子內部存在不利的相互作用。
另外值得注意的是改造后的R2基團。NO2和SO2CH3都是吸電子的帶負電基團，與靜電場等勢圖藍色標塊的改造意見不符，由此造成Ra命中結果貢獻值不高。分析其原因，可能是因為建模時Ra包含R1、R2、苯環的這種切割方式并未單獨考慮結構骨架中R1、R2基團的構效關系，從而增大了Topomer Search基于相似性搜索R1和R2的難度，降低了命中率，未能搜索出能夠提高貢獻值的R2基團。搜索結果的Topomer 距離(a1=181,a2=181,a3=125)較大也可以看出這一點。為了搜索更優的R1和R2，可以在進行Topomer Search的時候適當增大Topomer距離，即降低命中分子與提問結構的相似程度來提高命中率，進而有機會篩選出活性貢獻值高的基團。
2.5 新設計分子與蛋白受體作用機制
分子對接是一種能夠結合配體和受體雙方相互作用模式對藥物的作用機制進行全面合理解釋的研究方法,并能提供多種小分子結合構象，輔助3D-QSAR分析。Surflex-dock 的打分函數[31]采用-lg(kd)為單位表示結合親和力，主要對以下幾方面性質進行考量：疏水性、極性、排斥力、熵、溶劑化效應、碰撞效應。Total score是對以上性質的總體打分，一般地，分值越高，表示產生的構象越好，大于4即有一定活性，大于6則活性能夠達到微摩爾級別。不同評分函數都具有某種程度的局限性，然而研究證明，C-score能夠反映5個不同經驗函數Total score、G-score、D-score、Chem-score和PMF-score得分的一致性：每一項達標得1分，滿分為5分。因此，與單一評分函數比較，C-score的使用可以大大提高篩選的成功率[32-33]。本文以C-score滿分為標準篩選出對接評價優秀的設計分子，并在每個分子的構象中選擇Total score最高的確定為最佳構象。
表 2中顯示了滿足篩選條件C-score=5的1、3、5、8、13號分子以及與之對應的最佳構象打分。可以認為在C-score相同的情況下，Total score越高，分子構象越符合實際。因此選出Total score打分最高的3號分子構象進行研究。3號分子與2VJ6活性位點對接如圖 4所示。圖 5為簡化后3號分子與活性位點殘基形成的氫鍵示意圖。可以發現3號分子與活性中心一共形成了9個氫鍵，其中與Asn294和Thr293形成的氫鍵與文獻[34]報道一致。圖 6為Sybyl MOLCAD顯示的2VJ6結合口袋表面的氫鍵密度分布圖。R1、R2周圍多為紅色區域，代表適合引入氫鍵配體；R3周圍多為藍色區域，代表適合引入氫鍵受體。通過觀察氫鍵密度分布圖在活性口袋內部的情況，結構骨架上C=O對應的是紅色區域，OH對應的是藍色區域。結合圖 5發現，這都與3號分子和活性位點間的氫鍵對接模式相匹配。
表2
篩選后分子對接的打分結果
Table2.
Surflex-dock scores of newly-designed molecules after screening
圖4
3號分子與2VJ6活性位點對接圖
黃色虛線表示氫鍵
Figure4. Docking results of compound 3 with 2VJ6dotted yellow lines indicate hydrogen bonds
圖5
3號分子與2VJ6作用的氫鍵、疏水性示意圖
虛線表示氫鍵，前后數字分別表示氫鍵長度及角度，單位分別為埃和度
Figure5. Docking results of compound 3 with 2VJ6dotted lines indicate hydrogen bonds,the former and latter numbers indicate the length and the angle of hydrogen bonds,the unit are angstrom and angle,respectively
圖6
2VJ6結合口袋處MOLCAD氫鍵密度表面圖
紅色和藍色區域分別表示氫鍵供體和受體區域
Figure6. The MOLCAD of hydrogen bond density in binding pocket of 2VJ6the red and blue areas are fit for donors and acceptors,respectively
圖 7為2VJ6結合口袋表面的疏水性分布圖。棕色和藍色分別代表活性位點的疏水性區域和親水性區域。不難看出，在結構骨架的OH附近有一片藍色區域，在R3位有一片棕色區域。進一步研究發現，藍色區域附近主要以親水性殘基為主，如Ser290、Gly291。棕色區域附近主要以疏水性殘基為主，如Trp258、Tyr259、Tyr260等。因此，在R3位用疏水性基團取代有利于親和活性的提高。如樣本3的R3三元環，用1-methyl-2,3-dihydro-1H-indene取代后得到樣本41，其pIC50提高了0.55。又如17的R3位增加一條支鏈1-methoxy-3-methylbenzene得到樣本27，其pIC50提高了1.59。可以推斷疏水性是影響親和活性的重要因素。
圖7
2VJ6結合口袋處MOLCAD親脂性表面圖
棕色和藍色區域分別代表疏水和親水區域
Figure7. The MOLCAD of lipophilicity in binding pocket of 2VJ6the brown and blue areas indicate hydrophobic and hydrophilic areas,respectively
3 結論
采用Topomer CoMFA 與Topomer Search聯用進行虛擬篩選，設計了新型BACE1抑制劑分子，采用Surflex-dock分子對接直觀地給出新設計分子配體與蛋白受體結合模式和構效關系，進一步對Topomer CoMFA模型進行了驗證。本研究為BACE1抑制劑分子設計和結構改造提供了較實用的方法學和理論依據。本文采用R基團搜索技術設計了抑制劑新分子，從基于配體和受體兩方面證明了所設計分子的有效性。但是，我們在采用R基團組合新分子時，并沒有考慮R之間或者R與核心基團的作用特征，今后的工作可以繼續研究這兩方面的作用，深入解釋新分子的作用機制與模式；另一方面，基于Topomer CoMFA 、Topomer Search和Surflex-dock分子對接的分子設計技術，還需要與藥物化學實驗以及生物學實驗緊密合作，以從理論和實踐兩方面實現有效的分子設計。
引言
阿爾茨海默病 (Alzheimer’s disease，AD)，俗稱老年癡呆癥，是老齡人口中發病率最高的疾病之一[1]，可造成患者記憶、思維、情緒等方面的障礙。然而，目前尚無特效藥物治療AD。AD的特征性病理變化為大腦皮層萎縮，并伴有β淀粉樣蛋白(amyloid-β,Aβ)沉積、神經原纖維纏結，以及老年斑的形成等[2]。近年來，隨著對AD病理機制研究的逐漸深入，越來越多的證據表明Aβ蛋白在AD致病分子機制中扮演著重要角色[3-5]。Aβ的異常增多可致細胞損傷乃至死亡。因此，抑制Aβ的產生，是具有前景的AD治療藥物研發策略。
Aβ是由細胞膜上淀粉樣蛋白前體(amyloid precursor protein,APP)經γ分泌酶和β分泌酶水解產生的一系列長短不等的多肽。由于γ酶本身具促細胞生長等生理作用，抑制γ酶的活性是否會妨礙其它重要功能尚未可知[6]。相比而言，β分泌酶抑制劑發展相對成熟，且已有實驗表明抑制小鼠腦內β分泌酶活性不會產生嚴重的毒副作用[7]。BACEl和BACE2是β分泌酶的兩種形式[8]。BACE2在腦組織中表達很少，在Aβ的產生中不占重要地位[9]；而BACEl是腦內β分泌酶的主要存在形式，具有β分泌酶的所有特性，并控制催化Aβ產生的初始和控速步驟[8]。因此，BACE1成為近年來尋找AD治療藥物的重要靶標[10]。
BACE1屬于天門冬氨酸蛋白酶家族[8]，根據其催化裂解機制，目前發現的抑制劑有肽和小分子兩類[6]。肽類BACE1抑制劑具有毒副作用小、生物相容性好等特點，在BACE1抑制劑研究中占據重要地位[11]，按其作用的核心結構不同，可以分為Hydroxyethylene[12]、Hydroethylamine(HEA)[13]和Carbinamine[14]等幾類。這些抑制劑大部分是BACE1催化過渡態的類似物[8]，通過與活性中心兩個天冬氨酸殘基及其他基團相互作用而抑制其活性。
HEA[13]類BACE1抑制劑表現出高效和細胞膜滲透性良好等優勢。Freskos等[15]利用基于結構的藥物設計，改善了之前抑制劑代謝穩定性不足的缺點，并提高了對BACE1的選擇性。同屬HEA類的GSK188909是由Hussain等[16]研發的有效抑制劑，已成功使APP轉基因小鼠腦內Aβ減少，不僅具有更強的抑制效果，而且對同是天冬氨酸酶家族的BACE2以及CAT-D有較好的選擇性。隨后，根據電子等排原理和藥物的構效關系理論，Clarke等[17-20]也設計出類似的具有口服生物活性的高效抑制劑。近年來，越來越多的HEA類BACE1抑制劑被不同研究小組所發現，都具有相當不錯的抑制效果。因此，HEA類是最具前景的β分泌酶抑制劑之一。
HEA類抑制劑一般具有相似的結構骨架以及可變的R基團，因而可以根據基于片段的藥物設計方法對其進行研究。為了闡明HEA抑制劑中R基團的結構與活性的關系以及獲取在理論上更優的化合物，本文采用Topomer CoMFA對HEA類分子進行三維定量構效關系(three-dimensional quantitative structure-activity relationship，3D-QSAR)研究，在此基礎上利用Topomer Search進行虛擬篩選并設計新的分子。通過Surflex-dock對新設計的配體分子與BACE1的作用機制進行合理解釋。最終為HEA類β分泌酶抑制劑的分子設計及新型AD藥物的研發提供參考。
1 原理與方法
1.1 數據集
所選的HEA類抑制劑的結構及活性數據來自文獻[17-20]，該數據集集合了Clark等同一研究小組的4篇研究成果，在相同的實驗條件下測定活性。52個分子具有共同的結構骨架以及影響分子活性的3個R基團(圖 1)。該類抑制劑在體外試驗已達到納摩爾水平，采用其IC50的負對數值(pIC50)作為活性標度，其中IC50值為化合物對BACE1的體外半數最大抑制濃度 (單位：nmol)。將pIC50按照從小到大排序，按照隔二取一法，得到訓練集分子35個，測試集分子17個。分子的三維結構由Sybyl 8.1創建，分子的能量優化采用最陡梯度下降法，力場為Tripos力場，加載Gasteiger-Huckel電荷，迭代1 000次，能量收斂限定為0.001 kJ/mol，優化后得到穩定的構象。除特殊說明外，其它參數均采用Sybyl 8.1默認設置。
圖1
HEA 類抑制劑的結構骨架
Figure1.
Sketch of structure of HEA derivatives
1.2 Topomer CoMFA建模
采用Topomer CoMFA[21]對訓練集建模。為了保證核心基團的一次性疊合以得到較好的建模結果，把分子切割為一個公共骨架(結構骨架中除去靠近酰胺鍵的苯環)和Ra、Rb兩個基團。以訓練集中的5號樣本為模板，切割樣式如圖 2所示。切割完成后Topomer CoMFA自動計算靜電和立體場參數，經分子疊合及偏最小二乘構建3D-QSAR模型。采用留一法(leave one out,LOO)交互檢驗和外部驗證評價模型的內外部預測能力。
圖2
樣本5的Topomer CoMFA切割式樣
Figure2.
Cutting style of sample 5 the green,blue and red represent the core group,Ra and Rb in sample 5,respectively
1.3 Topomer Search虛擬篩選
基于Topomer CoMFA模型進行Topomer Search研究。分別以Ra、Rb為提問結構，設定評價結構相似性程度的Topomer距離值在185以內命中，搜索對象為ZINC(2012)[22]數據庫中的Leadlike standard類[23]，原數據庫共有4 554 059個分子。該類片段源分子理化特征滿足以下標準：2.5≤logP (脂水分配系數對數)＜3.5、250≤mol.wt. (分子量)＜350、5≤rotatable bond count (可旋轉鍵計數) ≤7。選擇原數據庫中編號0-9的文件共534 597個化合物，作為本文的目標搜索數據庫，進行虛擬篩選。
1.4 Surflex-dock分子對接
采用Sybyl8.1中的Surflex-dock方法研究新設計的配體分子與BACEl的作用模式和機制。對接研究采用的BACEl晶體結構下載于PDB數據庫(PDB code：2VJ6)，自帶配體VG5_1447。對蛋白質結構進行預處理：刪除所有結晶水，抽取出配體分子。修正蛋白質晶體結構的化學鍵：選擇末端質子化狀態為charge，并為受體加全氫，氫取向設為H-bonding，以便形成氫鍵。以自帶配體生成結合口袋。其他均為默認設置。小分子結構構建方法同訓練集。將新設計的分子與準備好的蛋白質受體對接，每個分子都能輸出結合后20個構象的經驗打分。打分包括Total score、G-score[24]、D-score[25]、Chem-score[26]、PMF-score[27]以及反映前五者打分一致性的C-score(consensus score)。通過打分函數高低結合MOLCAD作圖模塊，對分子與配體結合構象進行研究。
2 結果與討論
2.1 Topomer CoMFA建模結果
作為第二代的CoMFA[28]技術，Topomer CoMFA可以不需人工疊合配體，減少了3D-QSAR分析的準備需求，避免了CoMFA模塊中分子疊合的弊端[21]。Topomer[29]能夠自動根據三維構象生成碎片，由一定經驗快速完成片段的構象選擇和疊合。再結合CoMFA原理，以分子場參數為自變量，分子生物活性(pIC50)為因變量，采用偏最小二乘構建3D-QSAR模型。
以35個訓練集樣本建立Topomer CoMFA模型，得到其主成分數(components)為8，截距(intercept)為5.66，建模的復相關系數r2=0.928，LOO的復相關系數qloo2=0.605，估計標準誤差SEE為0.255，F值為41.943。其中，截距為一個參考值，其與每個基團活性貢獻值之和即為分子的活性預測值。模型對17個測試集樣本預測的復相關系數qpred2=0.626。可以看出該模型具有良好的擬合及預測能力。
2.2 Topomer CoMFA 模型等勢圖分析
圖 3為Topomer CoMFA模型的立體場和靜電場的三維等勢圖，為了清晰地顯示，Ra、Rb等勢圖只選取了部分活性較高的分子進行疊合。可以根據等勢圖提出相應的改造意見以提高活性：立體場中黃色區域表示適宜引入體積小或空間位阻小的取代基，綠色區域表示適宜引入體積大或空間位阻大的取代基。靜電場中藍色區域表示適宜引入帶正電性的取代基，紅色區域表示適宜引入帶負電的取代基。
圖3
Topomer CoMFA模型的等勢圖
(a)Ra立體場等勢圖；(b)Ra靜電場等勢圖；(c)Rb立體場等勢圖；(d)Rb靜電場等勢圖
Figure3. 3D contour plot of Topomer CoMFA model(a) steric field map of Ra; (b) electrostatic field map of Ra; (c) steric field map of Rb; (d) electrostatic field map of Rb
結合圖 2切割樣式分析，在Ra立體場等勢圖a中，R1苯環下方及R2遠離公共骨架一側各有一片較小和較大的黃色區域；R1苯環上方及R2靠近公共骨架一側各有一片較大和較小綠色區域。在Ra靜電場等勢圖b中，R1一側有大片的紅色區域；R2一側有兩片藍色區域。在Rb立體場等勢圖c中，接公共骨架的遠端有一片大面積的綠色區域，近端一側有幾片較小的黃色區域。在Rb立體場等勢圖c中，接公共骨架的近端有幾片黃色區域，遠端有一片大面積的綠色區域。Rb靜電場等勢圖d中，接公共骨架遠側各有幾片紅色區域，近側有一片藍色區域。由上述分析可知，R1適宜引入體積大、帶負電性的基團；R2適宜引入帶正電性的基團；R3近端適宜引入體積小、帶正電的基團，遠端適宜引入體積大、帶負電的基團。
疊合圖中不同性質基團對應的色塊位置符合等勢圖分析。R1為大體積的環狀基團，對應大塊的綠色區域，且電負性大的氧原子所對方向正好為較大的紅色區域。R2基團中大部分為脂肪胺或包含脂肪烴，為斥電子基團，該區域帶正電荷，對應藍色區域。R3遠端三氟化碳體積大，帶負電，分別被綠色與紅色色塊包裹。
通過對比也可驗證上述分子場分析與基團活性貢獻值的關系。對比樣本5與樣本45發現，R2結構相同，樣本45的R1取代基增加了一個側鏈，增大了體積。因此樣本45的 Ra貢獻值更大。對比樣本40與樣本51可知，樣本51的三氟甲基使遠端取代基體積增大，且帶負電，都有利于活性提升，因此樣本51的Rb貢獻值明顯大于樣本21。
2.3 Topomer Search
Topomer Search[30]是一種基于3D配體的虛擬篩選工具，能夠基于Topomer CoMFA模型，根據結構相似性從含有大量片段源的數據庫中預測識別出有利于生物活性提高的R基團。Topomer Search 的結果中主要包含兩項：Topomer距離和R基團的活性貢獻值。Topomer距離通過計算值評價化合物與提問結構基團的結構相似性程度；而R基團活性貢獻值則是基于Topomer CoMFA模型對R基團活性值的預測打分。一般情況下，在設定的結構相似性程度范圍內都會優先考慮貢獻值的大小，以用于基于配體的分子優化設計以及虛擬篩選。
基于建立的Topomer CoMFA模型對指定的分子數據庫進行R基團的搜索，分別得到2 039個Ra結構和5 000個Rb結構。以活性最高的樣本50的Ra、Rb貢獻值對命中分子進行過濾。分析Ra命中結果，發現命中結果中的貢獻值并未明顯優于樣本50的Ra貢獻值：所有結果中僅有5個命中分子貢獻值大于樣本50 的Ra貢獻值，且最高結果僅比樣本50提高12.7%。為了還原樣本的結構骨架，篩選結果還應包含苯環且苯環上含有兩個互為間位的取代基。在滿足條件的結果中，按貢獻值大小，篩選得到滿足條件的前三個基團a1、a2和a3。對于Rb命中的結果，發現5 000個命中結果貢獻值均比樣本50的Rb貢獻值高，最低結果提高29%，最高結果提高135%。根據搜索結果，Rb有很大的優化改造空間。按貢獻值大小，篩選得到排名前5的基團b1、b2、b3、b4和b5。篩選出的Ra、Rb的結構，Topomer距離，Ra及Rb貢獻值如表 1所示。
表1
新設計分子的結構、Topomer距離、R貢獻值和活性預測值
Table1.
Structures of newly-designed molecules,Topomer distances,contribution values of R groups,and predicted pIC50
2.4 分子設計
以切割時的公共骨架為設計模板，用篩選后的Ra、Rb基團分別對應模板分子的Ra和Rb基團。得到15個新分子的結構及預測活性標度pIC50(表 1)，可看出新分子pIC50均明顯大于文獻中活性最高的樣本50分子的活性。
經觀察可以發現,表 1中新設計分子的R基團與樣本分子的R基團有明顯變化。R1基團引入苯環或環狀化合物使體積明顯增大，引入了高電負性氧原子使取代基帶負電。R3都是大體積的取代基，且遠端都包含N、O、Cl等電負性高的原子形成帶負電取代基。這些與Topomer CoMFA模型等勢圖的分析相一致。
表 1中部分新化合物Rb的貢獻值，搜索完成時的結果與組合后分子活性預測時的結果不一致，體現在b2、b3、b4、b5的Rb均小于模板分子。這說明雖然篩選時著重考慮基團貢獻值大小，但改造后整體分子的活性預測值降低，表明修改R1這個基團對于整個分子活性的提高是有風險的，可能是它與核心基團不匹配的原因，或者是組合后的分子內部存在不利的相互作用。
另外值得注意的是改造后的R2基團。NO2和SO2CH3都是吸電子的帶負電基團，與靜電場等勢圖藍色標塊的改造意見不符，由此造成Ra命中結果貢獻值不高。分析其原因，可能是因為建模時Ra包含R1、R2、苯環的這種切割方式并未單獨考慮結構骨架中R1、R2基團的構效關系，從而增大了Topomer Search基于相似性搜索R1和R2的難度，降低了命中率，未能搜索出能夠提高貢獻值的R2基團。搜索結果的Topomer 距離(a1=181,a2=181,a3=125)較大也可以看出這一點。為了搜索更優的R1和R2，可以在進行Topomer Search的時候適當增大Topomer距離，即降低命中分子與提問結構的相似程度來提高命中率，進而有機會篩選出活性貢獻值高的基團。
2.5 新設計分子與蛋白受體作用機制
分子對接是一種能夠結合配體和受體雙方相互作用模式對藥物的作用機制進行全面合理解釋的研究方法,并能提供多種小分子結合構象，輔助3D-QSAR分析。Surflex-dock 的打分函數[31]采用-lg(kd)為單位表示結合親和力，主要對以下幾方面性質進行考量：疏水性、極性、排斥力、熵、溶劑化效應、碰撞效應。Total score是對以上性質的總體打分，一般地，分值越高，表示產生的構象越好，大于4即有一定活性，大于6則活性能夠達到微摩爾級別。不同評分函數都具有某種程度的局限性，然而研究證明，C-score能夠反映5個不同經驗函數Total score、G-score、D-score、Chem-score和PMF-score得分的一致性：每一項達標得1分，滿分為5分。因此，與單一評分函數比較，C-score的使用可以大大提高篩選的成功率[32-33]。本文以C-score滿分為標準篩選出對接評價優秀的設計分子，并在每個分子的構象中選擇Total score最高的確定為最佳構象。
表 2中顯示了滿足篩選條件C-score=5的1、3、5、8、13號分子以及與之對應的最佳構象打分。可以認為在C-score相同的情況下，Total score越高，分子構象越符合實際。因此選出Total score打分最高的3號分子構象進行研究。3號分子與2VJ6活性位點對接如圖 4所示。圖 5為簡化后3號分子與活性位點殘基形成的氫鍵示意圖。可以發現3號分子與活性中心一共形成了9個氫鍵，其中與Asn294和Thr293形成的氫鍵與文獻[34]報道一致。圖 6為Sybyl MOLCAD顯示的2VJ6結合口袋表面的氫鍵密度分布圖。R1、R2周圍多為紅色區域，代表適合引入氫鍵配體；R3周圍多為藍色區域，代表適合引入氫鍵受體。通過觀察氫鍵密度分布圖在活性口袋內部的情況，結構骨架上C=O對應的是紅色區域，OH對應的是藍色區域。結合圖 5發現，這都與3號分子和活性位點間的氫鍵對接模式相匹配。
表2
篩選后分子對接的打分結果
Table2.
Surflex-dock scores of newly-designed molecules after screening
圖4
3號分子與2VJ6活性位點對接圖
黃色虛線表示氫鍵
Figure4. Docking results of compound 3 with 2VJ6dotted yellow lines indicate hydrogen bonds
圖5
3號分子與2VJ6作用的氫鍵、疏水性示意圖
虛線表示氫鍵，前后數字分別表示氫鍵長度及角度，單位分別為埃和度
Figure5. Docking results of compound 3 with 2VJ6dotted lines indicate hydrogen bonds,the former and latter numbers indicate the length and the angle of hydrogen bonds,the unit are angstrom and angle,respectively
圖6
2VJ6結合口袋處MOLCAD氫鍵密度表面圖
紅色和藍色區域分別表示氫鍵供體和受體區域
Figure6. The MOLCAD of hydrogen bond density in binding pocket of 2VJ6the red and blue areas are fit for donors and acceptors,respectively
圖 7為2VJ6結合口袋表面的疏水性分布圖。棕色和藍色分別代表活性位點的疏水性區域和親水性區域。不難看出，在結構骨架的OH附近有一片藍色區域，在R3位有一片棕色區域。進一步研究發現，藍色區域附近主要以親水性殘基為主，如Ser290、Gly291。棕色區域附近主要以疏水性殘基為主，如Trp258、Tyr259、Tyr260等。因此，在R3位用疏水性基團取代有利于親和活性的提高。如樣本3的R3三元環，用1-methyl-2,3-dihydro-1H-indene取代后得到樣本41，其pIC50提高了0.55。又如17的R3位增加一條支鏈1-methoxy-3-methylbenzene得到樣本27，其pIC50提高了1.59。可以推斷疏水性是影響親和活性的重要因素。
圖7
2VJ6結合口袋處MOLCAD親脂性表面圖
棕色和藍色區域分別代表疏水和親水區域
Figure7. The MOLCAD of lipophilicity in binding pocket of 2VJ6the brown and blue areas indicate hydrophobic and hydrophilic areas,respectively
3 結論
采用Topomer CoMFA 與Topomer Search聯用進行虛擬篩選，設計了新型BACE1抑制劑分子，采用Surflex-dock分子對接直觀地給出新設計分子配體與蛋白受體結合模式和構效關系，進一步對Topomer CoMFA模型進行了驗證。本研究為BACE1抑制劑分子設計和結構改造提供了較實用的方法學和理論依據。本文采用R基團搜索技術設計了抑制劑新分子，從基于配體和受體兩方面證明了所設計分子的有效性。但是，我們在采用R基團組合新分子時，并沒有考慮R之間或者R與核心基團的作用特征，今后的工作可以繼續研究這兩方面的作用，深入解釋新分子的作用機制與模式；另一方面，基于Topomer CoMFA 、Topomer Search和Surflex-dock分子對接的分子設計技術，還需要與藥物化學實驗以及生物學實驗緊密合作，以從理論和實踐兩方面實現有效的分子設計。
