基于計(jì)算機(jī)輔助水解的中藥大豆寡肽的ETA拮抗活性預(yù)測(cè)

喬連生+蔣蘆荻+雒剛剛+路芳+陳艷昆+王靈芝+李貢宇+張燕玲

[摘要] 中藥寡肽是中藥發(fā)揮藥效的關(guān)鍵組分之一，系統(tǒng)地研究中藥寡肽的組成及其藥效是中藥物質(zhì)基礎(chǔ)及作用機(jī)制研究的關(guān)鍵。該研究擬基于計(jì)算機(jī)輔助水解和分子對(duì)接等模擬技術(shù)，以中藥大豆為研究載體，解析其降壓寡肽成分，并預(yù)測(cè)其潛在的內(nèi)皮素受體A（endothelin receptor A，ETA）拮抗活性。該文通過(guò)收集大豆中的儲(chǔ)存蛋白序列，基于計(jì)算機(jī)輔助水解方法進(jìn)行蛋白的模擬消化水解，并將水解獲得的寡肽構(gòu)建虛擬結(jié)構(gòu)數(shù)據(jù)庫(kù)。隨后該研究構(gòu)建了ETA肽類拮抗劑的藥效團(tuán)模型，包括1個(gè)疏水特征、1個(gè)負(fù)電中心、1個(gè)芳環(huán)特征和5個(gè)排除體積。同時(shí)，利用同源模建方法構(gòu)建ETA的蛋白三維結(jié)構(gòu)模型，并將其用于分子對(duì)接研究，采用一致性打分評(píng)價(jià)化合物的ETA拮抗活性。通過(guò)構(gòu)建ETA肽類拮抗劑的藥效團(tuán)模型、同源模建的三維蛋白結(jié)構(gòu)以及分子對(duì)接模型，共預(yù)測(cè)獲得27條可能具有潛在ETA拮抗活性的寡肽分子。該文進(jìn)一步分析關(guān)鍵氨基酸GLN165，并結(jié)合文獻(xiàn)，說(shuō)明其對(duì)拮抗劑活性的重要性。計(jì)算機(jī)輔助水解方法可以高效地對(duì)已知序列結(jié)構(gòu)的中藥蛋白進(jìn)行模擬水解，結(jié)合分子模擬模型，能進(jìn)一步辨識(shí)中藥寡肽潛在的生物學(xué)活性。該研究為快速、高效開(kāi)展中藥來(lái)源的肽類物質(zhì)的活性機(jī)制研究提供了方法學(xué)依據(jù)。

[關(guān)鍵詞] 計(jì)算機(jī)輔助水解； 寡肽； 大豆； 內(nèi)皮素受體A； 藥效團(tuán)； 同源模建； 分子對(duì)接

[Abstract] Oligopeptides are one of the the key pharmaceutical effective constituents of traditional Chinese medicine（TCM）. Systematic study on composition and efficacy of TCM oligopeptides is essential for the analysis of material basis and mechanism of TCM. In this study， the potential anti-hypertensive oligopeptides from Glycine max and their endothelin receptor A （ETA） antagonistic activity were discovered and predicted based on in silico technologies.Main protein sequences of G. max were collected and oligopeptides were obtained using in silico gastrointestinal tract proteolysis. Then， the pharmacophore of ETA antagonistic peptides was constructed and included one hydrophobic feature， one ionizable negative feature， one ring aromatic feature and five excluded volumes. Meanwhile， three-dimensional structure of ETA was developed by homology modeling methods for further docking studies. According to docking analysis and consensus score， the key amino acid of GLN165 was identified for ETA antagonistic activity. And 27 oligopeptides from G. max were predicted as the potential ETA antagonists by pharmacophore and docking studies.In silico proteolysis could be used to analyze the protein sequences from TCM. According to combination of in silico proteolysis and molecular simulation， the biological activities of oligopeptides could be predicted rapidly based on the known TCM protein sequence. It might provide the methodology basis for rapidly and efficiently implementing the mechanism analysis of TCM oligopeptides.

[Key words] in silico proteolysis； oligopeptides； Glycine max； ETA； pharmacophore； homology modeling； docking

