Prediction and assessment of xenoestrogens mixture effects using the in vitro ERα-CALUX assay

Introduction: Many natural or synthetic compounds used in foods, dietary supplements, and food contact materials (FCMs) are suspected endocrine disruptors (EDs). Currently, scientific evidence to predict the impacts on biological systems of ED mixtures is lacking. In this study, three classes of substances were considered: i) phytoestrogens, ii) plant protection products (PPP) and iii) substances related to FCMs. Fourteen compounds were selected based on their potential endocrine activity and their presence in food and FCMs.

Methods: These compounds were evaluated using an in vitro gene expression assay, the ERα-CALUX, to characterize their responses on the estrogen receptor alpha. Cells were exposed to fixed ratio mixtures and non-equipotent mixtures of full and partial agonists. The concentration-response curves measured for the three classes of compounds were characterized by variable geometric parameters in terms of maximum response (efficacy), sensitivity (slope) and potency (median effective concentration EC50). To account for these variations, a generic response addition (GRA) model was derived from mass action kinetics.

Results: Although GRA does not allow us to clearly separate the concentration addition (CA) and independent action (IA) models, it was possible to determine in a statistically robust way whether the combined action of the chemicals in the mixture acted by interaction (synergy and antagonism) or by additive behavior. This distinction is crucial for assessing the risks associated with exposure to xenoestrogens. A benchmark dose approach was used to compare the response of phytoestrogen blends in the presence and absence of the hormone estradiol (E2). At the same time, 12 mixtures of 2–5 constituents including phytoestrogens, phthalates and PPPs in proportions close to those found in food products were tested. In 95% of cases, the response pattern observed showed a joint and independent effect of the chemicals on ER.

Discussion: Overall, these results validate a risk assessment approach based on an additive effects model modulated by intrinsic toxicity factors. Here, the CA and IA approaches cannot be distinguished solely based on the shape of the concentration response curves. However, the optimized GRA model is more robust than CA when the efficacy, potency, and sensitivity of individual chemical agonists show large variations.

1 Introduction

Endocrine disruptors (EDs) are substances that affect the endocrine system and have adverse effects on the health of an organism and/or its offspring. EDs affect the hormonal balance and impact human health by modifying functions such as growth, development, sexual and reproductive function, behavior, obesity, and blood circulation (Diamanti-Kandarakis et al., 2009; Guarnotta et al., 2022). The increasing incidence of endocrine diseases and associated disabilities is assumed related to the increasing exposure of the population to EDs (Ahn and Jeung, 2023). These chemicals disrupt the body at different levels, often in discrete ways and in doses as low as the order of physiological hormone concentrations (Marlatt et al., 2022). The impact on human health is even more worrying since the effect of these substances on fetal development may be irreversible. It was estimated that the annual cost of EDs for European society is of the order of hundreds of billions of euros per year (Trasande et al., 2016).

According to REACH (updated list in April 2022), 24 chemicals were identified and legally adopted as endocrine disruptors for humans. These include phthalates, and the bisphenol and phenol groups. However, many other compounds are considered potential EDs including antioxidants, pesticides, biocides, detergents, metals, plasticizers, and organometallic compounds (Mukherjee et al., 2021; Wiesinger et al., 2021; Imparato et al., 2022). Diet is an important exposure route for these known or potential EDs, which enter the food chain in different ways. They can occur naturally as phytoestrogens, deliberately added as additives, or can contaminate food through production, processing, transportation, and migration from food contact materials (FCM) (Ong et al., 2020; Wiesinger et al., 2021). One of the major challenges with current testing and regulation of EDs is that they are assessed on an individual basis, whereas in real life the population is exposed to multiple compounds (“mixtures”) simultaneously. EDs enter the human body as complex mixtures, usually through chronic low-dose exposure. This implies that the organism may be exposed to the same compound from several sources or that a consumer product may contain several known or potential EDs (Kumar et al., 2020; Darbre, 2022).

Current practice is to determine the toxicity of individual compounds, after which the mixture is considered “tolerable” when the toxicity or concentration of individual compounds is below experimentally defined toxicity threshold values (Heys et al., 2016). However, the fact that the toxicity of a single compound is known does not guarantee its behavior when combined with other substances. Compounds that share a common mechanism of action may be able to cooperate with each other and induce a measurable effect even if they are individually present below their No Observed Adverse Effect Level—NOAEL (Silva et al., 2002). These joint effects therefore complicate risk assessment, hence the need to extend current methods of toxicological assessment (Perkins et al., 2019).

The two most used models to predict joint effects in mixtures are concentration addition (CA) and independent action (IA), both having various names in the literature (Cedergreen et al., 2008; Medlock Kakaley et al., 2019; Hsieh et al., 2021). The theoretical principle of CA is that the chemicals in the mixture share the same target molecule (same mode of action) and dose-response curves differ only in potency (Gottschalk and Dunn, 2005). They can, therefore, be considered as dilutions of each other. In the IA model, the chemicals in the mixture are assumed to have different modes of action, thus affecting different potential biological target sites. Since they act independently, the overall mixture might show no response when the compounds are present below their NOAEL (Martin et al., 2009). A common assumption between CA and IA is that chemicals in a mixture do not interact physically, chemically, or biologically. However, interactions between different chemicals are possible and can have a substantial impact on the overall potency or toxicity of the mixture. Either the result of the mixture exceeds that predicted by the sum of the individual effects, which is called synergy, or the result of the mixture is less than that expected, which is called antagonism (Cedergreen et al., 2017; Martin et al., 2021). A thorough study of the literature shows that cocktail effects for EDs with interactions are only observed in 11% of the cases studied. Overall, the deviations were small, and very often additivity could be used as the default concept to predict mixture outcomes (Martin et al., 2021).

