Quantification of Conceptual Model Uncertainty in the Modeling of Wet Deposited Atmospheric Pollutants

Conceptual model uncertainty and parameter uncertainty are dominant contributors to the total uncertainty of a radioecological model output. In the present study the focus is on conceptual model uncertainty, which is often not acknowledged. Conceptual model uncertainty is assessed by subtracting from the total uncertainty of the model output the propagated parameter uncertainty, obtained by means of Bayesian inference analysis. The conceptual model uncertainty is quantified for two process‐based models, which describe the interception of wet deposited pollutants under equilibrium and kinetic conditions, respectively. The natural variability due the chemical valence of the elements considered is accounted for in both models. Quantitative evidence has been obtained that the conceptual model uncertainty can contribute to the total uncertainty budget of the models for interception of wet deposited pollutants at least as much as, if not more than, parameter uncertainty.


Environmental Models and Uncertainty of Their Output
Environmental models deliver outputs that are characterized by an inevitable degree of uncertainty.
Model parameters are constant quantities (Kirchner & Steiner, 2008) whereas model input variables are time-dependent quantities. These are usually measured during experimental activities.
Epistemic uncertainties may include: • parameter uncertainty, which arises from a limited knowledge of the values for the model parameters and the correlations between model parameters; • conceptual model uncertainty (sometimes also called structural or model uncertainty), which arises from an inadequate or oversimplified model structure that does not represent key characteristics of an ecosystem and key processes; • input uncertainty or the uncertainty on the values of the input variables, which arises principally from inaccuracies while measuring these alongside the process of concern. • scenario uncertainty, which arises from an incomplete knowledge of the scenario to be modeled, especially site-specific information and data; • modeler uncertainty, which arises from the (imperfect) way in which a modeler translates an assessment situation in a model (inadequate interpretation of the assessment situation, selection of an inappropriate model); • computational/numeric uncertainty, which arises from inaccuracy in numerical solutions.
Methods to evaluate conceptual model uncertainty are manifold and no standard recipe is available. However, many make use of the assumption that the residual between experimental output and modeled output quantifies the total uncertainty. For example, starting from probabilistic or Bayesian treatment of parameter uncertainty, the other con-tributions to the overall uncertainty can be obtained with decomposition of mean squared error as in Alderman & Stanfill (2017), Jin et al. (2010), Draper (1995), and Wallach et al. (2017). If several different conceptual models are available multimodel analysis is possible which requires high computational effort as in Symonds & Moussalli (2011), Walsh & Kaiser (2011), McAllister & Kirchner (2002.
In the context of radioecology, conceptual model uncertainty has not been quantified so far. The present study addresses this gap by considering a specific case study involving two process-based models developed for quantifying the process of interception of wet deposited pollutants including radionuclides.

Conceptual Model Uncertainty in Radioecology
Radioecological models have often a very simple mathematical structure, for example transfer factors or parametric equations (Urso, Hartmann, Diener, Steiner, & Vives i Batlle, 2015), that lacks of causal attributes and that represents in an aggregated form the processes under consideration. For example, conceptual model uncertainty in radioecology plays an important role when (Gonze & Sy, 2016), (Salbu, 2016): • empirical, highly sensitive parameters are used, for example the so-called concentration ratios and transfer factors. Concentration ratios quantify the transfer from one medium to another one, for example concentration ratio soil-plant, soil-mushroom, soil-cow milk, and so on. These factors can vary over orders of magnitude and are based on measurements of the total activity concentration in bulk samples (Bq kg −1 , Bq m −2 ); • relevant processes are deliberately excluded, for example resuspension. Resuspension (Pröhl, 2003) is the process by which contaminated particulate material in soil or sediments are remobilized into the air or the water column, respectively, by mechanisms such as wind action, sea currents and mechanical disturbances by human activity. Although generally not considered to be a significant pathway of human exposure it may become an important secondary source of airborne exposure for example due to mining activities because of the larger amounts of dust produced and the increased amounts of radionuclides that can be inhaled or following a nuclear accident; • modeling of atmospheric dispersion of the radionuclides is too simple and inadequate. For instance, when sites with complex topography are modeled with Gaussian model, which is appropriate for flat sites (Leelőssy, Lagzi, Kovács, & Mészáros, 2018). • the radioactive particle characteristics are not considered in the characterization of the source term, for example characteristics of radionuclide inventory of a specific environmental compartment (IAEA, 2011); • it is assumed that the system is in an equilibrium state when this is not the case (i.e., the dynamics of the process are neglected); Because of the too simple mathematical structure used in many cases, conceptual model uncertainty together with parameter uncertainty can be considered a major contributor to the total uncertainty budget of model output. Objectives of the present study are to apply a methodology for characterizing and quantifying conceptual model uncertainty to the modeling of interception of wet deposited atmospheric pollutants and to distinguish contribution of the conceptual model uncertainty to the total uncertainty budget of the model output from the contribution provided by parameter uncertainty.

Total Uncertainty of a Model Output
The present study builds upon the analysis and results of (Gonze & Sy, 2016), in which two processbased models for interception of wet deposited pollutants are developed and analyzed in terms of parameter uncertainty and model performance. In (Gonze & Sy, 2016) as well as in the present study, it is assumed that the measure of the total uncertainty of model output is given by the residual between modeled output and the experimentally obtained quantity. Under this assumption, it is also implicitly granted that the data used to carry out quantitative uncertainty analysis is representative of the studied quantity/process. The effect of variability also impacts the residuals and hence differentiation between uncertainty and variability needs to be addressed. In the present study, variability is assumed to be solely induced by the physicochemical characteristics of the pollutant. Not all the uncertainty contributions listed in Section 1 propagate directly to the model output. For example, scenario uncertainty is usually addressed by a "What if" type of analysis (Refsgaard, Van der Sluijs, Brown, & Van der Keur, 2006). However, parameter uncertainty and conceptual model uncertainty do propagate to the model output. The authors that produced datasets used in the present study did not provide any information about quality of measurement and/or of the sampling strategy and therefore input uncertainty is not accounted for. However, in (Gonze & Sy, 2016) a (not formal) sensitivity analysis based on the mathematical analysis of the models, estimation of the ranges of variation for input quantities and Monte Carlo calculations was carried out and impact of input variables and model parameters on model output was assessed.

Equilibrium Model and Kinetic Model for Quantifying Wet Interception Process
Wet deposition is the process for which airborne pollutants are entrained by falling rain droplets and deposit onto the ground surface. In the case of a vegetated surface, wet interception is the process by which a fraction of the deposit is intercepted by the aerial biomass (e.g., shoot). The process of wet interception is quantified in terms of wet interception factor f B , defined as: where f is the ratio between the amount of pollutant retained by plant and the total amount of pollutant deposited on the soil plant system. The f B is obtained by normalizing f to the standing biomass B [kg m −2 d.w.] in order to account for different growth phases of the plants. f depends upon the plant/foliage characteristics, precipitation characteristics and element characteristics. The amount of pollutant absorbed by the plant will be due to the balance between drainage and absorption mechanisms, water storage capacity and chemical valence of the pollutant. In fact, it is observed that inert particles interact more with leaf surfaces than divalent cations and these are more reactive than monovalent cations (Gonze & Sy, 2016). Anions (e.g. iodine) do not interact with plant surfaces. Two process-based models are presented in detail in (Gonze & Sy, 2016). These have been derived from mass balance equations with explicit parameterization of the hydrological, biological and chemical mechanisms involved.
The equations used in the analysis are the analytical solutions to the system of mass balance equations described in detail in (Gonze & Sy, 2016). The equilibrium model (EM) has been obtained by assuming that the interaction between pollutant and plant foliage is instantaneous and reversible (Gonze & Sy, 2016). The kinetic model (KM), instead, has been obtained by assuming that the interaction between pollutant and plant foliage is irreversible and depends on the kinetics of absorption and drainage processes. The EM and KM equations are provided here for completeness: • EM: • KM: and To use Equation (2), (3), (4), (5) and (6) • need to be calculated starting from the parameters CR, JK, SLA, L and input quantities B, H, I, p, E 0 listed in Table I.
In Table I, the concentration ratio CR (mm) mathematically represents the interaction between Absorption velocity pollutant and plant foliage in the EM whereas the absorption velocity J/K (mm h −1 ) mathematically represents the interaction between pollutant and plant foliage in the KM. The throughfall coefficient p provides the fraction of rainfall the reaches the ground without hitting the plants. The biomass density B is specific to each herbaceous species whereas the rainfall height H and the rainfall intensity I are specific to the rainfall event considered in the experiment. Different assumptions for the dependence between wet interception factor and chemical valence have been tested to quantify the impact of variability related to chemical valence. In (Gonze & Sy, 2016), it has been shown that both the models perform best when the chemical valence of the element interacting with plant is accounted for. For anionic particles (i.e. which undergo no chemical interaction with the plant, that is CR = J/K = 0) EM and KM reduce to the same equation. For water (H 2 O), the models simplify to a form similar to the Horton equation (Horton, 1919).

Available Data
Data used in the present study consist of 440 observations for wet interception factor f B obtained mostly with controlled experiments using monovalent cations (Cs), divalent and trivalent cations (Ba, Be, Cd, Cr, Pb, Sr), inert particles (Polystyrene with median aerodynamic diameter of 1 μm and 3 μm), anions (I), water (H 2 O) and different herbaceous species (e.g., grasses, clovers and weeds). Data considered originated from experiments that were conducted either in greenhouses or outdoor plots, in which at least the physicochemical form of the pollutants and the rainfall height were known and documented in the related publications (Chadwick & Chamberlain, 1970), (Angeletti & Levi, 1977b), (Angeletti & Levi, 1977a), (Hoffman, Thiessen, Frank, & Blaylock, 1992), (Hoffman, Thiessen, & Rael, 1995), (Kinnersley, Goddard, Minski, & Shaw, 1997), (Bengtsson, Gärdenäs, Eriksson, Vinichuk, & Rosén, 2014). In addition, three observations from Chernobyl fallout are included (Jacob et al., 1993). The set of collected observations covers a wide range of variation both for rainfall amount and intensity. As an example, in Fig. 1, the experimental wet interception factor f B is plotted against the rain intensity I (mm h -1 ). In the experiments considered, observations are obtained with either constant or intermittent precipitation. For an intermittent precipitation, an equivalent rainfall amount and duration were estimated by cumulating event-based quantities, and an equivalent intensity was estimated as the ratio of the cumulated amount and duration.
Of the 440 observations, (Gonze & Sy, 2016) use 363 data points for calibrating the free parameters CR, SLA, L for EM and J/K, SLA, L for KM as shown in Fig. 1 (top). For these datapoints all relevant input variables in the models are available. Data for model calibration consists of 27 datapoints for monovalent cations (Cs), 125 datapoints for di-valent cations (Be, Sr), 88 datapoints for inert particles (Polystyrene with median aerodynamic diameter of 1 μm and 3 μm), 108 datapoints for anions (I) and 15 datapoints for H 2 O. The data for model prediction instead consists of 1 datapoint only for monovalent cations (Cs), 24 datapoints for anions (I), 11 datapoints for inert particles (Polystyrene) and 41 datapoints for divalent and trivalent cations (Ba, Be, Cd, Cr, Pb, Sr). The 77 data points used for prediction analysis are shown in Fig. 1 (bottom). Although the experiments are conducted on different types of herbaceous species, no clear dependency of f B upon herbaceous species is observed in (Gonze & Sy, 2016).
For the prediction data, information either on the standing biomass B and rainfall height H is not available and is produced "artificially" by means of random sampling over the experimental distribution of H or B of the calibration data.
Although the present study considers the process of wet interception in context of radioecology, not only data for radionuclides are considered but also data available for inert particles (i.e., Polystyrene) and water (H 2 0) are used in order to account for the benefit of a larger database to the statistical analysis. The EM and KM, in fact, account for chemical valence of pollutants and do not depend on radioactive properties. Data used in the present study is included in Section 6. Also, it is made extensively available in (Gonze & Sy, 2016).

Bayesian Inference Analysis and Parameter Uncertainty
In (Gonze & Sy, 2016) Bayesian analysis is carried out for calibrating the model parameters (CR, SLA, L) from EM and (J/K, SLA, L) from KM and for determining their uncertainty. Bayesian inference (Box & Tiao, 2011) consists of calculating the posterior distribution functions for relevant model parameters by means of Monte Carlo Markov Chain simulations on the basis of assumed priors and available observations (stochastic realizations of the modeled process).
The input variables, which are model components for which the value is known are B, H, I. These are driving input quantities specific to each herbaceous plant (for B) or rainfall event (H, I) considered in the laboratory experiment and ideally measured for each realization. Being different among the experiments, they cannot be estimated by any statistical calibration technique.
The model parameters, for which the uncertainty is quantified using the Bayesian analysis are CR, JK, SLA, and L. In particular, SLA and L are driving parameters that moderately vary with herbaceous species and can be assumed constant among the realizations (although uncertain). For the specific foliage storage L, a uniform distribution between 0.1 and 0.5 mm was defined as prior in order to cover the typical range of values found in the literature for the maximum thickness of the water film that can be retained by leaves (Clark, 1940) (Merriam, 1961). For the specific foliage area, SLA the range from 12 to 35 m 2 kg -1 d.w., consistent with the values in (Jouven, Carrère, & Baumont, 2006), was adopted as such range of values corresponds not only the interception by the foliage but also by the other organs such as the stem. For CR and J/K, no preliminary accurate information could be used to appropriately define a prior distribution. However, to avoid implausible values the upper bounds of the uniform distributions were fixed at 10 mm for CR and 100 mm h -1 for J/K respectively whereas the lower bounds were fixed to 0 for both quantities to account for nonreactive (i.e., nonabsorbed) species such as iodine.
E 0 and p are quite insensitive quantities and therefore are fixed to a deterministic value. As climatic conditions may have been very specific for especially in-door experiments, it cannot be excluded that relatively high evaporation rates occurred. In (Gonze & Sy, 2016) two extreme situations were considered and tested, by setting E 0 value to either 0 mm h -1 (no evaporation) or 0.6 mm h -1 (high evaporation). The best agreement in (Gonze & Sy, 2016) between model and experimental observations was obtained with E 0 = 0.6 mm·h −1 . Therefore, in the present study, the evaporation rate is set to E 0 = 0.6 mm·h −1 . The throughfall coefficient p is set to 0, which is valid for a dense, nonsparse herbaceous layer because the selected experiments were not dealing with sparse vegetation.
The posterior distribution functions of the model parameters (CR, SLA, L) from EM and (J/K, SLA, L) from KM are obtained, with which the posterior distribution functions of the models' output f B is calculated.
In the Bayesian formulation, the quantity f B is of stochastic nature and is normally distributed (based on the principle of maximum entropy (Jaynes, 1957)) with mean value μ i from EM and KM output and residual variance S i 2 obtained for each data point i: In (Gonze & Sy, 2016), 15.000 simulations were used to obtain Bayesian statistics for each of the relevant parameters. The probability distribution functions (pdfs) are obtained for all the relevant parameters and for residual uncertainty S i .

Quantification of the Conceptual Model Uncertainty
Given that in the analysis parameter uncertainty and variability are accounted for and given the formalism presented in Section 2.4, it is justified to consider that the residual uncertainty S fully represents the conceptual model uncertainty. Hence, S will be referred to as the conceptual model uncertainty throughout the rest of the study. Predictive calculations for all 77 datapoints (i = 1,…,77) are carried out with and without propagating the conceptual model uncertainty S, that is, the posterior distributions of the models' outputs have been calculated from the posterior distributions of the model parameters by using and The conceptual model uncertainty S is quantified by comparing the two pdfs obtained with Equation (8) and (9), which account for propagated S, p and only propagated p, respectively. The pdfs are then compared by means of different statistical metrics in order to evaluate the relative impact of parameter uncertainty and conceptual model uncertainty on the total uncertainty.

Metrics Used to Compare pdfs and to Quantify Conceptual Model Uncertainty
The conceptual model uncertainty for the EM and KM wet deposition models is quantified by comparing: • the error bars obtained with the 95th percentile for model output with propagated S and without propagated S; • coefficients of variations (Özkaynak, Frey, Burke, & Pinder, 2009).
The error bars for each of the datapoints considered are used to deduce from 95% credibility interval of model output and 95% credibility interval of parameter uncertainty, the 95% credibility interval of conceptual model uncertainty. From the width of the error bars obtained with and without propagating S, it is possible to visualize how the conceptual model uncertainty S affects the total uncertainty compared to parameter uncertainty for each datapoint.
The coefficients of variations (CV) are calculated using the following mathematical expression: (1) quantify the difference between the pdf that accounts for the contributions of the conceptual model uncertainty and parameter uncertainty and the pdf that includes only the contribution of the propagated parameter uncertainty; (2) quantify the effect of the variability on the model output by comparing f B obtained by calibrating S over the whole dataset (obtaining a unique range of S for all pollutants, independently of their chemical valence) and f B obtained by calibrating S over the radionuclidespecific subsets S i (obtaining single ranges of S for each type of pollutant). The latter is the standard approach in the present study.

Application of Data and Use of Output from Performed Parameter Uncertainty Analysis
Predictive f B was obtained by using posterior distributions for the parameters inferred in (Gonze & Sy, 2016). Note that, 15.000 values of f B were obtained by considering the 15.000 values of the posterior distributions of each of the relevant parameters and input values. The posterior distributions of f B are shown in Fig. 2 for one datapoint. The distribution of f B were compared to those recovered in (Gonze & Sy, 2016) for each of the datapoints. The values were very similar but not always identical, since for predictive calculations the missing input values were obtained randomly in both studies using different seeds for sampling over the available distributions. Data for iodine were used to check the consistency of predictions with EM and KM obtained in (Gonze & Sy, 2016). As mentioned in Section 2.2 for iodine the EM and KM provide identical model outputs. This was confirmed in the present study, indicating also the correct implementation of the models in the programming environment R (Grunsky, 2002) Version 3.6.2.

Graphical Quantification of Conceptual Model Uncertainty Via Error Bars
Error bars for all 77 datapoints and both KM and EM were produced. An example is displayed in Fig. 3, which shows the difference between EM and KM with and without considering S for obtaining the posterior distributions of f B for one datapoint for divalent cation (barium and lead), for monovalent cation (cesium) and for anion (iodine). In Section 6 (Supporting Information) the results obtained for 77 datapoints with and without S are presented in tabulated form. The difference between error bars with and without conceptual uncertainty S is in general larger for EM than KM. It is more or less the same in the case of data points for which the two models perform similarly. The difference between error bars with S and without S shows that parameter uncertainty contributes less to the overall uncertainty of the model output than the uncertainty of the model structure. An exception is Cs-137, for which the posterior distribution of f B without S for KM (due to parameter uncertainty) is already rather large. In fact, the posterior distribution function obtained for J/K for the monovalent cation Cs-137 is larger than for other particles and makes most of the parameter uncertainty.

Coefficient of Variations (CV 95 )
The conceptual model uncertainty S increases the dispersion of the model output for the KM case on average by a factor 13 (from the ratio of the CVs obtained) and by a factor 22 for the EM case. This is shown in Fig. 4, in which CV 95 from pdfs obtained with (right) and without propagating S (left) are plotted. For Cs-137 and KM, the dispersion of the model output around the median for Cs-137 without S is rather large compared to the other elements. From  Fig. 5, it can be observed, as expected, that the posterior distributions obtained with single S have smaller dispersion around the median value compared to these obtained with one unique S. This is due to the variability arising from the chemical valence of the pollutant. The dispersion of the f B values around the median increases by a factor 4 for both KM and for EM.

DISCUSSION
A methodology based on the analysis of residuals within a Bayesian framework has been applied to quantify conceptual model uncertainty of two process-based models, EM and KM, for modeling the process of interception of wet deposited pollutants. In particular, this methodology discriminates between conceptual model uncertainty and parameter uncertainty. In this way, the impact of conceptual model uncertainty is explicitly considered and not implicitly attributed to parameter uncertainty. The applied methodology combines the use of residual analysis and the Bayesian approach, which allows for quantification of conceptual model uncertainty efficiently and includes expert judgement in a formalized way. The coefficient of variation CV 95 is the metric applied to quantify impact of conceptual model uncertainty and it includes 95th percentile in its estimation. The 95th percentile is important in radioecological assessments, especially for regulatory purposes (ICRP, 2006). This metric proves to be adequate not only to quantify how conceptual model uncertainty affects the outputs of EM and KM, respectively, but also to quantify how natural variability contributes to modeling outputs.
Nevertheless, there are fundamental assumptions that underlie the study, namely that the Bayesian residual is a suitable measure of total uncertainty of EM and KM and that the data used are representative for the scientific question under study.
The Bayesian methodology allows for straightforwardly propagating the uncertainties related to the parameters and model structure to the model output. This has the advantage that typical assumptions on additivity and correlations between the different uncertainty contributions need not to be carried out.
If no adequate measurement and sampling information is available, the uncertainty of input variables that depend on the experimental set-up (such as H, I, B in the present study) cannot be properly quantified.
Dedicated effort already during the experimental procedure to gather information about uncertainty related to measurements and sampling is necessary. However, such information is not always gathered because data are often collected for specific purposes and are not targeted to modeling activities. In several studies-from various environmental disciplines -the same issue is reported for quantities like evapotranspiration, temperature, soil water content, and so on (Radwan et al., 2004) (Højberg & Refsgaard, 2005) (Jin et al., 2010) (Alderman & Stanfill, 2017). Nevertheless, methods to quantify uncertainties related to sampling activities and/or measurements exist (Kirkup & Frenkel, 2006). In particular, assessing the performance of measurement sensors (Chen, Chen, & Chen, 2018) (Damgaard, 2020) is crucial to determine the uncertainty due to measurements. In studies where input uncertainty is considered and discussed at least qualitatively, more exhaustiveness in the determination of the total uncertainty budget of a model is obtained (see Section 3.2 in Gondwe, Merediz-Alonso, & Bauer-Gottwein, 2011;Radwan et al., 2004).
In environmental models, data available not only should provide as much information as possible re-garding input quantities but also need to be representative of the situation being analyzed that is cover the entire validity domain of the model. Instead, data may be rather limited: for the present radioecological study the 440 observations available with large range of variation of quantities such as I or H and for different types of radionuclides are more the exception rather than the rule.
For Bayesian analysis the choice of priors requires information on the possible range of variation of the quantities involved and the choice of most sensitive parameters needs to be outweighed based also on the possibility to build a probability distribution function. A more formal sensitivity analysis on all input quantities (Saltelli, 2002) is an interesting perspective when the aim to characterize the total uncertainty budget.
In the field of radioecology not so many different conceptual models for a given transfer process may be considered or are available (Urso et al., 2015) as can be the case instead in other disciplines (e.g., Gondwe et al., 2011;Højberg & Refsgaard, 2005). Indeed, the availability of two different processbased models for describing the process of interception of wet deposited pollutants on herbaceous species is rather unique and motivated the present study. The results from various studies (Gondwe et al., 2011;Jin et al., 2010;Refsgaard et al., 2006) indicate that conceptual model uncertainty contributes even more than parameter uncertainty to the total uncertainty budget of a model in agreement with the present findings. Also, the analysis of conceptual model uncertainty improves model understanding and confidence in model predictability. It is especially when data are not available and hence models are needed for prognostic assessment that improved understanding of the uncertainty contributions is necessary. The present case is very specific and more modeling examples and research is required before generalizations valid for radioecological models can be made. However, the developed methodology can be certainly applied to other radioecological models.

CONCLUSION
The main finding is that KM describes experimental data better than EM in agreement with (Gonze & Sy, 2016). However, in (Gonze & Sy, 2016) the predictive model performance is carried out by means of posterior predictive loss criterion (PPLC). The present study suggests model performance testing by looking at how the model structure affects the gap between the observations and the corresponding predictions. In fact, the substantial differences in the EM and KM mathematical structures and the fact that the same set of data was used for both models justifies the consideration of conceptual model uncertainty as being the dominant contribution to residuals. This analysis shows that for the majority of data points the process-based approach, which models the interception process as a kinetic process, provides reduced structural uncertainty (and hence better agreement) than the modeling of the interception process under equilibrium conditions.
The main conclusions are: • the 95% credibility intervals are larger for EM than KM for the majority of the datapoints. • CV 95 values obtained with and without single S show that, on average, conceptual model uncertainty is larger for EM than KM by at least a factor 10. • CV 95 values obtained using single S and unique S show that the effect of variability related to chemical valence on EM and KM output affects the dispersion of the model outputs around the mean value by a factor 4.

ACKNOWLEDGMENTS
This work has been supported by the TERRITO-RIES project, part of the CONCERT European Joint Programme funded from the Euratom research and training programme 2014-2018 under grant agreement No. 662287. The present publication reflects only the authors' view. Responsibility for the information and views expressed therein lies entirely with the authors. The European Commission is not responsible for any use that may be made of the information it contains. Open Access funding enabled and organized by Projekt DEAL.