雙酚A類似物雌激素干擾效應的定量結構-活性關系模型

楊先海，劉會會，楊倩,3，劉濟寧,

1. 環境保護部南京環境科學研究所，南京 210042 2. 南京理工大學 環境與生物工程學院江蘇省化工污染控制與資源化高校重點實驗室，南京 210094 3. 南京工業大學生物與制藥工程學院，南京 210009

?

雙酚A類似物雌激素干擾效應的定量結構-活性關系模型

楊先海1，劉會會2，楊倩1,3，劉濟寧1,

1. 環境保護部南京環境科學研究所，南京 210042 2. 南京理工大學 環境與生物工程學院江蘇省化工污染控制與資源化高校重點實驗室，南京 210094 3. 南京工業大學生物與制藥工程學院，南京 210009

雙酚A類似物；內分泌干擾物；雌激素效應；內分泌干擾效應；定量結構-活性關系

Received 31 March 2016 accepted 5 May 2016

雙酚A(bisphenol A, BPA)是全球生產量最大的化學原料之一，被廣泛應用于日用品的生產[1]。然而，大量的研究表明，BPA具有較強的內分泌干擾等毒性效應[2-3]，可干擾激素受體[4-5]，進而影響胚胎/胎兒發育[6]，損害雄性生殖健康[7]等。考慮到BPA的潛在危害效應，部分國家已經開始限制使用含BPA的塑料制品，尤其是對嬰幼兒等敏感人群使用的制品，許多國家/地區直接采取了禁用政策[8-9]。為了應對管制措施，一些雙酚A類似物，如雙酚AF(bisphenol AF, BPAF)、雙酚S(bisphenol S, BPS)等被開發用來替代含BPA的工業化學品的生產。由于其與雙酚A結構相似，雙酚A類似物是否也具有潛在內分泌干擾效應就受到越來越多的關注。

研究表明，雙酚A類似物不僅具有急、慢性毒性效應[10]，也具有內分泌干擾效應[11-13]，且部分雙酚A類似物的雌激素活性和抗雄性激素活性比BPA更強。例如，Matsushima等[14]的體外細胞實驗結果顯示，BPAF與雌激素受體α(ERα)的親和力相當于BPA的20倍，與雌激素受體β(ERβ)的親和力相當于BPA的48倍。Kitamura等[11]也通過熒光素酶報告基因研究了19種BPA類似物在MCF-7細胞中的雌激素活性，發現這類物質在10-8～10-4mol·L-1就能刺激MCF-7細胞的增殖；此外，還發現除四溴雙酚A和四氯雙酚A外，其他幾種物質均具有較強的抗雄性激素活性。Ji等[15]的研究表明，當把雄性和雌性斑馬魚暴露于0.5 μg·L-1～ 50 μg·L-1BPS的水樣中21 d后，雌性和雄性斑馬魚血漿中的雌二醇含量顯著升高，而雄性斑馬魚血漿中的睪酮含量則顯著降低。因此，雙酚A類似物的潛在內分干擾效應也應得到足夠重視。目前，盡管已有研究表明部分雙酚A類似物具有雌激素活性和抗雄性激素活性等內分泌干擾效應，但是主要是體外試驗結果為主，考慮到雙酚A類似物種類較多，難以對雙酚A類似物一一開展整體動物試驗研究。為解決該困境，一方面需要新的測試策略，以快速填補缺失的數據，另一方面需要通過設定優先級，篩查高優先級的雙酚A類似物，減少進行整體動物測試的化學品數量。