In this study, to better appreciate these possible cocktail effects, the estrogenic activity of mixtures consisting of several known and/or potential EDs related to food and FCM was studied using the bio- ERα-CALUX assay. The experimental results were compared to the activity predicted by models derived from CA and IA concepts, to test whether the mixtures followed the additivity principle or showed signs of interaction. Mixtures have also been prepared containing compounds below their activity threshold to determine if the overall mixture activity was significant, illustrating the so-called “something out of nothing” paradigm (Silva et al., 2002).

2 Material and methods

2.1 Prioritization of potential estrogenic disruptors

Potential endocrine active substances present in foods for study come from the categories: i) phytoestrogens, ii) biocides and phytosanitary products and iii) substances present in the form of additives or present in FCMs. For each of these categories, a list based on the presence on the Belgium food market was established (see Bel et al., 2016), and a selection based on the potential for endocrine disruption was made after literature review (ex. Brouwers et al., 2009; EFSA, 2019). The selection of these compounds of interest guided the selection of foods and FCMs, which were analyzed by LC/GC-MS/MS as part of the FPS Public Health project ENDFOODTOX (RF 18/6326).

2.2 Analytical standards

The following analytical standards were purchased: i) from Sigma Aldrich: butylated hydroxytoluene (BHT), Benzophenone, Bis (2-ethylhexyl) phthalate (DEHP), Dibutyl phthalate (DBP), Diisobutyl phthalate (DIBP), Kaempferol, and Equol, ii) from Phytilab: Dadzein, Genistein, 6-Prenylnaringenin (6-PN), 8-Prenylnaringenin (8-PN) and Enterolactone and iii) from LGC: Benzyl butyl phthalate (BBP) and Triadimenol.

Each reference material was solubilized in appropriate solvents to make stock solutions at a concentration of 1 mg/mL. Most compounds were solubilized in MeOH, except daidzin and genistin, which required DMSO. BHT is soluble in MeOH, but anhydrous ethanol or hexane is required to prevent its degradation. Working standards at 20 μg/mL and 10 μg/mL were prepared in MeOH and ACN. Stock and working standards were stored at −20°C and retained for a maximum of 6 months.

2.3 Assessment of endocrine activities using the ERα-CALUX assay

The recombinant human breast cancer cell line VM7Luc4E2 (variant MCF7, formerly known as BG1Luc4E2) was used to determine estrogenic activities. These cells express ERα endogenously but lack any functional ERβ (Rogers and Denison, 2000; Brennan et al., 2016). The bioanalytical procedure and data analysis are described in (Elskens et al., 2023).

2.4 Interpretation of concentration-response curves for mixtures

A growing body of experimental evidence indicates that the in vitro activities of chemical mixtures can be predicted from the overall response of their individual components using the concept of the concentration addition model-CA (Evans et al., 2012). Mathematically, CA models assume that all dose-response curves are characterized by similar threshold values and similar efficacy and that their slope-values are not significantly different from 1.

${\mathit{E}\mathit{C}}_{{\mathit{X}}_{\mathit{M}\mathit{I}\mathit{X}}}={\left(\sum _{\mathbf{1}=\mathbf{1}}^{\mathit{n}}\frac{{\mathit{f}}_{\mathit{i}}}{{\mathit{E}\mathit{C}}_{{\mathit{X}}_{\mathit{i}}}}\right)}^{-\mathbf{1}}(1)$

where ${\mathit{E}\mathit{C}}_{{\mathit{X}}_{\mathit{i}}}$ is the effective concentration of compound i alone causing an effect X, f_i is the fraction of i present in the mixture and ${\mathit{E}\mathit{C}}_{{\mathit{X}}_{\mathit{M}\mathit{I}\mathit{X}}}$ is the effective concentration of the mixture causing the same effect X. To apply Eq. 1, the individual potencies of the compounds present in the mixture are required (Supplementary Table S1).

Assumptions that all concentration-response curves are characterized by similar threshold, efficacy and slope-values are not valid for phytoestrogen response on ER (Supplementary Table S1). The Hill slope are frequently different from 1, and the maximum response fluctuate depending on the presence of full-, partial-, or supra-agonists. To correct for these inconsistencies, a logistic function was derived from the mass action kinetics using chemical equilibrium reactions:

$\mathit{R}+{\sum }_{\mathit{i}=\mathbf{1}}^{\mathit{n}}{\mathit{n}}_{\mathit{i}}\bullet {\mathit{L}}_{\mathit{i}}\iff {\sum }_{\mathit{i}=\mathbf{1}}^{\mathit{n}}{\mathit{R}\mathit{L}}_{{\mathit{n}}_{\mathit{i}}}(2)$

where R is the receptor protein concentration, L_i the free unbound ligand concentration with n_i binding sites. From the law of mass action, the apparent dissociation constant K_i for ligand L_i is given by:

