貝類雌激素受體的同源建模與分子對接研究?

劉力銣， 苗晶晶， 趙安冉， 潘魯青

(海水養(yǎng)殖教育部重點實驗室(中國海洋大學)，山東 青島 266003)

在過去的一個世紀中，有超過80 000種化學物質被引入到環(huán)境中[1]。其中一些物質如壬基酚、多環(huán)芳烴或溴系阻燃劑等具有模擬和改變?nèi)祟惢蚱渌吧锛に睾铣膳c代謝的效應，被稱為環(huán)境內(nèi)分泌干擾物(Endocrine disrupting chemicals, EDCs)[2-6]。EDCs暴露可對水生生物的內(nèi)分泌系統(tǒng)產(chǎn)生有害影響，如性畸變、性別比例改變、雄性個體生殖障礙和發(fā)育異常、個體性早熟和腫瘤發(fā)生等[6]。因此，及早發(fā)現(xiàn)EDCs對生殖內(nèi)分泌系統(tǒng)的潛在干擾效應至關重要?？紤]到篩選與受體存在結合可能性的多種EDCs所耗費的時間、費用與工作量等因素，人們對開發(fā)預測化合物與受體親和力計算方法的興趣日益增加[7-8]。目前，同源建模與分子對接已逐漸成為基于化學結構的計算模擬分析中不可或缺的一部分。同源建模技術是使用確定的X衍射晶體學數(shù)據(jù)來預測具有相似氨基酸序列的另一種蛋白質構象的方法[9]，通常涉及四個步驟，分別是序列比對、模板選擇、模型構建和模型評估。分子對接能夠通過計算結合親和力預測化合物的生物活性態(tài)，并確定結合構象的優(yōu)先級[10]。

迄今為止，圍繞EDCs對脊椎動物的內(nèi)分泌干擾機制已進行了許多研究，特別是卵生脊椎動物[11]。由于現(xiàn)有水處理系統(tǒng)無法將EDCs從廢水中完全去除[12-13]，因此海洋及淡水環(huán)境中仍然存在著大量EDCs。盡管在實驗室和田間研究中已經(jīng)以魚類為研究對象，圍繞各種化學物質或混合物的雌激素活性進行了廣泛的評估[14]，但EDCs對貝類的分子靶標和內(nèi)分泌干擾效應還知之甚少，仍然需要更進一步的探究[15]。如今已經(jīng)確定，水生脊椎動物斑馬魚的基因組編碼三種雌激素受體，包括ERa、ER2a和ER2b。其中，斑馬魚的ERa在脊椎動物間具有最大的同源性，而ER2a和ER2b可能由于基因組復制而引起了幾個關鍵氨基酸變化[16]。有大量的實驗室和現(xiàn)場證據(jù)表明，貝類ER對雌激素的反應具有與魚類相當?shù)拿舾行訹17]，但在紫貽貝(Mytilusedulis和M.galloprovincialis)[18]、櫛孔扇貝(Chlamysfarreri)[19]、團聚牡蠣(Saccostreaglomerata)[20]與長牡蠣(Crassostreagigas)[21]中都只鑒定到一種ER亞型，其與魚類ERa序列呈較高同源性。因此，本研究利用同源建模技術構建并評估了長牡蠣和菲律賓蛤仔(Ruditapesphilippinarum)雌激素受體配體結合區(qū)(Estrogen receptor-ligand binding domain, ER-LBD)的同源模型，并同時構建了模式生物斑馬魚(Daniorerio)ERa模型作為參照。針對典型的脊椎動物內(nèi)分泌干擾物苯并芘(B[a]P)、對壬基酚(4-NP)和四溴雙酚A(TBBPA)，利用分子對接模擬技術探究了氨基酸殘基與這些化合物的相互作用，解析了ER雌激素活性相關的活性位點特征和與配體的結合親和力，比較了貝類ER與魚類ERa的配體親和特征，為研究貝類ER結構功能和與配體結合的作用機制提供參考，同時為探究EDCs對貝類ER內(nèi)分泌干擾機制奠定理論基礎和技術支撐。

1 材料與方法

1.1 序列比對與同源建模