中藥蛋白是中藥中一類常見(jiàn)的活性物質(zhì)，具有明顯而特殊的藥理活性[1]。從兩千多年前阿膠的應(yīng)用，到水蛭素的制備和大豆多肽的上市，中藥蛋白和生物活性肽的應(yīng)用具有悠久的歷史，已有大量的藥性、藥理和臨床數(shù)據(jù)。而中藥蛋白具有物理化學(xué)性質(zhì)不穩(wěn)定、易水解和易氧化的特點(diǎn)，因此，基于傳統(tǒng)的中藥化學(xué)方法，通過(guò)對(duì)中藥蛋白提取分離酶解獲取中藥活性肽，往往具有分離困難、耗時(shí)費(fèi)力的特點(diǎn)。而計(jì)算機(jī)輔助水解是基于生物信息學(xué)方法，預(yù)測(cè)蛋白質(zhì)序列在特定蛋白酶或化學(xué)物質(zhì)在給定條件下的切割位點(diǎn)的技術(shù)，該技術(shù)具有快速高效可靠的寡肽預(yù)測(cè)能力[2]，已較為廣泛的應(yīng)用于食品工業(yè)研究[3-4]。

中藥大豆是豆科植物大豆Glycine max的成熟種子，現(xiàn)代研究表明大豆活性肽具有良好的降壓活性[5]。內(nèi)皮素受體A（endothelin receptor A，ETA）是治療高血壓和肺動(dòng)脈高壓的經(jīng)典靶標(biāo)，其第一代拮抗劑即為肽類化合物[6]。因此，基于大豆寡肽開(kāi)發(fā)ETA拮抗劑具有較高的研究?jī)r(jià)值。本研究擬利用計(jì)算機(jī)輔助水解技術(shù)，進(jìn)行中藥大豆的模擬消化水解，并結(jié)合基于藥效團(tuán)和分子對(duì)接的虛擬篩選技術(shù)，預(yù)測(cè)中藥大豆寡肽的ETA拮抗活性。

1 材料與方法

1.1 大豆蛋白的模擬水解 本研究在Uniport（http：//www.uniprot.org/）中檢索已解析的大豆種子儲(chǔ)存蛋白的多肽序列，共獲得5條大豆球蛋白（glycinin）和3條伴大豆球蛋白（conglycinin）序列，主要包括大豆球蛋白G，大豆球蛋白G1-G4，β-伴大豆球蛋白的α，α′和β鏈。其為大豆中的主要儲(chǔ)存蛋白，占分離蛋白的80%以上[7]。基于獲得的蛋白序列，本研究采用BIOPEP（http：//www.uwm.edu.pl/biochemia/index.php/pl/biopep）進(jìn)行計(jì)算機(jī)輔助模擬水解，以預(yù)測(cè)蛋白序列的切割位點(diǎn)。本研究采用模擬胃腸道的水解策略，選擇胃蛋白酶（pH 1.3）、胰蛋白酶和糜蛋白酶作為消化酶[8]。針對(duì)獲得的寡肽序列，利用Discovery Studio 4.0（DS）進(jìn)行寡肽三維結(jié)構(gòu)數(shù)據(jù)庫(kù)的構(gòu)建。

1.2 ETA肽類拮抗劑的藥效團(tuán)模型的構(gòu)建 為了更可靠的篩選ETA的肽類拮抗劑，本研究擬采用DS中的3D-QSAR pharmacophore（Hypogen）模塊，構(gòu)建ETA肽類拮抗劑的藥效團(tuán)模型。構(gòu)建藥效團(tuán)的訓(xùn)練集和測(cè)試集中的50個(gè)活性化合物均來(lái)源于文獻(xiàn)報(bào)道[9]，其包括直鏈肽和環(huán)肽2類化學(xué)結(jié)構(gòu)，其活性值分布5個(gè)數(shù)量級(jí)。

對(duì)訓(xùn)練集中20個(gè)化合物在CHARMm中進(jìn)行三維構(gòu)象生成和能量最小化處理。隨后，對(duì)其進(jìn)行藥效特征分析，發(fā)現(xiàn)活性化合物中包括的藥效特征有氫鍵受體（hydrogen bond acceptor，A）、親脂性氫鍵受體（hydrogen bond acceptor lipid，Ali）、氫鍵供體（hydrogen bond donor，D）、疏水基團(tuán)（hydrophobic region，H）、脂性疏水基團(tuán)（hydrophobic aliphatic region，Hal）、芳香疏水基團(tuán)（hydrophobic aromatic region，Har）、芳香環(huán)（ring aromatic，R）、負(fù)電基團(tuán)（ionizable negative group，N）。