${\mathit{K}}_{\mathit{i}}=\frac{\left[\mathit{R}\right]\bullet {\left[{\mathit{L}}_{\mathit{i}}\right]}^{{\mathit{n}}_{\mathit{i}}}}{\left[\mathit{R}{\mathit{L}}_{\mathit{i}{\mathit{n}}_{\mathit{i}}}\right]}(3)$

The ratio of occupied receptor to total receptor is then:

$\mathit{\Phi }=\frac{\sum _{\mathit{i}=\mathbf{1}}^{\mathit{n}}{\mathit{R}\mathit{L}}_{\mathit{i}{\mathit{n}}_{\mathit{i}}}}{\mathit{R}+\sum _{\mathit{i}=\mathbf{1}}^{\mathit{n}}{\mathit{R}\mathit{L}}_{\mathit{i}{\mathit{n}}_{\mathit{i}}}}(4)$

Solving for the receptor activity, one finally obtains:

${\mathit{y}}_{\mathit{i}}={\mathit{y}}_{\mathbf{0}}+\left(\mathit{m}-{\mathit{y}}_{\mathbf{0}}\right)\bullet \frac{{\sum }_{\mathit{i}=\mathbf{1}}^{\mathit{n}}{\mathit{a}}_{\mathit{i}}\bullet {\left(\frac{{\mathit{x}}_{\mathit{i}}}{{\mathit{c}}_{\mathit{i}}}\right)}^{{\mathit{n}}_{\mathit{i}}}}{\mathbf{1}+{\sum }_{\mathit{i}=\mathbf{1}}^{\mathit{n}}{\left(\frac{{\mathit{x}}_{\mathit{i}}}{{\mathit{c}}_{\mathit{i}}}\right)}^{{\mathit{n}}_{\mathit{i}}}}(5)$

which is a generic formulation of the Hill equation for several binding ligands. In Eq. 5 y₀ and m are the lower and upper asymptote for E2, a_i = efficacy, c_i = EC₅₀, n_i = Hill slope of ligand L_i.

The parameter n_i represents a “sensitivity coefficient” describing the slope of the tangent to the dose response curve at its inflection point for each ligand. Eq. 5, referred to as GRA, is hereafter considered equivalent to a mixed CA/IA approach, because it makes possible to consider the variability of the parameters defining the potency, sensitivity, and amplitude of the concentration-response curves for each compound tested in a mixture. Values for a_i, c_i, and n_i are summarized in (Supplementary Table S1). CA and IA models assumes the additivity principle between tested compounds. It is the default assumption for predicting and evaluating the mixing effect. The accuracy of these models was assessed by comparing the predicted curves against the observed values with the root mean square error (RMSE), which measures the average difference between a model’s predicted values and the actual values. In addition, the measured concentration-response curves of the mixture are compared with those obtained by modelling, and the equivalence hypothesis rejected if the measured curve is outside the 95% confidence interval of the theoretical models obtained by means of Monte-Carlo simulations carried out using XLSTAT-Sim (XLSTAT version 2023.1.5) according to normal and/or uniform distribution laws.

2.5 Quality control and quality assurance (QC/QA)

The best-fit parameters for the concentration response-curves were obtained by minimizing the least squares residuals using a Levenberg-Marquardt algorithm in R version 4.1.1-2021-08-10 (Elskens et al., 2011). For the standard E2 curves, the parameter-values (95% CI) under repeatability conditions and the goodness of fit criteria (min—max) are given in Supplementary Table S1. Sample activity below the threshold values of 20% was designated as < LoQ. Sample activity below the threshold values of 10% was designated as not detectable (<LoD). Quality control samples were systematically performed in triplicate on the 96-well plates. They consisted of standard at the half maximal effective concentration (EC₅₀). The recovery rates were between 89% and 120% for E2, results which can be considered as satisfactory.

3 Results and discussion

3.1 Rationale and assumption behind the derivation of the Ca and GRA models

Concentration-response curves were used to study the in vitro effects of mixtures, assuming that the bioactivity of these mixtures followed by default the principle of additivity according to Eq. 1. This model assumes that all the ligands act in the same way as dilutions of each other. While the binding affinity on the receptor may differ between compounds (i.e., relative potency), they all exhibit the same slope (i.e., Hill coefficient) and maximum response (i.e., efficacy), resulting in parallelism between sets of concentration response curves (Gottschalk and Dunn, 2005). However, for the chemicals investigated in this study, the parallelism hypothesis is not valid; the efficacy (a) varying between 0.5 and 2 and Hill coefficient (n) between 0.5 and 5 (Supplementary Table S1). Consequently, a generalisation of the Hill-Langmuir function (Eq. 5) was derived from the mass action kinetics (Eq. 2) assuming that each ligand is characterised by an intrinsic potency (c_i), sensitivity (n_i) and maximum response (a_i), providing a generic response addition model (GRA).