計算毒理學技術在設定優先級及填補數據缺失方面具有一定的優勢。在雙酚A類似物的內分泌干擾效應研究方面，針對雙酚A類似物分子結構與其雌激素干擾效應關系的研究還較少[16-19]。現有的研究主要涉及作用機制研究和構建雙酚A類似物的雌激素干擾效應預測模型。Delfosse等[17]分析了BPA及其2種類似物與雌激素受體的作用模式。Zhuang等[18]采用分子動力學模擬方法，研究了BPA及其7種類似物與雌激素受體alpha(ERα)的作用機制，發現靜電作用對雙酚A類似物的雌激素效應具有重要作用。Coleman等[16]測試了24種雙酚A類似物的雌激素效應，并采用比較分子力場法(CoMFA)、比較分子相似性指數(CoMSIA)方法和分子全息QSAR方法(HQSAR)構建了預測模型。該研究未基于構建的模型預測其他雙酚A類似物的雌激素效應數據。Ng等[19]基于18個化合物，采用分子對接方法，計算了小分子與雌激素受體的結合能力(得分)，建立了化合物雌激素效應數據與打分之間的關系，并用構建的關系填補了其他27種化合物(其中雙酚A類似物17種)的雌激素活性數據。本研究將采用定量結構-活性關系(QSAR)方法，基于24個有雌激素效應的數據，構建雙酚A類似物雌激素效應的預測模型，并填補其他雙酚A類似物缺失的雌激素效應數據。

1 材料與方法 (Materials and methods)

1.1 數據集

數據集I包含24個雙酚A類似物的雌激素效應值。實驗數據來源于Coleman等[16]的測定結果。實驗數據是采用人乳腺癌細胞(MCF-7細胞)增殖實驗獲取。實驗測定的指標是半數效應濃度(EC50)。此外，本文還收集了另外22個無試驗測試數據的雙酚A類似物(數據集II)，將用本研究構建的模型預測這部分化合物的雌激素效應值。所有化合物信息及其效應值列于表1。

QSAR建模的指標是取對數的相對效應勢(logRP)，其定義為：

(1)

其中，EC50,E2和EC50,BPA分別是17β-雌二醇(E2)和雙酚A類似物的半數效應濃度。

1.2 分子結構優化

首先采用ChemBioOffice 2010軟件生成初始的化合物分子結構，進而根據上述分子結構生成Gaussian輸入文件，然后運用Gaussian 09優化化合物分子結構[20]。具體而言，采用Gaussian 09程序包中的密度泛函理論(DFT)中的B3LYP方法，在6-31+g(d,p)基組水平上進行結構優化，獲得其穩定的分子構象。對優化好的幾何構型進行頻率分析，以確保體系無虛頻。在上述幾何優化和DFT的計算中都使用自洽場(SCRF)和積分方程形式極化連續介質模型(IEFPCM)來考慮溶劑(水)效應[21-22]。

1.3 分子描述符計算

其次，基于Gaussian 09優化的結構，使用Dragon 6.0軟件計算描述符[23]。Dragon 6.0包括29大類4 885種分子描述符，是一個全面且高效的分子描述符計算軟件。基于如下標準，對計算的4 885種描述符進行了初步篩選：去掉常數和近似常數的描述符，去掉至少有一個缺失值的描述符。最終，描述符集包含914種描述符。

表1 雙酚A及其類似物信息、雌激素效應數據

注：*驗證集化合物。

Note: * Compounds selected as the external validation set.

1.4 QSAR模型構建與表征

(2)

此外，還計算了均方根誤差(RMSE)和模型的平均絕對誤差(MAE)，以評估預測方法的可靠性。其中，MAE計算式如下[26]：

(3)

其中，logRP (實驗值)和logRP(預測值)分別是實驗測定和預測的logRP值，n是化合物數。