基于特征分析，進(jìn)一步構(gòu)建藥效團(tuán)模型。首先，對(duì)訓(xùn)練集中活性化合物進(jìn)行構(gòu)象分析，采用Best模式，設(shè)置能量收斂值為20 kcal·mol-1，最大構(gòu)象數(shù)目255個(gè)。選擇生成的5個(gè)藥效特征為氫鍵受體、氫鍵供體、疏水基團(tuán)、芳香環(huán)和負(fù)電基團(tuán)，同時(shí)，模型中通過(guò)增加5個(gè)排除體積（maximum excluded volumes，Ev）以提高模型的特異性，并設(shè)置化合物的活性不確定度為1.5，利用Hypogen模塊構(gòu)建10個(gè)候選藥效團(tuán)模型。

針對(duì)生成的藥效團(tuán)模型需進(jìn)行進(jìn)一步的模型評(píng)價(jià)。利用剩余的30個(gè)活性化合物及從Binding Database（http：//www.bindingdb.org/bind/index.jsp）中隨機(jī)挑選的90個(gè)非活性化合物，以1∶3的比例生成測(cè)試集，并對(duì)測(cè)試集分子進(jìn)行同訓(xùn)練集分子相同的計(jì)算處理。隨后，利用生成的候選藥效團(tuán)模型與測(cè)試集化合物進(jìn)行匹配篩選，計(jì)算相應(yīng)的有效命中率（hit rate of active compounds，HRA）、辨識(shí)有效性指數(shù)（identify effective index，IEI）、綜合評(píng)價(jià)指數(shù)（comprehensive appraisal index，CAI）等指標(biāo)對(duì)藥效團(tuán)優(yōu)劣進(jìn)行評(píng)價(jià)[10]。3個(gè)外部評(píng)價(jià)指標(biāo)計(jì)算公式如下所示（1）～（3），其中TD代表測(cè)試集中化合物總數(shù)，TA代表測(cè)試集中活性化合物數(shù)目。Ha代表命中的活性化合物數(shù)目，Ht代表命中化合物總數(shù)。本文綜合利用藥效團(tuán)模型的多種評(píng)價(jià)參數(shù)，選出最優(yōu)的藥效團(tuán)模型。利用最優(yōu)的藥效團(tuán)模型進(jìn)行大豆寡肽庫(kù)的虛擬篩選，其中搜索模式為Best Search、匹配模式為柔性匹配，初步獲得潛在的活性寡肽。

HRA=Ha/TA×100%（1）

IEI=（Ha/Ht）/（TA/TD）（2）

CAI=IEI×HRA（3）

1.3 ETA的同源模建及分子對(duì)接 ETA（Uniprot/P25101）由427個(gè)氨基酸殘基組成，是典型的G蛋白偶聯(lián)受體（G protein-coupled receptor，GPCR）中A家族的7次跨膜蛋白。通過(guò)對(duì)已解析的同家族蛋白晶體進(jìn)行序列比對(duì)，分析其同源性和相似性的高低，選擇相似性較高的蛋白作為模板蛋白。本研究同源模建使用Modeller 9.13程序，采用多模板建模方法，生成5個(gè)候選模型進(jìn)行進(jìn)一步優(yōu)化。首先，基于逐級(jí)能量最小化方法，先約束所有重原子，優(yōu)化氫原子，再約束保守跨膜區(qū)，優(yōu)化側(cè)鏈及Loop區(qū)，最后優(yōu)化跨膜區(qū)。能量最小化過(guò)程采用最陡下降法和共軛梯度法進(jìn)行優(yōu)化。最后基于拉氏圖和ERRAT分析，利用SAVES服務(wù)器（http：//services.mbi.ucla.edu/SAVES/）對(duì)同源模建蛋白進(jìn)行評(píng)價(jià)。