Binding curves showing the characteristically sigmoid generated by CA and GRA for an equipotent binary mixture of fictitious agonists are illustrated in Figure 1. The effect of random noise (7%) on the shape of the response curve and the empirical adjustment of parameters was investigated (Supplementary Table S2). Both models give non-statistically different values for the fitted parameters when compared to actual values at the 5% level of significance. CA is a special case of GRA for (a_i = 1, n_i = 1) but as n_i or a_i increases, the saturation curves become steeper or higher for GRA (Figures 1A, B). The question arises as to the physiological significance of these parameters in the CALUX assay. The Hill coefficient plays a leading role in the study of ligand-receptor interactions, measuring the number of binding sites in cooperating systems (Abeliovich, 2005). For n_i >1, when a ligand is bound to the receptor, its affinity for other ligands increases. For n_i <1, once a ligand is bound to the receptor, its affinity for other ligands decreases. Finally, for n_i = 1, the affinity of the receptor for a ligand does not depend on whether other ligands are already bound to the receptor. In practice, however, n_i in the CALUX assay is an ‘apparent coefficient’ reflecting a series of processes. After ligand binding to the cytosolic receptor, translocation of the newly formed dimeric complex in the nucleus and binding to DNA, the CALUX cells induce different protein activities via transcription, translation, and enzymatic expression. While we cannot rule out positive or negative cooperative binding for n_i ≠1, neither can we rule out indirect or cumulative effects on either of the reaction sequences. The same reasoning applies to a_i. The concept of “super agonism” was described in the early 1980 s, in relation to peptide hormone analogues that produced greater functional responses than endogenous agonists. Whether these compounds can actually be more effective than endogenous agonists have long been debated, but relatively recent pharmacological evidence has indicated that super agonists may be more than an artefact (Langmead and Christopoulos, 2013). Consequently, it is risky to assume that the functional affinity operating at the level of the whole cell reflects the molecular ligand-receptor binding affinity; care must be taken not to over-interpret the significance or reliability of the fitted parameters, even if it can be shown empirically that a_i and n_i are statistically different from 1. Moreover, the GRA model allows the effect of antagonism to be considered. Whether a_i of a given compound tends towards 0, it can bind to the receptor but does not activate it (Figure 1C). This results in a decrease in agonist activity in the presence of increasing concentrations of antagonist, the inhibitory effect of which being estimated by the IC₅₀ (concentration of antagonist inhibiting 50% of the agonist response). Supplementary Table S2 shows how the parameters fitted with a background of 7% compare with the actual model outputs. As previously, the shape of the concentration-response curves is modulated by n_i. The inhibition reflects a reversible competitive antagonism. The latter can be overcome. At high concentrations, the effect of the antagonist is no longer observed but the sigmoid agonist curve shifts to the right and its “apparent” potency is, therefore, reduced (EC₅₀ increases). The higher the antagonist concentration, the greater the displacement (Figure 1D).

FIGURE 1

FIGURE 1. Concentration response-curves generated by CA and GRA for equipotent binary mixtures of fictitious agonists/antagonists (see Supplementary Table S2).

GRA, therefore, appears to be a generalisation of CA and, like it, assumes an additive pattern. The concentration-response curves being modulated by the intrinsic a_i and n_i values, allow us to consider both mixtures of competing agonists and agonists/antagonists. However, this approach does not allow us to assess the synergy between endogenous and exogenous agonists, i.e., the sigmoid curve of the mixture shifts to the left and its “apparent potency” increases (EC₅₀ decreases). Such synergies have been less frequently reported and might occur at lower concentration ratios (Cox et al., 2021). Interaction is, therefore, expressed as changes in amplitude (increase/decrease in efficacy), potency (increase/decrease in binding affinity to the receptor), sensitivity (increase/decrease in slope response) or all three together (Indurthi and Auerbach, 2021). It is essential to assess these interactions, for example, when antagonists or toxic compounds are present in a sample, the CALUX response may be partially masked or amplified. The aim is to determine i) which CA or GRA approaches give the best description in mixture of real life food samples, and ii) whether the investigated chemicals follow an additive pattern or show signs of synergistic and/or antagonistic activity. A random noise of 7% on the CALUX experimental data is close to reality for triplicate measurements. This results in RMSE values between 2.7 and 4.9 and RSQ values greater than 0.9903 (Supplementary Table S2). These goodness-of-fit parameters summarise the difference between the observed and predicted values that can be expected if the CA or GRA models describe the set of observations satisfactorily.

3.2 Mixture effects of phytoestrogens at EC₁₀ and EC₅₀ in the presence and absence of E2

In these experiments, a mixture of 4 phytoestrogens (8-PN, daidzein, kaempferol and equol) was studied using the benchmark dose (BMD) approach. The concentration-response curves measured at low dose (EC₁₀) and at their half-maximal responses (EC₅₀) were compared with those predicted by the CA and GRA models (Eqs 1, 5), respectively). All measured activities were expressed as % relative to the maximum response of the reference ligand E2; the negative control is represented by the DMSO blank.

Figure 2A depicts that joint effects were detected at EC₁₀ and there were no statistical differences between the GRA-predicted and observed curves. The two curves overlap, and the observations are within the interval defined by the two models. RSQ and RMSE values are reported in Supplementary Table S3. At EC₅₀, the observed curve is slightly shifted to the left and the fit is better with the CA model (Figure 2B). However, the observations are most always included in the interval defined by GRA and CA.

FIGURE 2

FIGURE 2. Luciferase activity in VM7Luc4E cells for phytoestrogens. #a: 8-PN, daidzein, kaempferol and equol at EC₁₀; #b: 8-PN, daidzein, kaempferol and equol at EC₅₀; #c: E2 at EC₅₀ + Mix daidzein, kaempferol and equol (see Supplementary Table S3).

