Exposure of the Irish population to PBDEs in food: consideration of parameter uncertainty and variability for risk assessment

Polybrominated diphenyl ethers (PBDEs) are brominated flame retardants used to retard the ignition and/or spread of fire. PBDEs are used in various consumer products, such as textiles, mattresses and TV screens. This study presents a chemical risk assessment for the Irish population based on exposure to PBDEs from food. Special regard is given to the influence of parameter uncertainty and variability on the margins of safety. To quantitatively model uncertainty and variability in concentration data and variability in consumer behavior, a hierarchical probabilistic model was constructed. This model was evaluated using a two-dimensional Monte Carlo simulation (2D-MCS) approach. By considering uncertainty and variability in concentration data, margins of safety (MOS) were derived that are lower by a factor of ∼2 compared to MOS based on dose estimates that only consider variability. The lowest MOS is 7.5 × 104 for BDE-99, with impaired spermatogenesis as toxic endpoint. Assuming an MOS of 104 as acceptable, we conclude that there is no significant risk for human health through intake of contaminated food. To investigate whether additional measurements could improve the quality of dose estimates, the statistic “uncertainty-to-variability (UVR)” was developed. By applying the UVR to our dose estimates, we show that, in our case, the datasets contain little uncertainty and additional measurements would not significantly improve the quality of dose estimates.


Introduction
Polybrominated diphenyl ethers (PBDEs) are brominated flame retardants that are added to household products, such as textiles, computers or television screens, to inhibit ignition and/or to retard the spread of fire. PBDEs have been produced in three different technical mixtures, pentabromo diphenyl ether (pentaBDE), octabromo diphenyl ether (octaBDE) and decabromo diphenyl ether (decaBDE) (De Wit 2002). In the European Union, the use of pentaBDE and octaBDE is prohibited (European Commission 2003), and the use of decaBDE is no longer allowed in electric and electronic (E&E) equipment (BSEF 2009). However, scientific discussion about the toxicological effects of decaBDE is ongoing (Hardy et al. 2009). Although no reference doses have been established by a governmental body so far, animal studies show toxicological effects for PBDEs such as endocrine disruption, neuro-developmental and behavioral alterations, as well as hepatic abnormalities and possibly cancer (McDonald 2002;Birnbaum and Staskal 2004;Darnerud 2008). In epidemiological studies, effects on male reproductive hormones and fertility, cryptorchidism, lower birth weight and thyroid hormone homeostasis have been found (Chao et al. 2007;Main et al. 2007;Akutsu et al. 2008;Turyk et al. 2008;Meeker et al. 2009).
Previous studies have shown that, beside dust, food intake is an important pathway for human exposure to PBDEs (Jones-Otazo et al. 2005;Darnerud et al. 2006;Lorber 2008;Fraser et al. 2009). Hence, it is of interest to accurately model the exposure of a human population to PBDEs in food and to provide ranges for the estimated doses. This can be achieved with probabilistic models that are used to generate distributions of the dose estimates based on distributions of the model input parameters (IPCS 2008). Probabilistic models make it possible to assess the percentage of a population that might be at risk by comparing the dose distributions, which reflect the variability among individuals, to toxicological thresholds of the respective substance. They also offer the opportunity to propagate parameter uncertainty through the model (Cullen and Frey 1999).
Exposure and risk assessments are often based on small datasets. These datasets reflect the variability of model parameters, such as PBDE concentrations in food samples and, additionally, are subject to measurement error. Even if the measurement error is small, the finite sample size of datasets leads to uncertainty in the estimated parameters of the variability distribution (Frey and Rhodes 1998;Frey and Burmaster 1999). The distinction of uncertainty and variability is important, because variability cannot be reduced by further measurements or better analytical equipment, whereas in cases of high uncertainty an exposure assessment would benefit from an extension of the dataset by adding reliable data points. For this reason, we want to assess how the uncertainty and variability in PBDE concentrations as well as the variability in consumer behavior influence the estimated doses.
An input variable of a probabilistic model is usually characterised by a parameterised probability distribution (e.g. a Gaussian distribution can be characterised by the parameters mean and standard deviation). Consequently, in the following, we distinguish model variables (e.g. concentrations) and distribution parameters. A mathematical model formulating the variability distribution with uncertain distribution parameters takes the form of a so-called hierarchical model, because these distribution parameters are not represented by precise numbers, but again by probability distributions. The parameters of these distributions are called hyperparameters. For fixed values of these parameters, the conditional probability distribution of the results represents the effects of the variability in the datasets on the results. The full, unconditional distribution considers also the uncertainty in the distribution parameters and represents the combined effect of uncertainty and variability on the results.
Numerically, estimates of the uncertainty of distribution parameters can in special cases be obtained from data by analytical estimators, in other cases by linearisation or by using bootstrapping techniques. Bootstrapping techniques are the most versatile approaches to characterise statistical properties of an estimator and are based on resampling (Efron 1982). Estimates of characteristics of conditional and unconditional distributions of doses can be obtained with a two-dimensional Monte Carlo simulation (2D-MCS) (Cullen and Frey 1999;Frey and Rhodes 1998).
To date, there are only few exposure studies that have attempted to consider parameter uncertainty and separate it from population variability (Ragas et al. 2009) and, to the best of our knowledge, this is the first study using this approach for human exposure to PBDEs from food. We used 2D-MCS to evaluate the extensive Irish database on PBDE concentrations in food (Fernandes et al. 2009) along with food consumption data for the Irish adult population (IUNA 2001). After accounting for uncertainty and variability in PBDE concentrations in food, we set up a statistic called uncertainty-to-variability ratio (UVR), to serve as a quality measure for the datasets used for dose estimation. Subsequently, we analysed the influence of parameter uncertainty on risk assessment.