結(jié)合對(duì)蛋白空腔的自動(dòng)識(shí)別，同時(shí)參考GPCR同家族活性口袋的位置，定義ETA蛋白的活性口袋。由于ETA活性口袋中可能存在多種結(jié)合模式[11]，為保證分子對(duì)接結(jié)果的準(zhǔn)確性，本研究選擇較大的活性口袋半徑。本研究利用DS中的Libdock模塊進(jìn)行分子對(duì)接。首先，以陽(yáng)性藥作為測(cè)試，分別對(duì)接進(jìn)入活性口袋中，并同時(shí)計(jì)算5種打分函數(shù)，采用一致性打分的方法，評(píng)價(jià)對(duì)接模型的可靠性，并分析ETA肽類拮抗劑可能存在的作用模式。隨后，利用上述的對(duì)接模型進(jìn)一步篩選藥效團(tuán)的命中結(jié)果，篩選潛在的ETA肽類拮抗劑。

2 結(jié)果與討論

2.1 大豆寡肽庫(kù)的構(gòu)建 通過(guò)模擬水解，本研究共獲得449條多肽，其中各段寡肽的數(shù)目分布見(jiàn)圖1，由于2～6肽可能具有較好的相對(duì)分子質(zhì)量與成藥性[2]，所以選擇2～6肽共299條作為本研究的主要研究對(duì)象。通過(guò)對(duì)模擬水解結(jié)果與已知的活性肽數(shù)據(jù)庫(kù)BIOPEP比對(duì)，共發(fā)現(xiàn)已知活性肽32條，其主要的生物學(xué)活性包括ACE抑制作用和抗氧化作用等，其結(jié)果見(jiàn)表1，表明大豆寡肽具有較為廣泛的藥理活性。

2.2 ETA肽類拮抗劑藥效團(tuán)模型的建立與評(píng)價(jià) 基于Hypogen的藥效團(tuán)建模方法主要包括2類內(nèi)部評(píng)價(jià)參數(shù)，即訓(xùn)練集活性化合物的預(yù)測(cè)結(jié)果與真實(shí)活性的相關(guān)系數(shù)Corrtrain，以及費(fèi)用函數(shù)（Cost）對(duì)所構(gòu)建的藥效團(tuán)模型的顯著性評(píng)價(jià)。本研究構(gòu)建的10個(gè)藥效團(tuán)內(nèi)部評(píng)價(jià)參數(shù)見(jiàn)表2。模型的Null Cost為425.75，ΔCost為Null Cost與Total Cost的差值，ΔCost大于60說(shuō)明模型有90%以上的顯著性。

利用測(cè)試集對(duì)10個(gè)候選藥效團(tuán)進(jìn)行評(píng)價(jià)，計(jì)算相應(yīng)外部評(píng)價(jià)參數(shù)，包括藥效團(tuán)的活性化合物命中率HRA、辨識(shí)有效性指數(shù)IEI和綜合評(píng)價(jià)指數(shù)CAI以及測(cè)試集中活性化合物的預(yù)測(cè)值與真實(shí)值之間的相關(guān)系數(shù)Corrtest，其結(jié)果見(jiàn)表2。本研究生成的藥效團(tuán)模型主要包括2類，即HNR型和DHR型。HNR型主要包括1個(gè)疏水特征，1個(gè)負(fù)電中心和1個(gè)芳環(huán)基團(tuán)，而DHR型主要包括1個(gè)氫鍵受體，1個(gè)疏水特征和1個(gè)芳香基團(tuán)。通過(guò)對(duì)比其CAI值差異，發(fā)現(xiàn)HNR型具有更好的測(cè)試集篩選能力。進(jìn)一步分析藥效團(tuán)對(duì)測(cè)試集的活性預(yù)測(cè)能力，發(fā)現(xiàn)1號(hào)和5號(hào)藥效團(tuán)具有相當(dāng)?shù)幕钚曰衔锉孀R(shí)能力。通過(guò)比較2個(gè)藥效團(tuán)的空間結(jié)構(gòu)差異，發(fā)現(xiàn)2個(gè)藥效團(tuán)中HNR 3個(gè)特征之間的空間距離及角度具有相似性，同時(shí)，5號(hào)藥效團(tuán)的D特征與其排除體積相重疊，因此最終選擇1號(hào)藥效團(tuán)作為最優(yōu)藥效團(tuán)對(duì)大豆寡肽庫(kù)進(jìn)行虛擬篩選，結(jié)果共獲得64條寡肽。其中，VVF與藥效團(tuán)具有較好的匹配，其匹配值為0.85，見(jiàn)圖2。