We then checked whether the presence of phytoestrogens could modify the response of E2 tested at its EC₅₀ (Figure 2C). At high dilutions of the phytoestrogen mixture (≤0.1 μM), the E2 response was not significantly modified. At higher concentrations (≥1 µM), the presence of phytoestrogens strengthens the E2 agonist response. Again, the observations are included in the interval defined by GRA and CA, and the RSQ & RMSE-values of both models are comparable (Supplementary Table S3). The mixture does not strongly modify the activity of E2. At high dilutions, there is a slightly higher overall activity of the mixture compared to the predictions, which decreases approaching saturation (Figure 2C). The fact that concentration-response curves can be steeper at lower exposure levels has been described for a variety of environmental toxicants, including lead, arsenic, and benzene (Hornung and Lanphear, 2014). Overall, current research on synergistic mixing effects is limited, and when discussed in the literature, their occurrence is low. A study on pesticide mixtures, for example, which examined the effects of mixing in realistic low-dose mixtures, reported that only 5% of the mixture combinations tested exhibited synergistic interactions (Cedergreen, 2014). Here, if there is a synergistic effect, it remains weak, and the GRA model provide overall a reasonable approximation of the behavior of the mixtures.

3.3 Mixtures at observed concentration ratios in food samples

Twelve mixtures of 2–5 compounds of phytoestrogens, phthalates, food additives and plant protection products were prepared according to concentration ratios found in real food samples, and their estrogenic activities were tested using ERα-CALUX. Concentration ratios were determined by LC/GC-MS/MS. as part of the ENDFOODTOX project funded by FPS Public Health (RF 18/6326). More information about how the food samples were selected can be found in the Supplementary Material section. The concentration-response curves obtained for the 12 mixtures and the corresponding model predictions are shown in Figure 3. Overall, the GRA model gave better results for these mixtures than the CA model (Supplementary Table S3). The observed experimental results are within the 95% CI of the prediction model. Consequently, no statistically significant differences between observation and prediction could be reported. There were some deviations from the expected mean activity with mixtures #1, #2, and #3 (Daidzein, Genistein) showing slightly below-average activities and mixtures #5 (Daidzein, Genistein, 6-PN, Kaempferol) and #6 (Daidzein, Genistein, 6-PN, Enterolactone) showing slightly above-average activities. Mixtures #10 (Daidzein, Genistein, 6-PN, Enterolactone, BHT) and #11 (Daidzein, Genistein, DIBP, Benzophenone, Triadimenol) show observations that fall on the leftmost edge of the 95 CI %, which could suggest the possibility of synergy in these mixtures.

FIGURE 3

FIGURE 3. Luciferase activity in VM7Luc4E cells for Mix #1 to #6 (see Supplementary Table S3). Luciferase activity in VM7Luc4E cells for Mix #7 to #12 (see Supplementary Table S3).