Materials and methods Modeling
Calculation of dose estimates The calculations were performed for nine important PBDE congeners: . Only foods of animal origin and vegetable oil were modeled, because all data points from other sample matrices (such as fruits) were below the analytical limit of detection (LOD).
Assuming chronic exposure conditions, the doses (representing daily intake of a yearly average) for food exposure to PBDEs were calculated according to the following equation: with D ki (ng/kg bw /day) as intake dose rate (i.e. before uptake into the body and herein referred to as dose) for congener k and individual i, C kj (ng/g fat ) for the concentration of congener k in the sample matrix (generally fat) of food group j of r ¼ 47 food groups, f j (g fat /g food ) as fat fraction of food group j and a ji (g food / kg bw /day) as the amount of daily consumed food per kg bw for food group j and individual i. For fish, the fat fractions are equal to 1, because their concentrations are already given on fresh-weight basis (and not on fatweight basis as for the other food matrices). Three food groups (mayonnaise, mayonnaise light and salad cream) consist of two fat types and their concentrations are, therefore, a composition of two concentrations measured in two different sample matrices multiplied by the according fat fractions (for details see Table S1). The dose of individual i for AEBDE is then calculated by Equation (2): with q ¼ 9 congeners. The food consumption was modeled using an individual-based approach. The daily averaged food consumption pattern was available for 468 male and 481 female adults (IUNA 2001). The distribution of food consumption patterns in the whole female and male adult population, respectively, was approximated by random re-sampling from these two groups of persons. This procedure likely overestimates the variance in food intake rates (especially for foods that are not consumed on a regular basis, such as ice cream), because the consumption of each person was only recorded for a period of 7 days. The available methods to adjust the questionnaire data to a so-called ''usual intake'' (Nusser et al. 1996;Dodd et al. 2006), however, have been developed for parametric models only and could, therefore, not be applied in this study.
The variation of fat fractions and concentrations was modeled by lognormal distributions (see also the section on PBDE concentration data). The variation in concentrations in food is caused by different factors, such as origin, manufacturing, or age of fish, for example. As it is unclear whether individuals mainly consume food items with similar concentrations (as probably a fisherman would do when eating his own fish) or consume food of an average over different concentrations (as probably a consumer in a city would do by eating different fish species of different origin, manufacturing, age, etc.), we calculate the results for two extreme scenarios. In the first scenario, the so-called high-variability scenario, we assume that an individual always consumes a food item of the same concentration. As discussed above, the concentrations used to calculate the dose of different individuals are assumed to be distributed according to the estimated concentration distribution, but for a single individual the concentration in this food item is randomly drawn and then constant. This leads to a high variability in the dose between consumers, as the variability is caused by variability in food intake and variability in concentration. In the second, so-called low-variability scenario, it is assumed that each individual averages over the concentration distribution of each food type by randomly buying food of different origin, manufacturing, etc. Consequently, we use the mean of the concentration distribution for all individuals, but still consider the uncertainty of this mean. This leads to a lower variability among individuals as the variability in this scenario is only caused by differences in food intake. Figure 1 shows the probabilistic model used to describe uncertainty and variability in all model parameters that are required to calculate the doses by Equations (1) and (2). The model is based on the following statistical assumptions:

Consideration of parameter uncertainty and variability
1. The concentrations C kj of the different congeners k in the fat fractions f j of all food groups j are independent of each other. The congener and food item specific distributions describing the variability of concentrations can be described by a lognormal distribution with parameters p kj ¼ ð C,kj , C,kj Þ, where C,kj represents the mean and C,kj the standard deviation of the concentrations. 2. High-variability (high-var) scenario: The variability in the concentrations C kj is assumed to characterise the variability of the concentrations in the foods consumed by different consumers on a chronic basis. This means that we assume that one consumer always eats food items with the same concentration (randomly drawn from the concentration distribution of this food item). If not indicated otherwise, results presented and discussed in  (1) and (2)). The fat lines illustrate how the dose for one congener and one food group is calculated. The model is operated on a distribution of individuals approximated by the empirical distribution of a sub-sample of the Irish population.
this study are derived with assumptions of the high-var scenario. Low-variability (low-var) scenario: To explore the influence of assigning one concentration pattern to each individual on the uncertainty/variability ratio, we constructed a second scenario where the concentrations of food items with different concentrations are averaged for each individual. Therefore, for each food item, we use the arithmetic mean of the concentrations. 3. The uncertainty of the parameters p kj is characterised by independent lognormal marginal distributions with hyperparameters determined from the fit of the distribution to measured data. 4. The variability in fat fractions f j is assumed to be characterised by independent lognormal distributions as described in the section on PBDE concentration data. 5. The consumption parameters a ji are modeled by the empirical distribution of the data for the population sample, as described in the section on consumption data.
These assumptions are motivated as follows. By assumption 1, uncertainty due to small datasets was only considered for the food concentration data. As the Irish datasets for food consumption and fat fractions in different food items are comparably large, they are assumed to have a negligible random sampling error and, hence, should reflect variability only. Other uncertainties, e.g. those associated with analytical measurements or the collection methodology for food consumption data, are assumed to be small and are, therefore, not taken into account. Although some dependence in the congener pattern of the PBDE concentrations measured in foods might exist due to the use of technical mixtures (e.g. pentaBDE) with a fixed congener pattern, assumption 1 is justified by the unknown primary source of the PBDEs in food and the transformation processes the different congeners are subject to. Assumption 2 (for high-var scenario) is based on the concept that the same consumer always consumes food with the same concentration during the modelled time period. This is a conservative assumption, as it likely overestimates the true variability of PBDE concentrations found in food on a chronic basis. Therefore, we relaxed assumption 2 for a second, less conservative low-var scenario by modelling the concentrations in fat matrices and fish with their arithmetic means and according standard errors. The fat fractions are modeled as point estimates (arithmetic mean of available fat fractions). As mentioned before, the high-var scenario likely overestimates the variability of the PBDE doses. This leads to higher values of large quantiles of the variability distribution, thus generating more conservative results in the risk assessment of food exposure to PBDEs. As this risk assessment is a major goal of our study, the focus in the main text is set on the high-var scenario. Detailed results of the low-var scenario are given in the Supporting Information (Section 8).
Assumption 3 was made for computational convenience and because it has been shown that lognormal distributions reasonably well characterise the data constructed by bootstrapping. Assumptions 4 and 5 are not critical, as they are justified by large datasets. Figure 2 illustrates for the case of one food group, one congener, and one individual how parametric bootstrapping and 2D-MCS were used in this study. Realizations of random variables representing the different parameters during the 2D-MCS are denoted with superscripts in brackets. The bootstrap simulations were performed in R (R Development Core Team 2009) and the 2D-MCS in Crystal Ball TM (Oracle 2009a).
In a first step (1), data were prepared for the bootstrap simulations, i.e. nondetects were replaced by extrapolated values with a procedure described in more detail in the section on PBDE concentration data. Second (2), we fitted a lognormal distribution to the concentration data with w data points (if at least one data point was above the LOD in the original dataset, for details see the section on PBDE concentration data). Third (3), we generated w random data points from this distribution to form one bootstrap sample. Fourth (4), the mean and standard deviation (SD) of the log-transformed bootstrap sample were calculated.
Steps (3) and (4) were repeated B ¼ 1000 times, so that altogether 1000 bootstrap samples were generated and their individual means and SDs determined. In a fifth step (5), the means and SDs of the B calculated means and SDs were determined and used to parameterise the two lognormal distributions of the concentration parameters. With the Anderson-Darling goodness-offit statistic (Stephens 1974), we confirmed that the parameters are best modeled with lognormal distributions. Dependencies between parameters were not modeled, as we assume that they have only a minor influence on the result of our simulation (Frey and Rhodes 1998). In a sixth step (6), we randomly generated one realisation of each parameter to parameterise the lognormal concentration distribution (within Crystal Ball TM these realisations were converted to log-scale to correctly parameterise the lognormal distribution of the concentrations) and performed a classical 1D-MCS (steps 7-9) that is called an inner loop in 2D-MCS. In step (9), Equation (1) is applied to calculate a dose estimate with the respective realisations of the concentration, the fat fraction distribution and the food intake distribution. The inner loop is repeated m times and steps 6-9 (outer loop) n times: In this study, the 2D-MCS was performed with m ¼ 1000 for the inner loop and n ¼ 100 runs for the outer loop. Thus, we get n samples from the dose distribution conditional on the n values of the parameters. These were arranged in a m Â n sample matrix D instead of only a vector of sample values in 1D-MCS. All samples were generated using Latin hypercube sampling to increase the quality of the results in the tails of the distribution (Cullen and Frey 1999).
Uncertainty-to-variability ratio (UVR) To quantify the relation between uncertainty and variability a statistic called uncertainty-to-variability ratio (UVR) is defined: with U as a measure for uncertainty and V as a measure for variability. To put U and V on a sound conceptual basis, they are derived from the law of total variance (Weiss 2006). The total variance of the output, i.e. the distribution of dose D (for clarity of presentation, we omit indices in this section): can be decomposed into the expected value of the variance of the conditional distribution of D given p, D|p, and the variance of the expected value of D|p.
If p denotes the parameters of the concentration distribution, this decomposition can be interpreted as a decomposition of the variance of the dose distribution into the contributions due to variability and uncertainty. The term takes the variance of D for each given p and averages these variances over the marginal distribution of p. It thus represents the average contribution to the variance of the result due to variability. The term averages over D for each given p and then takes the variance of these outputs with respect to the uncertainty in p. As the distribution of p represents the uncertain knowledge about the concentration distribution, this expression represents the contribution of uncertainty to the variance of dose distribution D. We have thus split the total variance of the result into the contributions due to variability and uncertainty. MCS offers straightforward ways to numerically estimate U and V. The expected value of D|p in U could be estimated by applying a smoothing algorithm to the sample points of D as a function of p. U could then be calculated by averaging these values over the marginal distribution of p and V as the difference Var[D]-U. The 2D Monte Carlo approach chosen for this paper makes it even easier to get estimates of U and V. As the columns of our sample matrix D represent samples of the conditional distribution D|p for samples of the marginal distribution of p, we can directly take the average of their empirical variance to get an estimator:V Similarly, we get The redundancy Var(D) ¼ V þ U allows us to test if the empirical estimates sum up to the empirical estimate of the total variance. The estimator for the UVR used in this study is then: The UVR consists of measures of uncertainty and variability that are easy to interpret. The UVR is larger than 1 if the uncertainty outweighs variability and smaller than 1 in the opposite case. This statistic is designed to help evaluating datasets that are used in 2D-MCS. If the UVR is large (i.e. above 1), welldistributed additional data points will considerably improve the accuracy and precision of dose estimates. Instead, when the UVR is small (i.e. below 1), additional measurements will have only a small impact on the quality of dose estimates.