2.3 ETA同源模建及分子對(duì)接研究 通過(guò)同源比對(duì)，選擇3個(gè)蛋白晶體結(jié)構(gòu)作為模板進(jìn)行同源模建，包括Delta Opioid Receptor（PDB：4N6H_A），Neurotensin Receptor 1（PDB：4BV0_A）和orphanin FQ receptor（PDB：4EA3_A）。3個(gè)受體均來(lái)源于GPCR的A家族，且其與ETA的序列相似度 （similarity）分別為40.0%，39.8%，38.8%，其一致性（identity）分別為26%，24%，25%。利用獲得的3個(gè)模板進(jìn)行多模板建模，通過(guò)能量最小化處理，利用拉氏圖和ERRAT評(píng)價(jià)模型的優(yōu)劣，最終模型和評(píng)價(jià)參數(shù)見(jiàn)圖3。其中，拉氏圖的最適殘基量為96.96%，ERRAT的值為91.41，顯示模建結(jié)構(gòu)具有良好的立體化學(xué)參數(shù)和空間結(jié)構(gòu)。

ETA分子對(duì)接結(jié)果：本研究通過(guò)分析ETA的空腔結(jié)構(gòu)，自動(dòng)識(shí)別了8個(gè)活性口袋，通過(guò)文獻(xiàn)中關(guān)鍵氨基酸的描述和參考GPCR同家族的位點(diǎn)情況，最終確定對(duì)接口袋，其半徑為15 。先將訓(xùn)練集中20個(gè)陽(yáng)性寡肽與活性口袋對(duì)接，20個(gè)化合物均能對(duì)接進(jìn)入活性口袋中，進(jìn)一步計(jì)算Jain，-PLP2，-PMF，-PMF04和-PLP1 5種打分函數(shù)，綜合已有的LibDockScore打分，進(jìn)行一致性打分評(píng)價(jià)[12]。其中，一致性打分大于3的化合物為95.0%，說(shuō)明了對(duì)接模型和一致性打分的可靠性[13]。隨后，分析ETA與陽(yáng)性藥發(fā)生相互作用的高頻氨基酸，產(chǎn)生氫鍵作用的氨基酸主要包括：LYS329，ARG326，ASP351和GLN165，產(chǎn)生疏水作用的氨基酸主要有LEU347，LYS140，LEU141，ARG145，ARG326，LYS329，與文獻(xiàn)報(bào)道基本一致[14-15]。

進(jìn)一步將測(cè)試集中30個(gè)活性分子與活性口袋對(duì)接，一致性打分大于3的化合物為96.7%，說(shuō)明對(duì)接模型和一致性打分評(píng)價(jià)對(duì)化合物具有較好的預(yù)測(cè)能力。發(fā)生氫鍵作用的關(guān)鍵氨基酸主要有GLN165，LYS329和ASP351，產(chǎn)生疏水作用的高頻氨基酸有ARG326，LEU322和LYS329等，與文獻(xiàn)報(bào)道及訓(xùn)練集化合物基本一致。其中高活性化合物Com1385與ETA之間的相互作用見(jiàn)圖4。

利用構(gòu)建的分子對(duì)接模型和一致性打分方法，篩選被藥效團(tuán)命中的64條寡肽。27條寡肽的一致性打分為6，提示大豆虛擬水解獲得的寡肽可能具有較好的ETA拮抗活性。其中，VVF能與關(guān)鍵氨基酸具有良好的匹配，VVF中苯丙氨酸的羧基與藥效團(tuán)的負(fù)電中心匹配，而其羰基氧與GLN165發(fā)生氫鍵相互作用，而其羥基氫能與ASP351發(fā)生氫鍵相互作用，因此，VVF的苯丙氨酸羧基可能是發(fā)揮活性的重要基團(tuán)[16]，其結(jié)果見(jiàn)圖5。

3 結(jié)論