通過NCBI網(wǎng)站查詢靶基因LBD的蛋白氨基酸序列(見表1)，并利用Editseq(DNAStar, Madison, WI, USA)保存為seq格式，隨后利用Modeller 9.22(Andrej Sali Lab, UCSF,USA)軟件進行靶序列的同源建模。從Modeller序列數(shù)據(jù)庫與PDB蛋白數(shù)據(jù)庫(http://www.rcsb.org/)中讀取模板reads，利用salign命令進行比對與計算質量得分。為了在輸出文件中獲得有更多的模板檢索結果，將參數(shù)max_aln_evalue修改為1。對靶蛋白與模板序列結構的同源性檢查后對齊，利用離散優(yōu)化蛋白能量(Discrete optimized protein energy, DOPE)評估新的模型坐標，選擇打分絕對值更高的模型作為更為準確的構象模板。隨后依據(jù)Modeller軟件教程生成并提取評分最高的10個三維模型。選擇軟件自帶的打分函數(shù)中DOPE分數(shù)最低且GA341分數(shù)最高的模型作為最優(yōu)模型。將模板模型與同源模型原子位置對齊后，利用PyMOL 2.3(PyMOL Molecular Graphics System, Schr?dinger)軟件計算其均方根偏差(Root mean square deviations, RMSDs)以確定二者結構的差異。將最大循環(huán)數(shù)設定為10，浮點數(shù)設置為2.0。

表1 同源建模研究所用的靶蛋白氨基酸序列Table 1 Amino acid sequences of target proteins used in homology modeling research

1.2 模型評價

使用PROCHECK、Verify3D和ERRAT程序對同源模型的質量進行了全面評價。PROCHECK(http://services.mbi.ucla.edu/PROCHECK/)常用于檢測蛋白質三維結構中氨基酸殘基的空間構象合理性，通過檢查蛋白質結構的殘基-殘基立體化學質量來進行拉氏構象分析[22]。ERRAT(https://servicesn.mbi.ucla.edu/ERRAT/)通過比較分子間非氫鍵性質以確定模型的置信度和整體質量[23]。Verify3D(https://servicesn.mbi.ucla.edu/Verify3D/)通過檢查建構模型(3D)與其本身氨基酸序列(1D)的兼容性以評估整體質量[24-25]。

1.3 對接配體制備

對接的配體分子為內(nèi)分泌干擾物苯并芘(B[a]P)、四溴雙酚A(TBBPA)、對壬基酚(4-NP)，所用配體結構如表2所示。利用Chemdraw 16.0(Cambridge Soft, Cambridge, USA)軟件畫出上述配體分子的平面結構，經(jīng)由Chem3D(Cambridge Soft, Cambridge, USA)轉化為三級結構。隨后利用Open Babel(http://openbabel.org/wiki/Main_Page)將其轉化為pdb格式，用于后續(xù)對接研究。

表2 分子對接研究所用配體結構

1.4分子對接

使用AutoDock 4.0(The Scripps Research Institute, USA)程序對ERs-LBD與配體之間的相互作用潛力進行了模擬分析。利用AutoGrid 4.0(The Scripps Research Institute, USA)軟件進行能量格點計算，將靶點蛋白的活性位點設置為中心坐標，使Grid box網(wǎng)格區(qū)域將可能活性區(qū)域完全覆蓋。然后采用Lamarc-kian GA 4.2進行區(qū)域搜索，初始種群數(shù)設置為150，能量評定最大次數(shù)為2 500 000，每個小分子配體都設定得到20個構象，其他參數(shù)選擇默認數(shù)值。為了進一步探究受體與配體的結合，運行AutoDock進行半柔性對接并對結果進行聚類分析。根據(jù)單個構象的相互作用能信息，選擇最優(yōu)聚類中的有利構象，并導入至Autodock/Vina plugin for PyMOL插件程序。基于最低的結合能(kJ/mol)確定三個模型最穩(wěn)定的配體結合模式。使用PyMOL插件程序查看了生成的對接形式與關鍵作用氨基酸殘基，并分析配置/分數(shù)關系。最后，使用PyMOL軟件進行了圖像優(yōu)化。

2 結果

2.1 同源建模分析與模型評價

分別將斑馬魚、長牡蠣與菲律賓蛤仔ER-LBD蛋白序列與Modeller數(shù)據(jù)庫和PDB數(shù)據(jù)庫中進行相似性搜索后，得到一系列與靶基因同源性較高的蛋白質X衍射晶體結構。除去相近受體與突變體的晶體結構(如過氧化物酶體增殖劑激活受體PPAR，維甲酸X受體RXR等)后，比較各模板晶體結構的同源性與E值以篩選合適模板。該研究涉及靶基因與內(nèi)分泌干擾物的結合，因此優(yōu)先選擇與化學品結合的復合結構?；谏鲜龊Y選條件，選擇1A52(人ERa-雌二醇復合物)、5KRH(人ERa-16-亞芐基雌酮復合物)、3UU7(人ERa-雙酚A復合物)和4N1Y(長牡蠣ER-LBD)等相似性極顯著(E值=0)的4個PDB序列作為模板進行模型建構。建模后查看了同源模型與模板模型的差異，發(fā)現(xiàn)斑馬魚同源模型與同樣為脊椎動物的人類ERa模板序列相似性較高，蛋白結構的RMSD較低；長牡蠣與菲律賓蛤仔同源模型與同樣為貝類的長牡蠣模板序列相似性較高，蛋白結構的RMSD較低。證明了本研究所構建的三個同源模型的原子位置相對于未知的實際結構誤差較小(見表3)。理論上，規(guī)范的二級結構(例如形成配體結合位點的(螺旋)應該是同源模型中預測最準確的部分[26]。因此本研究所構建ERs-LBD同源模型計算的RMSD值可能是與實際結構形成誤差的上限值。

本研究構建的斑馬魚ER-LBD三維結構包含氨基酸殘基245個，長牡蠣三維結構包含227個氨基酸殘基，菲律賓蛤仔三維結構包含226個氨基酸。分別利用PROCHECK、ERRAT和Verify3D軟件程序對同源模型的可信度進行了評價(見表4)。PROCHECK產(chǎn)生的拉氏構象分析結果表明，斑馬魚、長牡蠣與菲律賓蛤仔ER-LBD三維結構在核心區(qū)域的氨基酸殘基所占比率分別為95.5%、96.7%和94.8%，證明本研究所構建的三個同源模型是質量較高的模型，具有良好的立體化學特征。ERRAT評估了模型中不同類型原子的隨機分布情況，結果表明三個模型的ERRAT值分別為98.3%、99.5%和100.0%，均高于95%，說明模型具有良好的高分辨率結構和合理的原子分布。Verify3D評價表明三個同源模型均沒有產(chǎn)生構象錯誤，且由超過83.2%、100.0%和99.6%的氨基酸殘基得分>0.2判斷上述三個模型的三級結構得到了合理折疊。

表3 斑馬魚、長牡蠣與菲律賓蛤仔ER-LBD同源模型與模板蛋白的結構差異Table 3 Structural differences between templates and homology models of zebrafish, oyster and clam ER-LBD

表4 斑馬魚、長牡蠣與菲律賓蛤仔ER-LBD模型質量檢查的參數(shù)結果Table 4 Parameters results for the model quality of zebrafish, oyster and clam ER-LBD /%

2.2 對接能量分析

為了進一步從分子水平上探究ERs氨基酸殘基與配體的結合敏感性與結合模式，我們評估了上述ERs-LBD模型與三種EDCs(B[a]P、TBBPA和4-NP)的對接結果。在綜合考慮對接打分和模型結構性的情況下，選擇的最優(yōu)對接構象與相應能量如表5所示。斑馬魚ER-LBD活性中心為x=12.538，y=6.108，z=-7.844(間距=0.508)，長牡蠣ER-LBD活性中心為x=16.779，y=-13.759，z=4.022(間距=0.447)，菲律賓蛤仔ER-LBD活性中心為x=16.442，y=-15.143，z=3.807(間距=0.447)。其中，結合自由能是分子間能、分子內(nèi)能與扭轉能之和，再減去未結合能量；對接能為分子間能與分子內(nèi)能的總和[27]。

表5 斑馬魚、長牡蠣與菲律賓蛤仔ER-LBD三維結構與配體分子的對接結果更改圖層

2.3 對接構象分析

使用Autodock軟件預測了魚類和兩種貝類ER-LBD與EDCs的相互作用，三維對接模型顯示了B[a]P(見圖1)、4-NP(見圖2)和TBBPA(見圖3)與結合位點殘基產(chǎn)生疏水相互作用的口袋。配體以球棒模型顯示，并被云表面包圍。當配體與受體接觸時，接觸的配體分子表面呈灰色；配體與受體之間沒有接觸的區(qū)域呈透明或白色。受體與配體間產(chǎn)生氫鍵的部分以空間網(wǎng)狀填充球形式表示，與配體相互作用的主要氨基酸殘基以紅色字體標明。對接結果表明， 三個ERs-LBD同源模型與EDCs形成相互作用網(wǎng)絡的主要氨基酸殘基基本相似。

在設置比例因子為1時，斑馬魚ER-LBD同源模型與B[a]P對接構象中參與的氨基酸殘基共計11個(見表5)。對接時未形成氫鍵，其結合能、結合自由能、

圖1 斑馬魚、長牡蠣與菲律賓蛤仔ER-LBD三維結構與B[a]P的對接結果

對接能都為-35.37 kJ/mol。與4-NP對接的構象中參與的氨基酸殘基為11個，對接時的結合能和結合自由能都為-25.41 kJ/mol，對接能為-39.89 kJ/mol。其中，Glu321和Arg362分別與4-NP羥基的H和O原子形成氫鍵，氫鍵的距離為2.1 ?和1.8 ?。與TBBPA對接的構象中參與的氨基酸殘基同樣為11個，對接時的結合能為-31.26 kJ/mol，結合自由能為-31.31 kJ/mol，對接能為-40.10 kJ/mol。其中，Glu321和Val295分別與TBBPA羥基的H和O原子形成氫鍵，氫鍵的距離為2.2 ?和2.1 ?。

設置比例因子為1時，長牡蠣ER-LBD同源模型與B[a]P對接構象中參與的氨基酸殘基共計9個(見表5)。其中對接時未形成氫鍵，其結合能、結合自由能、對接能都為-37.88 kJ/mol。與4-NP對接的構象中參與的氨基酸殘基有15個，對接時的結合能和結合自由能都為-27.29 kJ/mol，對接能為-40.77 kJ/mol。Glu290和Arg331分別與4-NP羥基的H和O原子形成氫鍵，氫鍵的距離為2.2 ?和2.0 ?。與TBBPA對接的構象中參與的氨基酸殘基為13個，對接時的結合能為-33.95 kJ/mol，結合自由能為-33.98 kJ/mol，對接能為-42.86 kJ/mol。其中，Glu290與TBBPA羥基的H形成氫鍵，氫鍵距離為2.2 ?。

(形成氫鍵的氨基酸殘基以紅色字體加方框標明。The amino acid residues forming hydrogen bonds are marked in red and boxed.)

(形成氫鍵的氨基酸殘基以紅色字體加方框標明。The amino acid residues forming hydrogen bonds are marked in red and boxed.)

同樣設置比例因子為1時，菲律賓蛤仔ER-LBD同源模型與B[a]P對接構象中參與的氨基酸殘基共計10個(見表5)。對接時未形成氫鍵，其結合能、結合自由能、對接能都為-38.09 kJ/mol。與4-NP對接的構象中參與的氨基酸殘基有12個，對接時的結合能為-28.80 kJ/mol，結合自由能為-28.84 kJ/mol，對接能為-42.15 kJ/mol。其中，Glu242和Arg283分別與4-NP羥基的H和O原子形成氫鍵，氫鍵的距離為1.9 ?和2.2 ?。與TBBPA對接的構象中參與的氨基酸殘基同樣為12個，對接時的結合能和結合自由能都為-31.23 kJ/mol，對接能為-39.72 kJ/mol。其中，Glu242與TBBPA羥基的H原子形成氫鍵，氫鍵距離為2.2 ?。

基于預測的配體效率與結合能，兩種貝類ERs-LBD同源模型與EDCs結合的能力強于斑馬魚(見表5)。三個物種同源模型在與EDCs對接中大多形成氫鍵以激活最大的結合活性，除了B[a]P的結合能來自于分子間能外，其他種類EDCs(4-NP和TBBPA)的結合能由分子間能、分子內(nèi)能、扭轉能、未結合能和靜電勢能共同參與?；陬A測的結合能，EDCs與三個同源模型ERs-LBD的能量順序為B[a]P > TBBPA > 4-NP(見表5)。作為線性碳鏈，4-NP(11.22 kJ/mol)和TBBPA(4.98 kJ/mol)在受體的配體活性腔中常常需要消耗更多的能量進行構象變化。在不考慮扭轉能差異下，EDCs與三個同源模型ERs-LBD的對接能順序為TBBPA > 4-NP > B[a]P，除了菲律賓蛤仔同源模型表現(xiàn)為4-NP > TBBPA > B[a]P。

3 討論

研究人員普遍認為，脊椎動物ER蛋白的配體結合結構域(Ligand binding domain, LBD)是一個由激素控制的變構化開關[28-29]。當沒有配體存在時，LBDs通常處于一種不活躍的狀態(tài)；而當ERs與配體結合后，LBD構象將被重構為活性構象。在之前的一些研究中，曾采用計算機模擬在分子水平上觀察配體和受體之間的相互作用[30-31]。例如Celik等發(fā)現(xiàn)，多氯聯(lián)苯(PCBs)、增塑劑和農(nóng)藥都可以和人類ERa結合，且受體與這些類固醇結合位點的兩個親水端中至少一端相互作用[30]。另一方面，人類ERa與天然雌激素E2的復合物晶體結構表明，配體的羥基與Glu353、Arg394和受體結合位點中的保守水分子形成氫鍵以穩(wěn)定結合構象，并且與他莫昔芬和雷洛昔芬等EDCs絡合的晶體結構也顯示出相同的結合模式[32]。這些結果表明EDCs的羥基可能是與脊椎動物ERs緊密結合并進行轉錄激活的相關結構特征。與脊椎動物相比，關于EDCs對無脊椎動物(特別是軟體動物)的生殖內(nèi)分泌干擾作用仍存在著較大的爭議。

Bridgham等對長牡蠣ER的x-ray晶體結構的分析中指出，長牡蠣ER-LBD的結合空腔中沒有可容乃配體的電子密度，配體可能與類固醇受體結合的內(nèi)腔被幾條龐大的疏水性側鏈(例如F425和F525)占據(jù)[29]。因此，研究學者們認為長牡蠣ER-LBD結合空腔的總體積((168±8) ?3)比人ERa-LBD相應空腔(402 ?3)小得多，不足以容納17β-雌二醇(245 ?3)或其他EDCs。我們的對接結果表明，當以相同的標準(比例因子為1)的情況下查看配體-受體周圍的關鍵氨基酸殘基時，斑馬魚ER-LBD周圍呈8個氨基酸殘基，長牡蠣模型周圍呈11個氨基酸殘基，菲律賓蛤仔模型周圍呈12個氨基酸殘基。這表明相比于魚類的配體結合腔，長牡蠣和菲律賓蛤仔ER-LBD配體結合口袋確實小得多。但是，在隨后的分子對接研究中，我們發(fā)現(xiàn)長牡蠣和菲律賓蛤仔ER-LBD都能夠與EDCs結合形成相互作用網(wǎng)絡，其組成的主要氨基酸殘基與斑馬魚基本一致。除了B[a]P不能與ERs-LBD形成氫鍵，4-NP和TBBPA都能夠與斑馬魚、長牡蠣和菲律賓蛤仔的關鍵氨基酸殘基結合。通過對關鍵氨基酸殘基的多序列比對，我們發(fā)現(xiàn)EDCs的結合活性位點主要為谷氨酸Glu(斑馬魚：Glu321，牡蠣：Glu290，菲律賓蛤仔：Glu242)和精氨酸Arg(斑馬魚：Arg362，牡蠣：Arg331，菲律賓蛤仔：Arg283)。這表明上述兩種氨基酸殘基在與配體結合過程中，可能具有關鍵的功能。

自由結合能越低，一般表明配體與受體結合作用越穩(wěn)定[27]。本研究通過比較EDCs和三個ERs-LBD模型的自由結合能，探究了配體與受體的結合能力與穩(wěn)定程度。與EDCs模擬對接的結合能結果表明，B[a]P與ERs的親和力顯著高于TBBPA，而4-NP具有較低的結合能力(見表5)。我們注意到，B[a]P在與三個物種同源模型的結合中雖然不形成氫鍵，但由于多環(huán)類芳烴結構的特殊性，B[a]P在進入配體口袋中時不需要消耗扭轉能。而4-NP和TBBPA往往具有較強的分子內(nèi)能和分子間能，但改變構象時消耗的扭轉能降低了配體的結合效率。相比較于TBBPA(4.98 kJ/mol)，長線性碳鏈4-NP(11.22 kJ/mol)往往需要更高程度的扭轉以便在最低能量狀態(tài)下進入結合腔。因此根據(jù)對接的相互作用能與氫鍵形成的數(shù)量，理論上TBBPA應具有較之4-NP更強的結合潛力。根據(jù)Autodock的軟件說明，消耗的扭轉能在各種構象中不會發(fā)生改變[27]。為了分析EDCs與受體的結合能力，利用對接能分析EDCs與同源模型結合時消耗的能量。EDCs與三個同源模型ERs-LBD的對接能順序為TBBPA>4-NP>B[a]P，表明TBBPA的結合能力強于4-NP和B[a]P。這一推測在其他研究中得到了驗證：全球水質研究聯(lián)盟(Global water research coalition, GWRC)統(tǒng)計的環(huán)境雌激素活性檢測工具中發(fā)現(xiàn)，利用化學激活的螢光素酶基因表達分析法、報告基因檢測法和E-Screen 細胞增殖試驗法等方法，都能夠檢測到TBBPA對于人類ERa的雌激素活性顯著高于4-NP[33]。

4 結論

本研究利用同源建模與分子對接模擬技術探究了斑馬魚、長牡蠣與菲律賓蛤仔ER與三種典型EDCs的相互作用與結合模式，實驗結論如下：

(1)在激活構象中，EDCs作為配體分子存在于兩種貝類與斑馬魚ER的活性空腔中。相比于魚類的配體結合腔，長牡蠣和菲律賓蛤仔ER-LBD配體結合口袋要小得多。因此，在于相同的EDCs結合時，貝類ER往往需要消耗更多的能量以容納與結合EDCs。三個物種與EDCs結合的穩(wěn)定性順序為斑馬魚>長牡蠣>菲律賓蛤仔。

(2)形成三個同源模型配體結合口袋的氨基酸殘基高度保守，其中谷氨酸Glu(斑馬魚：Glu321，牡蠣：Glu290，菲律賓蛤仔：Glu242)和精氨酸Arg(斑馬魚：Arg362，牡蠣：Arg331，菲律賓蛤仔：Arg283)是與EDCs結合的主要活性位點，可能行使著與配體結合的關鍵的功能。

(3)EDCs與三個同源模型的對接能順序為TBBPA>4-NP>B[a]P。其中，TBBPA和4-NP能夠與關鍵氨基酸殘基形成氫鍵，從而激活最大的雌激素活性。而由于多環(huán)類芳烴結構的特殊性，B[a]P與同源模型結合將不會形成氫鍵。

本研究在分子水平上的對接模擬能夠為貝類ER與EDCs的作用位點和結合機制提供參考，同時為EDCs的優(yōu)選和親和程度的研究奠定理論基礎。但不可避免的是，分子對接模擬與真實結合情況必然存在一定的差異。因此，仍需進一步的實驗探究真實環(huán)境中貝類ER與天然雌激素和EDCs的親和情況與結合程度、實際結合的中間態(tài)與能量釋放，通過對各種情況的綜合分析以探究EDCs對貝類的內(nèi)分泌干擾機制并完善EDCs毒性評價體系。