However, the experimental uncertainties (particularly for #9), do not allow us to conclude that there is a significant difference amongst the observed results. Triadimenol is present only in #11, while #10 and #12 are the only ones to contain BHT. Mixture #12, however, contains a lower BHT concentration, which might suggest that the potentially increased potency compared to the expected activity could be due to the antioxidant BHT in #10 and the fungicide triadimenol in #11. Synergistic effects with fungicides have already been described (Rondeau and Raine, 2022). This has been documented for triadimenol, which exhibits synergy effects when used in combination with other triazole-based fungicides (Li et al., 2018). For BHT, no description of synergistic activity in mixtures has been reported in the literature. Overall, our results show that the additive model is a sufficiently robust default model to describe mixture effects. Similarly, Yu et al. (2019) determinized the estrogenic activity of environmental EDs using the yeast estrogen screening test (YES) which suggested additive effects best described the observed results. Furthermore, more recent work on the combined toxicity of EDs revealed that the CA approach provided reliable estimates for describing joint effects when the compounds tested shared a common mode of action (Hamid et al., 2021).

The CA model assumes, however, that all the chemicals in the mixture behave as a dilution of each other. Mathematically, this implies that the concentration-response curves for each compound follow a parallel-line logistic regression model (Gottschalk and Dunn, 2005). It also means that the Hill slope is not significantly different from 1 and that all the concentration-response curves of the compounds in the mixture have equivalent background and maximum responses (Cedergreen et al., 2008). These assumptions do not apply to our CALUX measurements, where the Hill coefficient may differ from 1 and the maximum responses may fluctuate greatly, particularly in the case of phytoestrogens (Supplementary Table S1). Daidzein and genistein, for example, are supra-inductive agonists (180%–200% efficacy), while kaempferol (60% efficacy) is a partial agonist. An additive model was therefore adapted, to simultaneously consider variations in the efficacy, potency, and sensitivity of the different compounds. Monte Carlo simulations were then used to test the robustness of the additivity hypothesis. This made it possible to reconstruct concentration response curves for highly variable experimental conditions and to accurately predict the activity of the mixtures in most cases. The GRA approach therefore appears to be a suitable complement to the CA approach and can be used when the assumptions of parallel response curves are not met.

4 Conclusion

Overall, the GRA model developed as part of this research successfully predicts the effects of the mixtures, with the % of variance explained (RSQ-values) ranging from 84% to 99%. Furthermore, GRA outperforms the CA model applied by default (% variance explained 23%–98%) in describing the responses of non-equipotent mixtures of partial and complete agonists. The assumption that each component of the mixture can be considered as a dilution of the other is therefore not validated. This does not necessarily exclude the CA hypothesis, nor does it validate the IA ones. Recombinant cells used in the ERα-CALUX bioassay contain a stably transfected firefly luciferase reporter gene that responds to chemicals that can bind to the ER and activate it, resulting in the induction of luciferase gene expression. This induction occurs in a time-, concentration- and chemical-specific manner, which could lead to dissimilarities depending on the mixtures of agonists tested (Rogers and Denison, 2000; Besselink, 2015).

These results highlight the importance of applying an additive model for xenoestrogen risk assessment, but that care should be taken at low doses, as mixtures may exhibit joint behaviour if the chemicals are combined below their individual level with no observable adverse effects (NOAEL). Additionally, the compounds Triadimenol and BHT should be studied more intensely, as they show potential signs of synergy when tested in combination. Taken together, these results also underscore the need to broaden the approach to ED risk assessment, going beyond the treatment of individual compounds, for example, using adverse reaction pathways (AOPs) as tools to support cumulative risk assessment of co-exposure in food and food contact materials.

Data availability statement

The datasets presented in this study can be found in online repositories. The names of the repository/repositories and accession number(s) can be found in the article/Supplementary Material.

Ethics statement

Ethical approval was not required for the studies on animals in accordance with the local legislation and institutional requirements because only commercially available established cell lines were used.

Author contributions

ME made substantial contributions to the design of the work, analysis of data, and co-writing and editing the manuscript. IB made substantial contributions to conducting the analyses of samples and data collection. SE contributed to the analysis of data collected and writing and editing the manuscript. All authors contributed to the article and approved the submitted version.

Funding

This work was funded by the VUB Strategic Research Program (SRP68) and was supported by a Federal Public Service, Health, Food, Chain Safety and Environment grant to ME (RF 18/6326).

Acknowledgments

The authors would like to thank Mélanie Di Mario for providing the concentration ratios analysed in food samples as part of the ENDFOODTOX project (RF 18/6326).

Conflict of interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

The handling editor RD declared a past co-authorship with the authors ME and IB.

Publisher’s note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Supplementary material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/ftox.2023.1252847/full#supplementary-material

References

Abeliovich, H. (2005). An empirical extremum principle for the Hill coefficient in ligand-protein interactions showing negative cooperativity. Biophysical J. 89, 76–79. doi:10.1529/biophysj.105.060194

PubMed Abstract | CrossRef Full Text | Google Scholar

Ahn, C., and Jeung, E. B. (2023). Endocrine-Disrupting chemicals and disease endpoints. Int. J. Mol. Sci. 24 (6), 5342. doi:10.3390/ijms24065342

PubMed Abstract | CrossRef Full Text | Google Scholar

Bel, S., Van den Abeele, S., Lebacq, T., Ost, C., Brocatus, L., Stievenart, C., et al. (2016). Protocol of the Belgian food consumption survey 2014: objectives, design and methods. Arch. Public Health 74, 20. doi:10.1186/s13690-016-0131-2

PubMed Abstract | CrossRef Full Text | Google Scholar

Besselink, H. T. (2015). Validation of the (anti-) ERα CALUX bioassay. U2-OS cells Transcriptional ERalpha CALUX-assay for the detection of estrogenic and anti-estrogenic chemicals for inclusion in TG455. OECD ERα CALUX validation – v5.

Google Scholar

Brennan, J. C., Bassal, A., He, G. C., and Denison, M. S. (2016). Development of a recombinant human ovarian (bg1) cell line containing estrogen receptor alpha and beta for improved detection of estrogenic/antiestrogenic chemicals. Environ. Toxicol. Chem. 35 (1), 91–100. doi:10.1002/etc.3146

PubMed Abstract | CrossRef Full Text | Google Scholar

Brouwers, M. M., van Tongeren, M., Hirst, A. A., Bretveld, R. W., and Roeleveld, N. (2009). Occupational exposure to potential endocrine disruptors: further development of a job exposure matrix. Occup. Environ. Med. 66 (9), 607–614. doi:10.1136/oem.2008.042184

PubMed Abstract | CrossRef Full Text | Google Scholar

Cedergreen, N., Christensen, A. M., Kamper, A., Kudsk, P. S., Mathiassen, K., Streibig, C., et al. (2008). A review of independent action compared to concentration addition as reference models for mixtures of compounds with different molecular target sites. Environ. Toxicol. Chem. 27, 1621–1632. doi:10.1897/07-474.1

PubMed Abstract | CrossRef Full Text | Google Scholar

Cedergreen, N., Dalhoff, K., Li, D., Gottardi, M., and Kretschmann, A. C. (2017). Can toxicokinetic and toxicodynamic modeling Be used to understand and predict synergistic interactions between chemicals? Environ. Sci. Technol. 51 (24), 14379–14389. doi:10.1021/acs.est.7b02723

PubMed Abstract | CrossRef Full Text | Google Scholar

Cedergreen, N. (2014). Quantifying synergy: a systematic review of mixture toxicity studies within environmental toxicology. PLoS One 9 (5), e96580. doi:10.1371/journal.pone.0096580

PubMed Abstract | CrossRef Full Text | Google Scholar

Cox, L. A., Ketelslegers, H. B., and Lewis, R. J. (2021). The shape of low-concentration dose-response functions for benzene: implications for human health risk assessment. Crit. Rev. Toxicol. 51 (2), 95–116. doi:10.1080/10408444.2020.1860903

PubMed Abstract | CrossRef Full Text | Google Scholar

Darbre, P. D. (2022). “Exposure to mixtures of EDCs and long-term effects,” in Endocrine disruption and human health, 165–182. doi:10.1016/B978-0-12-821985-0.00020-7

CrossRef Full Text | Google Scholar

Diamanti-Kandarakis, E., Bourguignon, J. P., Giudice, L. C., Hauser, R., Prins, G. S., Soto, A. M., et al. (2009). Endocrine-disrupting chemicals: an Endocrine Society scientific statement. Endocr. Rev. 30 (4), 293–342. doi:10.1210/er.2009-0002

PubMed Abstract | CrossRef Full Text | Google Scholar

Elskens, M., Baston, D. S., Stumpf, C., Haedrich, J., Keupers, I., Croes, K., et al. (2011). CALUX measurements: statistical inferences for the dose-response curve. Talanta 85 (4), 1966–1973. doi:10.1016/j.talanta.2011.07.014

PubMed Abstract | CrossRef Full Text | Google Scholar

Elskens, M., Van Langenhove, K., Carbonnel, V., Brion, N., and Eisenreich, S. J. (2023). Dynamics of estrogenic activity in an urban river receiving wastewater effluents: effect-based measurements with CALUX. Water Emerg. Contam. Nanoplastics 2 (2), 9. doi:10.20517/wecn.2023.15

CrossRef Full Text | Google Scholar

Evans, R. M., Scholze, M., and Kortenkamp, A. (2012). Additive mixture effects of estrogenic chemicals in human cell-based assays can be influenced by inclusion of chemicals with differing effect profiles. PLoS One. 7 (8), e43606. doi:10.1371/journal.pone.0043606

PubMed Abstract | CrossRef Full Text | Google Scholar

Gottschalk, P. G., and Dunn, J. R. (2005). Measuring parallelism, linearity, and relative potency in bioassay and immunoassay data. J. Biopharm. Stat. 15 (3), 437–463. doi:10.1081/BIP-200056532

PubMed Abstract | CrossRef Full Text | Google Scholar

Guarnotta, V., Amodei, R., Frasca, F., Aversa, A., and Giordano, C. (2022). Impact of chemical endocrine disruptors and hormone modulators on the endocrine system. Int. J. Mol. Sci. 23 (10), 5710. doi:10.3390/ijms23105710

PubMed Abstract | CrossRef Full Text | Google Scholar

Hamid, N., Junaid, M., and Pei, D. S. (2021). Combined toxicity of endocrine-disrupting chemicals: a review. Ecotoxicol. Environ. Saf. 215, 112136. doi:10.1016/j.ecoenv.2021.112136

PubMed Abstract | CrossRef Full Text | Google Scholar

Heys, K. A., Shore, R. F., Pereira, M. G., Jones, K. C., and Martin, F. L. (2016). Risk assessment of environmental mixture effects. RSC Adv. 6 (53), 47844–47857. doi:10.1039/c6ra05406d

CrossRef Full Text | Google Scholar

Hornung, R. W., and Lanphear, B. P. (2014). The supralinear dose-response for environmental toxicants: a statistical artifact? Clin. Toxicol. (Phila) 52 (2), 88–90. doi:10.3109/15563650.2013.878946

PubMed Abstract | CrossRef Full Text | Google Scholar

Hsieh, N. H., Chen, Z., Rusyn, I., and Chiu, W. A. (2021). Risk characterization and probabilistic concentration–response modeling of complex environmental mixtures using new approach methodologies (NAMs) data from organotypic in vitro human stem cell assays. Environmental Health Perspectives. doi:10.1289/EHP7600

CrossRef Full Text | Google Scholar

Imparato, C., Bifulco, A., Silvestri, B., and Vitiello, G. (2022). Recent advances in endocrine disrupting compounds degradation through metal oxide-based nanomaterials. Catalysts 12 (3), 289. doi:10.3390/catal12030289

CrossRef Full Text | Google Scholar

Indurthi, D. C., and Auerbach, A. (2021). Agonist efficiency from concentration-response curves: structural implications and applications. Biophys. J. 120 (9), 1800–1813. doi:10.1016/j.bpj.2021.02.034

PubMed Abstract | CrossRef Full Text | Google Scholar

Kumar, M., Sarma, D. K., Shubham, S., Kumawat, M., Verma, V., Prakash, A., et al. (2020). Environmental endocrine-disrupting chemical exposure: role in non-communicable diseases. Front. Public Health 8, 553850. doi:10.3389/fpubh.2020.553850

PubMed Abstract | CrossRef Full Text | Google Scholar

Langmead, C. J., and Christopoulos, A. (2013). Supra-physiological efficacy at GPCRs: superstition or super agonists? Br. J. Pharmacol. 169, 353–356. doi:10.1111/bph.12142

PubMed Abstract | CrossRef Full Text | Google Scholar

Li, J., Liu, H., Guo, Z., Yang, M., Song, J., and Ma, H. (2018). Two Cu(ii)-triadimenol complexes as potential fungicides: synergistic actions and DFT calculations. RSC Adv. 8 (6), 2933–2940. doi:10.1039/c7ra10572j

PubMed Abstract | CrossRef Full Text | Google Scholar

Marlatt, V. L., Bayen, S., Castaneda-Cortes, D., Delbes, G., Grigorova, P., Langlois, V. S., et al. (2022). Impacts of endocrine disrupting chemicals on reproduction in wildlife and humans. Environ. Res. 208, 112584. doi:10.1016/j.envres.2021.112584

PubMed Abstract | CrossRef Full Text | Google Scholar

Martin, H. L., Svendsen, C., Lister, L. J., Gomez-Eyles, J. L., and Spurgeon, D. J. (2009). Measurement and modeling of the toxicity of binary mixtures in the nematode Caenorhabditis elegans - a test of independent action. Environ. Toxicol. Chem. Int. J. 28 (1), 97–104. doi:10.1897/07-215.1

PubMed Abstract | CrossRef Full Text | Google Scholar

Martin, O., Scholze, M., Ermler, S., McPhie, J., Bopp, S. K., KienzlerKortenkamp, A. A., et al. (2021). Ten years of research on synergisms and antagonisms in chemical mixtures: a systematic review and quantitative reappraisal of mixture studies. Environ. Int. 146, 106206. doi:10.1016/j.envint.2020.106206

PubMed Abstract | CrossRef Full Text | Google Scholar

Medlock Kakaley, E., Cardon, M. C., Gray, L. E., Hartig, P. C., and Wilson, V. S. (2019). Generalized concentration addition model predicts glucocorticoid activity bioassay responses to environmentally detected receptor-ligand mixtures. Toxicol. Sci. 168 (1), 252–263. doi:10.1093/toxsci/kfy290

PubMed Abstract | CrossRef Full Text | Google Scholar

Mukherjee, R., Pandya, P., Baxi, D., and Ramachandran, A. V. (2021). Endocrine Disruptors-'Food' for thought. Proc. Zool. Soc. 74 (4), 432–442. doi:10.1007/s12595-021-00414-1

PubMed Abstract | CrossRef Full Text | Google Scholar

Ong, H. T., Samsudin, H., and Soto-Valdez, H. (2020). Migration of endocrine-disrupting chemicals into food from plastic packaging materials: an overview of chemical risk assessment, techniques to monitor migration, and international regulations. Crit. Rev. Food Sci. Nutr. 62, 957–979. doi:10.1080/10408398.2020.1830747

PubMed Abstract | CrossRef Full Text | Google Scholar

Perkins, E. J., Ashauer, R., Burgoon, L., Conolly, R., Landesmann, B., Mackay, C., et al. (2019). Building and Applying Quantitative Adverse Outcome Pathway Models for Chemical Hazard and Risk Assessment. Environ. Toxicol. Chem. 38 (09), 1850–1865. doi:10.1002/etc.4505

PubMed Abstract | CrossRef Full Text | Google Scholar

Rogers, J. M., and Denison, M. S. (2000). Recombinant cell bioassays for endocrine disruptors: development of a stably transfected human ovarian cell line for the detection of estrogenic and anti-estrogenic chemicals. Vitr Mol Toxicol 13 (1), 67–82.

PubMed Abstract | Google Scholar

Rondeau, S., and Raine, N. E. (2022). Fungicides and bees: a review of exposure and risk. Environ. Int. 165, 107311. doi:10.1016/j.envint.2022.107311

PubMed Abstract | CrossRef Full Text | Google Scholar

EFSA Silano, V., Barat Baviera, J. M., Bolognesi, C., Chesson, A., Cocconcelli, P. S., et al. (2019). Update of the risk assessment of di-butylphthalate (DBP), butyl-benzyl-phthalate (BBP), bis(2-ethylhexyl)phthalate (DEHP), di-isononylphthalate (DINP) and di-isodecylphthalate (DIDP) for use in food contact materials. EFSA J. 2019, e05838. doi:10.2903/j.efsa.2019.5838

PubMed Abstract | CrossRef Full Text | Google Scholar

Silva, E., Rajapakse, N., and Kortenkamp, A. (2002). Something from “nothing” - eight weak estrogenic chemicals combined at concentrations below NOECs produce significant mixture effects. Environ. Sci. Technol. 36, 1751–1756. doi:10.1021/es0101227

PubMed Abstract | CrossRef Full Text | Google Scholar

Trasande, L., Zoeller, R. T., Hass, U., Kortenkamp, A., Grandjean, P., Myers, J. P., et al. (2016). Burden of disease and costs of exposure to endocrine disrupting chemicals in the European Union: an updated analysis. Andrology 4 (4), 565–572. doi:10.1111/andr.12178

PubMed Abstract | CrossRef Full Text | Google Scholar

Wiesinger, H., Wang, Z., and Hellweg, S. (2021). Deep dive into plastic monomers, additives, and processing aids. Environ. Sci. Technol. 55 (13), 9339–9351. doi:10.1021/acs.est.1c00976

PubMed Abstract | CrossRef Full Text | Google Scholar

Yu, H., Caldwell, D. J., and Suri, R. P. (2019). In vitro estrogenic activity of representative endocrine disrupting chemicals mixtures at environmentally relevant concentrations. Chemosphere 215, 396–403. doi:10.1016/j.chemosphere.2018.10.067

PubMed Abstract | CrossRef Full Text | Google Scholar