采用歐幾里德距離圖和基于杠桿值(leverage)的Williams圖定義模型的應用域。歐幾里德距離圖使用AMBIT Discovery v0.04 (http://ambit.sourceforge.net/download_ ambitdiscovery.html)軟件計算。

Williams圖是由標準殘差(δ*)和杠桿值(以hi表示，i代表不同的化合物)定義的一種QSAR模型應用域[24]。其中，δ*使用下式計算：

(4)

訓練集化合物的H值計算式為：

H = X(XTX)-1XT

(5)

其中，X為n × m的矩陣，n為訓練集化合物的個數，m為模型中預測變量的個數，X矩陣表征了模型中訓練集化合物的描述符空間。

訓練集和驗證集中化合物的hi值可由式6計算得到：

hi= xiT(XTX)-1xi

(6)

其中，xi為訓練集和驗證集中第i個化合物預測變量的行向量。

警戒值(h*)定義式為：

h*= 3(m + 1)/n

(7)

其中，m為預測變量的個數。

由于使用SPSS逐步回歸時，可能得到一系列模型，本研究基于上述的模型表征參數、方差膨脹因子(VIF)值及描述符個數要求等指標，選取最優模型。

2 結果與討論(Results and discussion)

2.1 最優模型及其表征結果

最優模型為：

logRP = -10.7-2.53Mor17u+64.1Vs,max

(8)

2.2 應用域表征

在OECD關于QSAR構建與驗證導則中[24]，明確指出要對所建模型進行應用域表征。模型應用域表征結果如圖1和圖2所示。模型定義的歐幾里德距離小于0.748。訓練集和驗證集中化合物都在模型應用域中，說明模型訓練集具有較好的代表性。此外，在Williams圖中(圖2)，若化合物的標準殘差δ*落在(-3.0, +3.0)以外時，認為該點是離群點。如圖可見，所有化合物的標準殘差δ*均落于(-3.0, +3.0)以內，且所有化合物的杠桿值均小于警戒值h*，說明模型沒有離群點。因此，該模型能夠用于預測應用域內其他雙酚A類似物的抗雄激素效應值。

圖1 基于歐幾里德距離方法表征的模型應用域Fig. 1 Characterization of application domain for developed model based on the Euclidean distance

圖2 模型Williams圖Fig. 2 Williams plot of developed model

圖3 靜電勢分布圖Fig. 3 Electrostatic potentials on molecular surface

圖4 數據集II中化合物的適用性評估Fig. 4 Applicability evaluation of the compounds in the data set II

2.3 描述符意義與機理解釋

所建模型包含2個描述符，即Mor17u和Vs,max，Mor17u是3D-MoRSE類描述符，Vs,max是量化描述符。其中，Mor17u是未加權的3D-MoRSE描述符；3D-MoRSE(3D molecule representation of structures based on electron diffraction)描述符是基于電子衍射表征的分子三維結構描述符[29]，表征了官能團的影響。與雙酚A相比，分子結構中引入吸電子取代基團會使Mor17u值增加，而引入供電子取代基團會使Mor17u值減小。例如，當雙酚A羥基鄰位被溴或氯原子取代時，四溴雙酚A和四氯雙酚A的Mor17u值從雙酚A的-1.072，增加至-0.559和-0.319；而當雙酚A羥基鄰位被1或2個甲基取代時，2,2-雙(4-羥基-3,5-二甲基苯基)丙烷和4,4'-(1-甲基亞乙基)雙(2-甲基苯酚)的Mor17u值從雙酚A的-1.072，減小至-1.273和-1.479。方程中Mor17u為負系數(r = -0.80, P< 0.0001)，說明化合物logRP值隨Mor17u值的減小而增加。因此，當雙酚A羥基鄰位引入吸電子取代基團會減弱分子的雌激素效應，而引入供電子取代基團則會增強分子的雌激素效應。

圖5 基于模型預測的數據集II中化合物logRP值Fig. 5 Predicted logRP values of compounds in data set II based on developed model

Vs,max是分子表面上的最正靜電勢。圖3顯示了分子靜電勢分布，對于含酚羥基的化合物(如圖3中雙酚A)，分子中靜電勢最正的區域是酚羥基對應的部分。前人的研究表明，與雌激素受體形成氫鍵相互作用、靜電作用是小分子表現雌激素效應的關鍵[18,21-22]。對雙酚A類似物而言，酚羥基是其與雌激素受體形成氫鍵相互作用、靜電作用的核心基團[18]。因此，該描述符表征了化合物分子與雌激素受體之間的氫鍵相互作用、靜電作用。從圖4可以看出，酚羥基鄰位取代基也會影響分子表面正靜電勢的分布。與雙酚A相比，分子結構中引入吸電子取代基團會使分子表面正靜電勢減小。在方程中Vs,max為正系數(r = 0.75, P< 0.0001)，說明化合物logRP值隨Vs,max值的增加而增加。因此，當雙酚A酚羥基鄰位引入吸電子取代基團會減弱分子的雌激素效應。綜上，模型選出的2個變量都解釋了取代基對雙酚A類似物雌激素效應的影響。

2.4 預測其他雙酚A類似物的雌激素效應值

本文采用構建的模型預測其他22個雙酚A類似物的雌激素效應數據。首先計算了模型對這22個化合物的適用性。計算結果如圖4所示，有3個化合物歐幾里德距離大于0.748，處于模型應用域外，3個化合物分別是2,2-二苯基丙烷、二苯基甲烷和雙酚A雙(二苯基磷酸酯)。由于模型訓練集化合物分子結構中都含有酚羥基，因而2,2-二苯基丙烷和二苯基甲烷分子結構中不含酚羥基基團是導致其處于模型應用域外的原因。對雙酚A雙(二苯基磷酸酯)而言，其分子結構中含2個二苯基磷酸酯基，與訓練集化合物分子結構差異較大。綜上，模型對其余19個化合物適用。

模型對數據集II中化合物的預測值列于表1。此外，前人的研究結果表明，在該實驗體系下，17β-雌二醇(E2)的半數效應濃度(EC50)約為10-11mol·L-1。根據E2的EC50值和式(1)，計算了數據集I中化合物實測的EC50值和數據集II中化合物預測的EC50值(表1)。圖5顯示了預測的效應值分布范圍及其對應的化合物數量。結果表明，數據集II中化合物logRP值介于-8 ～ -2之間，而數據集I(建模數據集)其logRP值介于-9 ～ -3之間。數據集II中僅有3個雙酚A類似物的預測logRP值介于-3 ～ -2，超出了建模數據集化合物的logRP值范圍。但是由于這3個化合物在模型應用域內，說明預測值可信。

從logRP值的定義可以看出，其值越大，表示對應化合物的雌激素效應越強。從表1可知，有14個雙酚A類似物的logRP值大于雙酚A。根據EC50值可以看出，數據集II中有7個化合物的EC50值僅是雙酚A類似物EC50值的二十分之一。在將來的研究中，需進一步試驗確認其雌激素效應。

對于不在模型應用域內的3個化合物，模型預測表明，2,2-二苯基丙烷、二苯基甲烷和雙酚A雙(二苯基磷酸酯)的logRP分別為-7.02、-7.56和-3.03，即表明預測的2,2-二苯基丙烷和二苯基甲烷雌激素效應較弱，而雙酚A雙(二苯基磷酸酯)雌激素效應較強。根據上述分析，分子結構中不含酚羥基基團是導致2,2-二苯基丙烷和二苯基甲烷處于模型應用域外的原因。由于分子結構中酚羥基基團的缺失，導致2,2-二苯基丙烷和二苯基甲烷不能與雌激素受體形成氫鍵相互作用。因此，2,2-二苯基丙烷和二苯基甲烷雌激素效應較弱的預測結果可信度較高。對于雙酚A雙(二苯基磷酸酯)而言，由于其分子結構與訓練集化合物差異最大，且分子中不含酚羥基等關鍵基團。在模型訓練集中，同為酯類的雙酚A雙(二甲基酯)的雌激素效應較弱(logRP = -7.36)，因此，雙酚A雙(二苯基磷酸酯)的雌激素效應預測數據不確定性較大。

2.5 模型比較

現有雙酚A類似物相關的預測模型情況列于表2。與表2中的模型相比，僅本研究所建預測模型嚴格按照OECD提出的QSAR模型構建與驗證準則進行了模型擬合優度、穩健性、預測能力及應用域的表征，且還進行了機理解釋。因此，本研究所建預測模型能夠較好地預測應用域內其他雙酚A類似物的雌激素效應數據。

表2 不同模型之間的比較

注：alogRP是取對數的相對效應勢，其定義為式(1)；logRBA是取對數的相對親和勢，其定義為：

其中，IC50,E2和IC50,BPAs分別是17β-雌二醇(E2)和雙酚A類似物的半數抑制濃度。

c原文未計算驗證集表征參數，表中數據是根據原文驗證集的實驗效應數據和預測效應數據計算(原文表5)。

Note:alogRP is the logarithm of the relative competing potency and is defined by equation 1; logRBA is logarithm of relative binding affinity, which is defined as:

where IC50,E2and IC50,BPAsare the concentrations of 17β-estradiol and bisphenols at 50% inhibition of 17β-estradiol binding from estrogen receptor, respectively.

cThe statisticsis calculated by employing the observed and predicted data of the compounds in validation set (Table 5 in the original paper) becuase of the statistics for validation set is not calculated in the original paper.

[1] Eladak S, Grisin T, Moison D, et al. A new chapter in the bisphenol A story: Bisphenol S and bisphenol F are not safe alternatives to this compound [J]. Fertility and Sterility, 2015, 103(1): 11-21

[2] Rochester J R. Bisphenol A and human health: A review of the literature [J]. Reproductive Toxicology, 2013, 42: 132-155

[3] Glausiusz J. Toxicology: The plastics puzzle [J]. Nature, 2014, 508(7496): 306-308

[4] Sui Y, Ai N, Park S H, et al. Bisphenol A and its analogues activate human pregnane X receptor [J]. Environmental Health Perspectives, 2012, 120(3): 399-405

[5] Delfosse V, Grimaldi M, le Maire A, et al. Nuclear receptor profiling of bisphenol-A and its halogenated analogues [J]. Vitamins & Hormones, 2014, 94: 229-251

[6] Virginie R F, René H, Gabriel L. Effects of endocrine disruptors on the human fetal testis [J]. Annales Dendocrinologie, 2014, 75(2): 54-57

[7] Li D K, Zhou Z J, Miao M, et al. Urine bisphenol-A (BPA) level in relation to semen quality [J]. Fertility and Sterility, 2011, 95(2): 625-630

[8] Commission Directive 2011/8/EU of 28 January 2011 amending Directive 2002/72/EC as regards the restriction of use of Bisphenol A in plastic infant feeding bottles [EB/OL]. (2011-01-28) [2015-12-17]. http://faolex.fao.org/docs/pdf/eur100741.pdf

[9] Health Canada. Archived — Order Amending Schedule I to the Hazardous Products Act (bisphenol A) [EB/OL]. (2010-03-31) [2015-12-17]. http://www.gazette.gc.ca/rp-pr/p2/2010/2010-03-31/html/sor-dors53-eng.html

[10] Chen M Y, Ike M, Fujita M. Acute toxicity, mutagenicity, and estrogenicity of bisphenol-A and other bisphenols [J]. Environmental Toxicology, 2002, 17(1): 80-86

[11] Kitamura S, Suzuki T, Sanoh S, et al. Comparative study of the endocrine-disrupting activity of bisphenol A and 19 related compounds [J]. Toxicological Sciences, 2005, 84(2): 249-259

[12] Riu A, Grimaldi M, le Maire A, et al. Peroxisome proliferator-activated receptor γ is a target for halogenated analogs of bisphenol A [J]. Environmental Health Perspectives, 2011, 119(9): 1227-1232

[13] Molina-Molina J M, Amaya E, Grimaldi M, et al. In vitro study on the agonistic and antagonistic activities of bisphenol-S and other bisphenol-A congeners and derivatives via nuclear receptors [J]. Toxicology and Applied Pharmacology, 2013, 272(1): 127-136

[14] Matsushima A, Liu X, Okada H, et al. Bisphenol AF is a full agonist for the estrogen receptor ERα but a highly specific antagonist for ERβ [J]. Environmental Health Perspectives, 2010, 118: 1267

[15] Ji K, Hong S, Kho Y, et al. Effects of bisphenol S exposure on endocrine functions and reproduction of zebrafish [J]. Environmental Science & Technology, 2013, 47(15): 8793-8800

[16] Coleman K P, Toscano Jr W A, Wiese T E. QSAR models of the in vitro estrogen activity of bisphenol A analogs [J]. QSAR & Combinatorial Science, 2003, 22: 78-88

[17] Delfosse V, Grimaldi M, Pons J L, et al. Structural and mechanistic insights into bisphenols action provide guidelines for risk assessment and discovery of bisphenol A substitutes [J]. Proceedings of the National Academy of Sciences of the United States of America, 2012, 109(37): 14930-14935

[18] Zhuang S L, Zhang C L, Liu W P. Atomic insights into distinct hormonal activities of bisphenol A analogues toward PPARγ and ERα receptors [J]. Chemical Research in Toxicology, 2014, 7(10): 1769-1779

[19] Ng H W, Shu M, Luo H, et al. Estrogenic activity data extraction and in silico prediction show the endocrine disruption potential of bisphenol A replacement compounds [J]. Chemical Research in Toxicology, 2015, 28(9): 1784-1795

[20] Frisch M J, Trucks G W, Schlegel H B, et al. Gaussian Inc. Wallingford CT [CP]. 2009

[21] Li F, Li X H, Shao J P, et al. Estrogenic activity of anthraquinone derivatives: In vitro and in silico studies [J]. Chemical Research in Toxicology, 2010, 23(8): 1349-1355

[22] 李斐. 部分有機污染物雌激素效應和甲狀腺激素效應的計算模擬與驗證[D]. 大連: 大連理工大學, 2010: 56-59

Li F. Computational simulation and validation on the estrogenic and thyroid hormone activities of some organic chemicals [D]. Dalian: Dalian University of Technology, 2010: 56-59 (in Chinese)

[23] Talete srl, DRAGON (Software for Molecular Descriptor Calculation) Version 6.0-2012 [CP]. http://www.talete.mi.it/

[24] Organisation for Economic Co-operation and Development (OECD). Guidance document on the validation of (quantitative) structure activity relationships ((Q)SAR) models, Technical Report for OECD Environment, Health and Safety Publications Series on Testing and Assessment No. 69 [R]. Paris: OECD, 2007

[25] Schüürmann G, Ebert R U, Chen J W, et al. External validation and prediction employing the predictive squared correlation coefficient — Test set activity mean vs training set activity mean [J]. Journal of Chemical Information and Modeling, 2008, 48(11): 2140-2145

[26] Hyndman R J, Koehler A B. Another look at measures of forecast accuracy [J]. International Journal of Forecasting, 2006, 22(4): 679-688

[27] Golbraikh A, Shen M, Xiao Z, et al. Rational selection of training and test sets for the development of validated QSAR models [J]. Journal of Computer-Aided Molecular Design, 2003, 17(2-4): 241-253

[28] Kaneko H, Funatsu K. Criterion for evaluating the predictive ability of nonlinear regression models without cross-validation [J]. Journal of Chemical Information and Modeling, 2013, 53(9): 2341-2348

[29] Schuur J, Gasteiger J. Infrared spectra simulation of substituted benzene derivatives on the basis of a 3D structure representation [J]. Analytical Chemistry, 1997, 69(13): 2398-2405

[30] Wiese T E, Kral L G, Dennis K E, et al. Optimization of estrogen growth response in MCF-7 cells [J]. In Vitro Cellular & Developmental Biology - Animal, 1992, 28A(9-10): 595-602

Predicting Estrogenic Activity of Bisphenols Using Quantitative Structure-Activity Relationships

Yang Xianhai1, Liu Huihui2, Yang Qian1,3, Liu Jining1,*

1. Nanjing Institute of Environmental Science, Ministry of Environmental Protection, Nanjing 210042, China 2. Jiangsu Key Laboratory of Chemical Pollution Control and Resources Reuse, School of Environmental and Biological Engineering, Nanjing University of Science and Technology, Nanjing 210094, China 3. The College of Biotechnology and Pharmaceutical Engineering, Nanjing Tech University, Nanjing 210009, China

bisphenols; endocrine disrupting chemicals; estrogenic activity; endocrine disrupting effects; quantitative structure-activity relationships

國家自然科學基金(No.201507038)；江蘇省自然科學基金(No.BK20151100)；中央級公益性科研院所基本科研業務專項(2015)

楊先海(1985-)，男，工學博士,助理研究員，研究方向為環境內分泌干擾效應的計算模擬研究，E-mail: xhyang@nies.org；

*通訊作者(Corresponding author)， E-mail: ljn@nies.org

10.7524/AJE.1673-5897.20160331015

2016-03-31 錄用日期：2016-05-05

1673-5897(2016)4-069-10

X171.5

A

簡介：劉濟寧(1977-)，男，綠色化學工學博士，副研究員，主要研究方向為化學品環境毒理學研究，發表學術論文70余篇。

楊先海, 劉會會, 楊倩, 等. 雙酚A類似物雌激素干擾效應的定量結構-活性關系模型[J]. 生態毒理學報，2016, 11(4): 69-78

Yang X H, Liu H H, Yang Q, et al. Predicting estrogenic activity of bisphenols using quantitative structure-activity relationships [J]. Asian Journal of Ecotoxicology, 2016, 11(4): 69-78 (in Chinese)