本研究基于計(jì)算機(jī)輔助水解方法，對(duì)大豆中的主要蛋白進(jìn)行模擬消化水解，獲得2～6肽的寡肽基本結(jié)構(gòu)共299條。在此基礎(chǔ)上，通過(guò)構(gòu)建ETA肽類拮抗劑的藥效團(tuán)，以及ETA蛋白同源模建與分子對(duì)接，篩選并預(yù)測(cè)大豆寡肽的ETA拮抗活性，共獲得27條可能具有ETA拮抗活性的寡肽分子，其結(jié)果有待進(jìn)一步的活性驗(yàn)證。本文通過(guò)基于生物信息學(xué)和計(jì)算機(jī)輔助藥物設(shè)計(jì)的方法，從大豆蛋白中篩選潛在的ETA拮抗活性寡肽，為探討中藥蛋白發(fā)揮藥效的作用機(jī)制探討提供了一條新的研究方法，為開(kāi)展中藥蛋白的酶解分離或藥理驗(yàn)證提供指導(dǎo)。

[參考文獻(xiàn)]

[1] 鄒吉利，徐南平.中藥活性多肽研究進(jìn)展[J].湖北中醫(yī)藥大學(xué)學(xué)報(bào)， 2012， 14（4）：66.

[2] 鄒平. 基于生物信息學(xué)與QSAR及分子對(duì)接的菜粕活性肽篩選及活性研究[D]. 杭州：浙江大學(xué)， 2014.

[3] Dziuba B， Dziuba M. New milk protein-derived peptides with potential antimicrobial activity： an approach based on bioinformatic studies[J]. Int J Mol Scis， 2014， 15（8）：14531.

[4] 孫國(guó)威，樂(lè)國(guó)偉，施用暉.模擬酶解大豆7S、11S蛋白及其抗氧化活性的研究[J].食品工業(yè)科技， 2010（7）：101.

[5] 于婷婷，韓飛，陳光.大豆降壓肽研究進(jìn)展[J].糧油食品科技， 2008， 16（2）：27.

[6] 楊菁.內(nèi)皮素受體拮抗劑的設(shè)計(jì)與合成[D].北京：中國(guó)人民解放軍軍事醫(yī)學(xué)科學(xué)院，2008.

[7] 田琨，管娟，邵正中，等.大豆分離蛋白結(jié)構(gòu)與性能[J].化學(xué)進(jìn)展， 2008， 20（4）：565.

[8] Lafarga T， O′Connor P， Hayes M. Identification of novel dipeptidyl peptidase-IV and angiotensin-I-converting enzyme inhibitory peptides from meat proteins using in silico analysis[J]. Peptides， 2014， 59（9）：53.

[9] Fukami T， Yamakawa T， Niiyama K， et al. Synthesis and structure-activity relationships of 2-substituted D-tryptophan-containing peptidic endothelin receptor antagonists： importance of the C-2 substituent of the D-tryptophan residue for endothelin A and B receptor subtype selectivity[J]. J Med Chem， 1996， 39（12）：2313.

[10] Jiang L， Zhang X， Chen X， et al. Virtual screening and molecular dynamics study of potential negative allosteric modulators of mGluR1 from Chinese herbs[J]. Molecules， 2015， 20：12769.

[11] Mey J G R D， Compeer M G， Pieter L， et al. ETA -receptor antagonists or allosteric modulators[J]. Trends Pharmacol Sci， 2011， 32（6）：345.

[12] 劉維國(guó).抗藥物依賴功效藥效團(tuán)模型的構(gòu)建及其虛擬篩選[D].北京：中央民族大學(xué)， 2013.

[13] Krovat E M， Langer T. Impact of scoring functions on enrichment in docking-based virtual screening： an application study on renin inhibitors[J]. J Chem Inf Comp Sci， 2004， 44（3）： 1123.

[14] Breu V， Hashido K， Broger C， et al. Separable binding sites for the natural agonist endothelin-1 and the non-peptide antagonist bosentan on human endothelin-A receptors[J]. Eur J Biochem， 1995， 231（1）：266.

[15] Jens L， Alexander O， Michael B， et al. Structural determinants for selective recognition of peptide ligands for endothelin receptor subtypes ETA and ETB[J]. J Pept Sci， 2009， 15（7）：479.

[16] Khuraijam Dhanachandra S， Karthikeyan M. Molecular modeling， quantum polarized ligand docking and structure-based 3D-QSAR analysis of the imidazole series as dual AT1 and ETA receptor antagonists[J]. Acta Pharm Sin， 2013， 34（12）：1592.