Sensitivity analysis
To determine the input parameters (in this study, these are model variables and their associated distribution parameters) with most influence on the output, the contribution to variance (CtV) was calculated (Oracle 2009b). The CtV is calculated by determining the Spearman rank correlation between the output and each of its input parameters. The correlation coefficients are then squared and normalised to total 100%. The higher the CtV the higher is the influence of the input parameter on the output (see Figures S12 and S13).

PBDE concentration data
Due to the lipophilic nature of PBDEs, they are mainly found in fats and oils. Consequently, in this study, only food of animal origin and vegetable oil was modelled. Almost all occurrence data of PBDEs in the investigated food groups stem from Ireland and were provided by the Food Safety Authority of Ireland (FSAI) (published in Fernandes et al. (2009) (Holland et al. 1989;Chan et al. 1994Chan et al. 1995Chan et al. 1996. For three food groups a combination of two sample matrices was needed to derive the fresh-weight concentration (for details, see Table S1 and Trudel et al. (2010)).
To prepare the concentration dataset of the sample matrices for bootstrapping, the robust regression on order statistics estimation (ROS) method was applied. This was to appropriately extrapolate values below LOD based on information gained from data points above LOD (Helsel 2005;Trudel et al. 2010). ROS assumes that the data are lognormally distributed. Trudel et al. (2010) could show that this assumption is reasonable for the Irish concentration data and, therefore, the concentrations were modeled with lognormal distributions. More details on how datasets with nondetects were treated are given in Section 2 of the Supporting Information.

Consumption data
Food consumption data for 949 persons were collected over a survey period of 7 consecutive days. The survey periods were distributed evenly over the whole year to reflect seasonal variation in food intake patterns. All consumption data are given on a per kilogram body weight basis (IUNA 2001) and are used directly in the modeling (see section on calculation of dose estimates).

Risk assessment
For risk assessments, it is important to have reliable dose estimates that can be compared to benchmarks or reference doses, e.g. a tolerable daily intake (TDI). To date, however, no intergovernmental expert body such as the European Food Safety Authority (EFSA) or the Joint FAO/WHO Expert Committee on Food Additives and Contaminants (JECFA) has determined TDIs for PBDEs, because there is as yet insufficient information available on their toxicity due to intake of contaminated food.
In the absence of TDIs for PBDEs, exposure can be compared to critical dose descriptors such as the NOAEL (no-observed-adverse-effect level) confirmed by expert bodies like the European Commission (EC), the UK Committee on Toxicity (COT) or the US EPA. Therefore, the doses derived in this study together with doses derived in a previous study (Trudel et al. 2010) are compared to the relevant NOAEL or LOAEL (lowest-observed-adverse-effect level) identified in the EU risk assessment reports by following the margin of safety (MOS) concept (van Leeuwen and Vermeire 2007). The risk assessment is performed for effects on the liver caused by the commercial mixtures of pentaBDE (product name: DE-71) and decaBDE, and for effects on reproduction caused by octaBDE: Individual congener data are compared to the relevant NOAEL, LOAEL or BMDL (benchmark dose lower confidence level) for potential effects on the liver, neurodevelopment, reproduction and impaired spermatogenesis as identified by a number of different researchers (for details, see Supporting Information, Section 4).
To judge whether the derived margins of safety are acceptable, the size of the margins needs to be considered. For chemicals with a toxicological threshold and a comprehensive toxicological database, health-based guidance values are usually derived by applying a default uncertainty factor (UF) of 100 to allow for interspecies differences (10) and human variability (10) (Dorne 2010). This factor of 100 applied to the threshold can be regarded as minimum MOS to ensure that the derived guidance value protects the health of exposed humans. Where no NOAEL has been identified, an additional UF of 3-10 is required for extrapolation from a LOAEL to a NOAEL. Further uncertainties, such as limitations in the database (e.g. study duration or data gaps) sometimes also need to be taken into account (up to another factor of 10). However, recently research efforts have been made to refine these UFs taking into account toxicokinetic and toxicodynamic aspects and replace them with pathway-related UFs, when metabolic routes are known, or, ideally, with chemical-specific adjustment factors (Dorne 2010;Dorne et al. 2005). In our case, chemical-specific or pathway-related adjustment factors are not available. Therefore, the combination of the default UFs suggests that a MOS of 41000 is acceptable if the point of departure (POD) for the risk assessment is a NOAEL or a BMDL, and a MOS of 3000-10,000 if the POD is a LOAEL.

Exposure assessment
The results of the 2D-MCS for AEBDE in the high-var scenario (Figure 3 and Figures S7 and S8) show that the spread of the distribution of dose D increases towards higher quantiles of the cumulative distribution function (CDF). From a risk assessment point of view, mainly the upper quantiles Q of D|p are of interest, therefore we show the 97.5th quantile and 99th quantile of D|p in Table 1 (for low-var scenario see Table S5). To also reflect the uncertainty in dose estimates the 2.5th, the 50th and the 97.5th quantile of the above-mentioned quantiles are given. This is denoted with Q [Q {D|p}] with ¼ {2.5, 50, 97.5} representing uncertainty and ¼ {97.5, 99} representing variability. For every quantile, BDE-47 contributes the highest doses followed by -99 and -209. Together, these three congeners contribute more than two thirds to the total dose. BDE-49 and -100 also contribute considerably to total doses (about a fifth, independent of the quantile). Of minor importance are BDE-28, -153 and -154, which contribute together about 10% to the total dose (for all quantiles). The doses for Q 97.5 [Q {D|p}] with ¼ {97.5th quantile and 99th quantile} are on average 1.5-fold higher than the medians with BDE-209 having the largest difference (97.5th quantile about two times higher than the median). To better illustrate the distribution of Q 97.5 {D|p}, we show the estimated probability density in Figure 4 for female and male adults for AEBDE (for information on congeners, see Figures S9-S11 in the Supporting Information). The above findings are also valid for dose estimates derived for the low-var scenario, except that the spread of the uncertainty distributions is smaller by a factor of up to 2 (for ¼ 97.5 and between 2.5 and 97.5).

UVRs for the Irish population
The results for the uncertainty, variability and the UVRs for the Irish dataset and AEBDE, BDE-47, -99, and -209 are presented in Table 2. The UVRs for AEBDE, -47, and -99 are quite similar (0.0032-0.0248), showing that the uncertainty is small compared to the population variability. The UVRs for BDE-209 are slightly higher (0.13 and 0.29), but still clearly below 1, indicating that also for exposure to BDE-209 variability is dominating uncertainty. The UVRs calculated for the low-var scenario are also all below 1 (Table S6).

Sensitivity analysis
The sensitivity analysis shows that dairy fat and lean fish are the most important food matrices for exposure to PBDEs (see Figures S12 and S13). The dominating congeners are BDE-47, -99, and -209. The model variables are generally more influential on the output than their distribution parameters (except for exposure to BDE-209). Table 3 shows the MOS derived for the individual congeners BDE-47, -99, -153, -183 and -209 and for the total dose of congeners BDE-47, -99, -100, -153 and -154, which are the major components of commercial penta mixtures (Konstantinov et al. 2008). All derived MOS for the 97.5th quantile of the doses as presented in Trudel et al. (2010) (for comparison we use the external doses presented therein, i.e. before uptake of the substance into the body) are larger than 1 Â 10 5 with the lowest MOS at 1.6 Â 10 5 for BDE-99 (endpoint impaired spermatogenesis). If also parameter uncertainty is considered in dose estimates (high quantile doses HQD), the safety margins are reduced by up to 60% with the lowest margin of safety identified for impaired spermatogenesis (reduced to 7.5 Â 10 4 ). Using BDE-99 dose estimates derived for the low-var scenario, the MOS equals 1.2 Â 10 5 and lies between the MOS based on 1D-MCS (1.6 Â 10 5 ) and high-var scenario (7.5 Â 10 4 ) dose estimates. Notes: a The penta-mix consists of BDE-47, -99, -100, -153, and -154 (note: this is not the sum of the single congeners due to nonlinearity of the model.) b Note: this is not the sum of the single congeners due to nonlinearity of the model.

Discussion
We show that the consideration of parameter uncertainty and variability can be influential for dose estimation and should, therefore, be considered in risk assessments if uncertainty in datasets is assumed to be large and dose estimates are close to reference values. This statistical concept has successfully been applied using a hierarchical, probabilistic model that was implemented numerically with parametric bootstrapping and 2D-MCS. Bootstrapping as used in this study, however, cannot reveal all the uncertainty contained in a (small) dataset, because the distribution fitted to the few given data points in step 2 of the applied modeling procedure (Figure 2) cannot exactly represent the ''true'' distribution of the examined model variable. Therefore, the less data points are available to perform bootstrapping and subsequent 2D-MCS, the more careful the derived dose estimates have to be interpreted. More research is also needed to improve the rather crude assumption that observed concentration variability either fully transfers to the calculation of chronic doses (high-var scenario) or not at all (lowvar scenario)). These assumptions might lead to an over-(high-var scenario) or underestimation (low-var scenario) of the variance in dose estimates. The highvar scenario might lead to an overestimation, because it is unlikely that an individual i will consume the same type of food at different points in time and that this food always possesses the same concentration. The low-var scenario might lead to an underestimation of the variance, because it is unlikely that concentrations   The penta-mix consists of BDE-47, -99, -100, -153, and -154 (note: this is not the sum of the single congeners due to nonlinearity of the model).
f Note: this is not the sum of the single congeners due to nonlinearity of the model. always average out on the long run and no variability can be observed. Hence, the true dose estimates likely lie between those of the high-and the low-var scenario. We need to keep in mind, however, that even in the low-var scenario the dose estimates are probably still on the conservative side, because the re-sampling technique used to model the food intake likely overestimates the variability in consumption patterns due to the short survey period of 7 consecutive days per investigated individual. The visual reflection of a small uncertainty in dose estimates for AEBDE compared to variability (Figure 3) is supported by the corresponding UVRs for AEBDE (Table 2). This finding is reasonable, as the most sensitive model variables for total PBDE exposure are concentrations of dairy products such as skimmed milk or cream, which are based on comparably large dairy fat datasets (see Table S2 for number of data points and Figures S12 and S13 for CtV). Further measurements will, therefore, not significantly reduce the uncertainty of dose estimates for AEBDE. The same is true for the dose estimates for the two influential congeners BDE-47 and -99 (for CtV see Figures S12 and S13). For BDE-209, however, the situation is slightly different as the UVR is as high as 0.1 and 0.3 for female and male adults, respectively. Although the UVRs are still clearly below one, it seems worthwhile to check whether additional measurements might reduce uncertainty. In particular, the dataset for salmon, which consists of three data points only, could be increased. However, as these three samples are composite samples and, hence, already better representations of the ''true'' concentration distribution of BDE-209 in salmon than single samples, it is likely that additional measurements will not lead to a large change in dose estimates. Additionally, the MOS for exposure to BDE-209 is comparably large (at least 2.6 Â 10 6 ), so that in this case more confidence in the exposure assessment is not relevant for the risk assessment.
The UVRs for the low-var scenario (Table S6) are also all below 1, even though the variability in concentration data and fat fractions is no longer included in V. This indicates that the remaining variability in food intake data is dominant and that also under the assumptions of the low-var scenario additional concentration measurements will not significantly improve the quality of dose estimates. Even if we consider that the variance in measured food intake data is possibly too high due to temporal variance (limited measurement period of 7 consecutive days), usual intakes (Dodd et al. 2006) will not alter the UVRs to values above one (the highest UVR is 0.095 for male adults and exposure to BDE-99). Therefore, our conclusion that variability dominates uncertainty is not a result of our model assumptions but reflects the good quality of the datasets.
Although the UVRs for all relevant congeners (Table 2) are below 1 in our study, we could show that the MOS based on considering uncertainty and variability (HQD in Table 3) are lower by about a factor of 2 compared to those MOS based on dose estimates that only consider variability (Q 97.5 in Table 3). Consequently, in the case of dose estimates close to a reference value or very small datasets (prone to contain significant uncertainty), 2D-MCS should be used to derive a factual conservative estimate. If the variability in relevant concentration data is large, it might be helpful to conduct specific surveys targeting those consumer groups that likely experience an increased risk (e.g. fish eaters) and to determine whether the variance in concentrations will average over time, i.e. to determine whether the low-or the high-var scenario better describes the true exposure situation.
For the actual risk assessment, we could show that all calculated MOS based on estimated high quantile doses are greater than 10,000 and, hence, exposure to PBDEs in food should not cause a risk to the health of Irish female and male adults. We need, however, to consider that the above exposure scenarios relate only to exposure of the Irish adult population to PBDEs in food. As it is known that younger consumers experience by trend higher doses than older consumers (Bakker et al. 2008;Lorber 2008) and might be more sensitive to certain toxic substances than adults (Bruckner 2000), it would be of great interest to also model the exposure of Irish infants, toddlers and children to PBDEs. Another aspect that needs to be considered is the influence of other pathways than food on the exposure of the Irish population. It is known from other studies (Jones-Otazo et al. 2005;Lorber 2008) that the inadvertent ingestion of house dust might considerably contribute to the total exposure to PBDE. It would be interesting to extend the present study to include this additional contribution, in particular, because uncertainty may be relevant for dust intake due to the lack of PBDE concentration data in dust for Ireland and the high variability and measurement uncertainty that are generally associated with dust ingestion rates.