Using a Bayesian Network Model to Predict Risk of Pesticides on Aquatic Community Endpoints in a Rice Field—A Southern European Case Study

Bayesian network (BN) models are increasingly used as tools to support probabilistic environmental risk assessments (ERAs), because they can better account for uncertainty compared with the simpler approaches commonly used in traditional ERA. We used BNs as metamodels to link various sources of information in a probabilistic framework, to predict the risk of pesticides to aquatic communities under given scenarios. The research focused on rice fields surrounding the Albufera Natural Park (Valencia, Spain), and considered three selected pesticides: acetamiprid (an insecticide), 2‐methyl‐4‐chlorophenoxyacetic acid (MCPA; a herbicide), and azoxystrobin (a fungicide). The developed BN linked the inputs and outputs of two pesticide models: a process‐based exposure model (Rice Water Quality [RICEWQ]), and a probabilistic effects model (Predicts the Ecological Risk of Pesticides [PERPEST]) using case‐based reasoning with data from microcosm and mesocosm experiments. The model characterized risk at three levels in a hierarchy: biological endpoints (e.g., molluscs, zooplankton, insects, etc.), endpoint groups (plants, invertebrates, vertebrates, and community processes), and community. The pesticide risk to a biological endpoint was characterized as the probability of an effect for a given pesticide concentration interval. The risk to an endpoint group was calculated as the joint probability of effect on any of the endpoints in the group. Likewise, community‐level risk was calculated as the joint probability of any of the endpoint groups being affected. This approach enabled comparison of risk to endpoint groups across different pesticide types. For example, in a scenario for the year 2050, the predicted risk of the insecticide to the community (40% probability of effect) was dominated by the risk to invertebrates (36% risk). In contrast, herbicide‐related risk to the community (63%) resulted from risk to both plants (35%) and invertebrates (38%); the latter might represent (in the present study) indirect effects of toxicity through the food chain. This novel approach combines the quantification of spatial variability of exposure with probabilistic risk prediction for different components of aquatic ecosystems. Environ Toxicol Chem 2024;43:182–196. © 2023 The Authors. Environmental Toxicology and Chemistry published by Wiley Periodicals LLC on behalf of SETAC.


INTRODUCTION
Today's environmental risk assessment (ERA) of chemicals is mainly based on deterministic approaches, which usually rely on a single-value risk characterization.The risk scores used for lower tier assessments are, for example, risk quotients (RQs) or environmental quality standards, which are often based on protective exposure and effect characterizations such as a worst-case estimation of exposure and the most sensitive effect concentration (toxicity test).However, pesticide fate simulation models (Pereira et al., 2017) already allow the calculation of exposure distributions that encompass spatial and temporal variability (Schmolke et al., 2010).Such process-based pesticide fate models can aid in obtaining realistic exposure characterization when monitoring data are scarce.They can also assist in the projection and analysis of future land-use and climate change impacts on pesticide exposure concentrations by integrating a wide diversity of scenario combinations such as agricultural practices, soil properties, crop types, and meteorological conditions in the prospective exposure assessment.Moreover, exposure models are relatively quick and costefficient tools for assessing the exposure of pesticides to the environment, compared with extensive environmental monitoring campaigns (Lammoglia et al., 2018).
In recent years, the influence of climate change on pesticide fate and transport has been the subject of increased concern (Bloomfield et al., 2006;Delpla et al., 2009;Lamon et al., 2009;Noyes et al., 2009).Changes in climate conditions (temperature and precipitation) and land-use practices can lead to shifts in ecosystem structure and function, as well as their hydrological processes.In turn, this may lead to changing responses of the affected species to contaminants (Landis et al., 2013).A realistic assessment of risks posed by these expected stressors is crucial for the future ecological sustainability of the Albufera National Park in Valencia, Spain (Figure 1), which is the site of the present case study.
Traditional effect assessment based on standard toxicity tests, which are commonly used in lower tier regulatory risk assessments, tend to ignore many aspects of ecological realism: indirect effects, the complexity of the population and population dynamics, and complex interactions occurring between populations in a community structure.Instead, to account for uncertainties related to effect assessment, assessment factors are commonly applied to the most sensitive species of toxicity tests or a given hazardous concentration (typically the hazardous concentration to 5% of a species [HC5]) from a species sensitivity distribution (SSD; Schmolke et al., 2010;Topping et al., 2020).This commonly used lower tier effect assessment neglects species interaction and lacks insights into the propagation of effects on single species to different trophic levels of the ecosystem ( Van den Brink et al., 2006).Predicted no-effect concentrations (PNECs) derived from micro-and mesocosm experiments are generally used in the prospective risk assessment of single compounds, and barely perform spatially explicit or scenario-based risk assessments, although they can offer such opportunities.The Predicts the Ecological Risk of Pesticides (PERPEST) model (Van den Brink et al., 2002) was designed to perform more realistic effect assessments by using case-based reasoning to search for analogous pesticide contamination cases based on type of experimental ecosystem, the exposure pattern, and the mode of action of the evaluated substance (Larras et al., 2022;Van den Brink et al., 2002, 2006).The model predicts pesticide effects on a series of biological structural and functional endpoints, which have been derived for species tested within micro-and mesocosm setups (Davis et al., 2013).The PERPEST model uses data created by micro-and mesocosm experiments to conduct more sophisticated risk assessments, by trying to predict specific population and community endpoints that may be impacted by pesticide contamination (Larras et al., 2022;Van den Brink et al., 2002, 2006).
In Europe, prospective ecological risk characterization is often calculated by comparing predicted exposure concentration (PEC) with PNEC following a so-called deterministic approach, that is, generating a single-value risk score (Di Guardo & Hermens, 2013;Schmolke et al., 2010).These deterministic approaches involve the use of uncertainty factors (e.g., assessment factors), or the exceedance and frequency of exceedance of safe thresholds to derive qualitative output that lacks indications of the level of certainty related to the input and output parameters (European Commission, 2006).In reality, however, pesticide exposure and effects have spatial and temporal variability.This is mainly influenced by environmental and biological characteristics, and by pesticide application patterns (European Commission, 2007).Improving prospective ERAs and considering and integrating future scenarios into the risk assessment of pesticides will aid in the prevention of further and future damage (Topping et al., 2020).Probabilistic approaches have been proposed as more suitable methods to account for variability and uncertainty of exposure, effect, and risk (Carriger & Newman, 2012;European Commission, 2006;Solomon et al., 2000;Verdonck, 2003).Previously proposed probabilistic methods are the joint probability curves (Giddings et al., 2000;Verdonck, 2003), the quantitative overlap (Hall et al., 2000;Solomon et al., 2000), and the RQ distribution approach (Campbell et al., 2000;Mentzel et al., 2021;Verdonck, 2003).Limitations of the existing probabilistic tools are that their outputs can be hard to interpret and/or that it can be difficult to communicate the risk to decision-makers (Dreier et al., 2020;Giddings et al., 2000).
Bayesian networks (BNs) can overcome some of these limitations and better communicate the quantified uncertainties to decision-makers and other stakeholders (Carriger et al., 2016;Carriger & Newman, 2012).Bayesian networks can be defined as directed acyclic graphs or probabilistic graphical models that contain nodes (variables) linked through arcs representing conditional probability tables (Aguilera et al., 2011;Kaikkonen et al., 2021).The nodes have assigned states (intervals or categories), which are quantified by probability distributions.Based on new evidence (e.g., scenarios), BNs use Bayes' rule to update the probability distributions throughout the network (Carriger et al., 2016;Kanes et al., 2017).While being used in situations where data are limited, BNs are able to incorporate various sources of information, for example, expert elicitation, other model outputs, or other information literature (Carriger et al., 2016;Carriger & Newman, 2012;Gibert et al., 2018;Hamilton & Pollino, 2012).Moreover, BNs can act as a metamodel (see Martínez-Megías et al., 2023;Mentzel et al., 2022) allowing the incorporation of input and output data and assumptions from different models, for example, process and case-based prediction models into a single network model.This approach can allow the display of uncertainty of all model compartments in a transparent way.Despite their extended use in disciplines such as human medicine or social sciences, their application for the ecological risk assessment of chemicals remains limited, although it is increasing (Kaikkonen et al., 2021).
Our overall purpose was to develop and evaluate a probabilistic (Bayesian) network model for calculating the risk of pesticides to aquatic communities, by combining information from different scenarios, a pesticide exposure model and a pesticide effect model.More specifically, the goals of the model and case study we present were to explore the following issues: (1) How will the risk of pesticide to individual biological endpoints vary across the environmental scenarios?
(2) How will the risk to biological endpoint groups vary across the different pesticide types?and (3) How will the variation in risk across scenarios and pesticide types manifest at the community level?

DATA AND METHODS
Our general approach was to integrate predicted outputs from exposure, quantified as probability distributions, and effect models quantified as probabilities, into a BN serving as a metamodel (Figure 2).We developed a BN metamodel structure incorporating spatial variability in the risk estimation of pesticides for various endpoints in the aquatic ecosystem, as explained in the next sections.
The first subsection describes the study area including the scenarios of climate and land-use changes and some of the related challenges.The following two subsections describe the derivation of inputs for the BN model according to the conceptual model (Figure 2): predicted exposure data for three selected pesticides, and predicted concentration-effect relationships for biological endpoints.The last subsection gives more information on the structure and parameterization of the BN model, including the approach for integrating risks from the individual endpoint level to endpoint groups and to the community level.
The model was run for six scenarios (Table 1) representing different management practices ("application") and FIGURE 2: Conceptual model for the effect estimation of a pesticide (acetamiprid) on an aquatic community.The pesticide exposure distribution derives input from the RICEWQ model and is determined by the associated future climate and application scenarios.The PERPEST model input derives the probability of effect on biological endpoints, which is aggregated to endpoint groups and in turn to the community.meteorological and hydrological conditions ("Year"), for each of the selected pesticides.

The case study region
The case study region was located in the Albufera Natural Park (a lake enclosed by rice fields) known for its diversity of bird and fish species (Soria, 2006; Figure 1).The study area is a coastal wetland approximately 5 km south of Valencia on the Mediterranean Spanish coast, with an area of approximately 210 km 2 (Figure 1).The Natural Park has ecological relevance because it is a nesting and transfer point for approximately 250 species of migratory birds; it has been designated a special protection area by the Birds Directive of the European Commission (2010), and is listed as a core breeding and resting site in Natura 2000 (European Commission, 2008) as well as the Ramsar Convention of wetlands (Calvo et al., 2021;Generalitat Valenciana, 2020;Vera-Herrera et al., 2021).Within its boundaries, 34% of Spanish rice is produced (Canet et al., 2003); 73% of the wetland is dedicated to rice cultivation (Vera-Herrera et al., 2021), and the remaining area is almost all an oligohaline lake (i.e., Albufera Lake).The rice paddies serve as a habitat for several species of aquatic plants, invertebrates, and fish, and provide food to the local and migratory birds.Rice production relies on the use of herbicides, insecticides, and fungicides, all of which pose potential hazards for the aquatic ecosystem, as residues of the pesticides applied in the rice paddies are transported downstream through the drainage channels and have been detected in Albufera Lake (Calvo et al., 2021;Vera-Herrera et al., 2021).Rice production and other anthropogenic stressors have been harming the lake's water quality and ecosystem over the last century (Calvo et al., 2021;Vera-Herrera et al., 2021).Moreover, future land use and climate change are expected to alter the distribution and fate of pesticides in aquatic environments (Bloomfield et al., 2006;Noyes et al., 2009).In the Mediterranean, it is expected that droughts will occur more frequently, and water will be less available, thereby resulting in lower dilution of pesticides.On the other hand, severe precipitation events are expected to occur more often, which may result in higher pesticide runoff.For this southern region, an expected temperature increase may facilitate microbial degradation of pesticides (Arenas-Sánchez et al., 2016;Balbus et al., 2013;García de Jalón et al., 2014;Noyes et al., 2009).Nevertheless, in general, pesticide impacts are expected to increase with higher temperatures (Op de Beeck et al., 2017).These authors have formulated three reasons for this increase: (1) increasing application rates at higher dosages related to a higher pest abundance; (2) increasing temperatures causing higher toxicity of the same concentrations; and (3) increasing pesticide runoff from fields related to increased precipitation (Op de Beeck et al., 2017).

Exposure prediction with RICEWQ-Prediction and settings
We have focused on pesticide applied to rice fields, and therefore the RICE Water Quality model (RICEWQ) was used to simulate pesticide exposure in the water of rice paddies (Karpouzas & Capri, 2006;Miao et al., 2004).The exposure data was obtained from Martínez-Megías et al. (2023), where more details about the exposure prediction can be found.The RICEWQ model is a process-based model that, at the field level, simulates pesticide runoff specific for use in rice paddies (Martínez-Megías et al., 2023;Williams et al., 1999).Thus far, it is considered to be the most suitable and reliable model for higher tier pesticide fate and exposure prediction (Daam et al., 2013;European Commission, 2003;Karpouzas & Capri, 2006;Pereira et al., 2017).In addition, it has been widely applied Also given, for each mean and SD, the selected pesticides and pesticide type, as well as details of the RICEWQ input scenarios by year (determining the climate conditions) and the application scenario (three significant digits were chosen).MCPA = 2-methyl-4-chlorophenoxyacetic acid; RICEWQ = RICE Water Quality model.
Probabilistic risk of pesticides to aquatic communities-Environmental Toxicology and Chemistry, 2024;43:182-196 (Karpouzas & Capri, 2006;Miao et al., 2004) to track the fate of both parent compounds and their metabolites (Christen et al., 2006).Various pesticide fate processes are included in RICEWQ modeling, such as biological, hydrological, and physicochemical processes (Wang et al., 2019).The model requires the following inputs: daily weather information, paddy soil properties, pesticide chemical properties, pesticide management information, and water management practices (Wang et al., 2019).The reader is referred to Williams et al. (1999) for more detailed information on the model function, assumptions, and description.Validation of the RICEWQ model for this specific study area has not yet been completed (Martínez-Megías et al., 2023), but the model has been validated in similar rice production scenarios and has provided satisfactory results (Karpouzas and Capri, 2006;Pereira et al., 2017).In the present study, we used predictions obtained from the latest version, RICEWQ 1.92 (Waterborne Environmental, 2022).We selected three active substances (Table 2) that are regularly applied by farmers in the rice paddies of the Albufera Natural Park.The insecticide acetamiprid, the fungicide azoxystrobin, and the herbicide MCPA were selected to show the applicability and integration of different pesticides types into the approach developed previously by Calvo et al. (2021) and Rodrigo et al. (2022).Acetamiprid and azoxystrobin are usually applied at the beginning of July, whereas azoxystrobin is applied at the end of July and the beginning of August (Martínez-Megías et al., 2023).
The derived pesticide application scenarios are based on the dosages recommended by the manufacturers, from which we derived two scenarios: one maximum recommended dosage (referred to as "baseline" application throughout the study) and one that is 150% of that baseline dosage ("baseline+50%").The baseline application is related to current practice according to the product description of the plant protection product applied.The baseline+50% can be considered a worst-case scenario for increased dosage, as indicated by farmers due to increasing pest resistance and occurrence.
Initially, we had aimed to have climate projections for at least three different emission scenarios, to include more variability and to use them as input for the RICEWQ model runs based on what had been previously used in a study by Pool et al. (2021).However, climate projection data for the nearby meteorological station ( 8416) was limited to one emission scenario: representative concentration pathway 8.5.Therefore, only one climate projection data set was collected from the Spanish State Meterological Agency (2021) at "Climate projections for the XXI Century-Daily data," derived with the model GCM MPI-ESM-LR.Based on this data set, three "climate conditions" scenarios for 2008, 2050, and 2100 were used to run the exposure prediction model (Table 1).
The exposure prediction model was run for 552 spatial units of rice crop clusters (Martínez-Megías et al., 2023), which constituted the spatial variation underlying the exposure distribution in the present study.A detailed description of the assumptions made to derive these clusters, and the automatization of the RICEWQ with a handy interface, are available in Martínez-Megías et al. (2023).The autoRICEWQ (open source under GPL-3.0License, programmed in Python 3) was run for the different pesticides, crop clusters, and scenario combinations and can be accessed from Fuentes-Edfuf and Martínez-Megías (2022).The peak exposure concentration from each cluster was assumed to follow a normal distribution, which was fitted to quantify the spatial variability of exposure within each scenario by mean and standard deviation (Table 1).

Effect prediction with PERPEST: Model assumptions and processing of predictions
The PERPEST model was developed to predict pesticide effects on various biological endpoints from aquatic community studies, such as phytoplankton, microcrustacean, or fish (Van den Brink et al., 2002).It can be used for risk assessment of the single or mixtures of pesticides (Rämö et al., 2018), but was only used for single pesticides separately in our study.It is considered more comprehensive than the traditional ERA that uses point estimate approaches such as risk or hazard quotient (Polidoro & Morra, 2016;Rämö et al., 2018).As mentioned, this effect prediction model applies a case-based reasoning approach based on empirical ecotoxicological data from microand mesocosm experiments (Davis et al., 2013;Rämö et al., 2018).It can represent indirect effects of contaminants, for example, via the food chain, while incorporating hydrological properties and acute and chronic exposure in the prediction ( Van den Brink et al., 2002).The PERPEST model compares environmental exposure concentrations with previous observations in mesocosm and microcosm toxicity experiments to estimate the probability of the pesticide having a toxic effect on various pesticide-type-dependent biological endpoints and endpoint groups.The PERPEST model was last Above max range-the value was higher than the maximum amount of days that could be entered, in which case the max of 1000 days was used.DT50 = degradation half-life; HC50 = hazard concentration for 50% of the species; Koc = soil adsorption coefficient; MCPA = 2-methyl-4-chlorophenoxyacetic acid.
updated in 2009 and has a limited number of studies in the database on which the predictions are based.The mode of action is a pesticide property that can be considered for some pesticides, but the selected pesticides used in our study do not have modes of action currently available for the model predictions.Therefore, the PERPEST model applies assumptions across the same pesticide class for the prediction of the effect on the different biological endpoints.
For each pesticide and biological endpoint, PERPEST predicts the probability of three effect classes along a pesticide concentration gradient.The selection of biological endpoints depends on the modeled pesticide type.Following Van den Brink et al. (2002), these three classes are defined as (1) No effect: no consistent adverse effects are observed as a result of the treatment; observed differences between treated test systems and controls do not show a clear causality; (2) Slight effect: confined responses of sensitive endpoints (e.g., partial reduction in abundance) are observed on individual sampling dates only and/or are of a concise duration directly after treatment; and (3) Clear effect: severe reductions of sensitive taxa over a sequence of sampling dates are demonstrated (Davis et al., 2013).
The reader is referred to Van den Brink et al. (2002) for more detailed information on the model function, assumptions, and description.We used the PERPEST model to predict the effect of a fungicide, herbicide, and insecticide on the biological endpoints associated with being affected by the different types of pesticides (Table 2).The pesticides we selected were not currently available in the PERPEST database.Therefore their physicochemical properties were collected from the literature and databases such as the Pesticide Properties Database (PPDB; Lewis et al., 2016), PubChem (Kim et al., 2020), and CompTox (Williams et al., 2017).The HC50 was calculated from an SSD for each of the pesticides using MOSAIC (Charles et al., 2017) with the effect concentration for 50% of the species (EC50) toxicity data collected from the ECOTOXicology Knowledgebase (Olker et al., 2022).The input information used in the toxicity data set by the PERPEST model is shown in Table 2.The PERPEST model uses a single HC50 derived from an SSD for the pesticide to estimate the effect on the different taxa.Because the model uses mesocosm studies to estimate the effects, no taxa-specific HC values were used.The HC50 is used by the model as an overall measure of toxic pressure for all species in the ecosystem.It is used to calculate hazard units for the different concentrations, and with those, make a search in the mesocosm database.For this reason, it is not taxon specific, because it is not used to make a quantitative risk assessment for the taxonomic group of concern, but rather to calculate overall toxic pressure and search for similar cases in the database.
The latest version of the PERPEST software (Van den Brink et al., 2002; Ver.4.0.0;www.perpest.wur.nl) was used to predict the probability of effects for the three selected pesticides.In general, the model was used with default settings, additional weight was put on "toxic unit" and "degradation half-life," and exposure was set to "not used." The output from PERPEST is a set of tables with probability distributions across the three effect classes along a gradient of pesticide concentration values, for each endpoint.These table were used to parameterize the conditional probability tables from the effect concentration node to each of the endpoint nodes, which are represented by arrows in the conceptual model (Figure 2).The intervals of this concentration gradient have exponentially increasing (doubling) width in PERPEST.All of the PERPEST output tables used to quantify the probabilistic links to the biological endpoint nodes can be found in the Supporting Information, Data III.An example of the PERPEST gradient output is illustrated in the Supporting Information,

BN conceptual model and assumptions
The developed BN metamodel is composed of three modules (example of concept for the herbicide MCPA in Figure 2): the scenarios and exposure (blue), the effect on the biological endpoint (green), and the predicted effect on the aquatic communities (gray).The first module, "scenario and exposure," is composed of the scenario combination (red node) that defines the exposure concentration distributions (blue node), which were obtained by fitting them to the RICEWQ model outputs (see Supporting Information, Data II).The second module is derived from the PERPEST model outputs that provide the effect concentration states and the associated probabilities of: no, slight, and clear effects of the biological endpoint nodes (green nodes; see Supporting Information, Data III).In the third module, "predicted effect on communities," each of the biological endpoint nodes are converted to Boolean nodes (true/false; dark gray) before being aggregated to their respective endpoint group nodes (Figure 2, light gray): plants, invertebrates, vertebrates, and ecosystem processes (Table 3).
For a simplified risk characterization at the community level, the probability of slight and clear effect was merged into the Boolean variable "effect" (true/false).The probability of effect on individual biological endpoints was first aggregated to endpoint groups, and further to the community level (Figure 2).The risk to an endpoint group was defined as the joint probability (P) of effect to any of the biological endpoints (n) in this group: × [ ] × ( − ( )) (1) This equation represents an "OR" expression and calculates the cumulative probability that any of the endpoints will be affected.For example, the node "effect on insect" is meant to quantify the probability of a pesticide effect to insects (true/ false); the node "effect on invertebrates" will quantify the cumulative probability of effects to any of the invertebrate endpoints (i.e., macrocrustacea and insects, zooplankton, and molluscs); and the node "effect on communities" will quantify the cumulative probability of effects to any of the endpoint groups.(This expression is equivalent to the equation used for the response multiplication for mixture toxicity of substances with an independent mode of action (see Boedeker et al., 1993); this is also sometimes referred to as "response addition" (see Gregorio et al., 2013).
In the following, an example for acetamiprid shows the aggregation from the PERPEST-defined categorized states to Boolean nodes of the biological endpoint.This example is based on the biological endpoint of insects, macrocrustacea, macroinvertebrate, microcrustacean, and rotifers (Figure 3A).The Boolean biological endpoint is summarized, and the probability of an effect on these endpoints (Figure 3B).
An assumption can be made for the effect on biological endpoints and the endpoint groups (Figure 3 and Supporting Information, Figure S.3).When some of the biological endpoints for the insecticide acetamiprid (with a scenario 4) are compared, it can be observed that nonarthropod macroinvertebrate (labeled "other macroinvertebrates" in PERPEST) was predicted to be affected by acetamiprid with a probability of 0%.Unlike the insects, macro-and microcrustaceans had a higher probability (25%-30%) of being affected.
Furthermore, the effect on the endpoint group can also be aggregated with the developed BN.The joint probability of the Boolean nodes derived the endpoint group for example invertebrates (Figure 3C).In this example, for a scenario 4, it can be concluded that the predicted probability of an effect on any of these biological endpoints of the invertebrate's endpoint group was false, with a probability of 55.1% (Figure 3C and Supporting Information, Figure S.4).
In our study, the BN was constructed with Netica software (Norsys Software).For each of the selected pesticides, the parameterized model was run by selecting a set of scenarios, for example, climate conditions at a specific time and application scenario, as evidence.The BN then calculated the predicted exposure concentration and subsequently the The biological endpoints were determined by the Predicts the Ecological Risks of Pesticides (PERPEST) model depending on the pesticide type.To enable comparison of predicted effect on the endpoint groups and community between the different selected pesticides, the biological endpoints were grouped into endpoint groups as shown in this table.probability of effect on the various biological endpoints and endpoint groups and the communities.The probability distributions were updated throughout the BN to predict the probability of effect classes on the output nodes.Model assumptions and a more detailed node description are provided in Table 4.

Interpretation of BN predictions
A parameterized example of the BN for the insecticide acetamiprid is shown in Figure 4.It displays the outcome of two scenarios with climate conditions for the year 2050: the baseline application (scenario 3, Figure 4A) and the baseline+50% application (scenario 4, Figure 4B).Considering scenario 4, the resulting exposure concentration shows the highest probability (17.4%) of values in the range 0.317 to 0.375 mg/L.In the effect concentration node, where the interval boundaries are given by the PERPEST model, the exposure concentration range is covered by the interval 0.29 to 0.58 mg/L (63.7% probability).The highest probability of clear effect, referred to as "clear risk," was predicted for the biological endpoint of macrocrustacea (12.3%).This is followed by microcrustacea with 9.19% clear risk and 9.96% slight risk.For this scenario, there is no risk to fish or to other (nonarthropod) macroinvertebrates.Furthermore, considering the endpoint groups, invertebrates show the highest risk (44.9% probability of effect).This represents the cumulative probability that any of the five invertebrate endpoints have shown a slight or clear effect.Finally, the aggregation node, summarizing the effects on the community, predicts a probability of 49.3% of the pesticide affecting at least one of the biological endpoint groups.Other examples of the parameterized BN models for the fungicide and insecticide can be found in the Supporting Information, Data I and Figures S.2, S.5, and S.8.
In the following (Figures 5 and 6), the risks predicted by the BN are summarized in a bar chart to facilitate comparison between the different scenarios, biological endpoints, and pesticide types.

Risk to biological endpoints: Comparison across scenarios
Focusing on the individual biological endpoints, the BN predicted the effect of the insecticide acetamiprid for eight biological endpoints (Supporting Information,Figure S.3.).Insects and macro-, and microcrustaceans had up to 20% risk (probability of slight or clear effect).Community metabolism, algae, and rotifers were mostly unaffected, with less than a 7% risk.Fish and nonarthropod macroinvertebrates were most likely not to be affected by the insecticide (Supporting Information, Figure S.3.).For the insecticide, the risk was similar for the scenarios representing years 2008 and 2100 (Supporting Information, Figures S1, S2, S5, and S6).For the year 2050 (scenarios 3 and 4), the risk was similar (13%-24%) for the biological endpoints of community metabolism, insects, macrocrustacea, microcrustacea, and rotifers.Scenario 4 (baseline+50% in 2050) caused the highest risk, of up to 16% to 21% for macrocrustaceans and insects, compared with the other scenarios.
For the fungicide azoxystrobin, 11 biological endpoints were considered by the PERPEST model (Supporting Information,Figure S.6).Nonarthropod macroinvertebrates and microcrustacea were predicted to have a risk of more than 50%, followed by other zooplankton taxa, phytoplankton, community metabolism, macrocrustacea, and insects.Fish and macrophytes were similarly affected, with a risk of 15% to 20%.Decomposition and periphytic algae were predicted to not be affected (4%-10% risk).The risk in general decreased slightly across the different climatic periods, being lowest in 2100, and increased with the application scenario baseline+50% (scenarios 2, 4, and 6) compared with the baseline scenarios (scenarios 1, 3, and 5).
The PERPEST model considered eight biological endpoints for the herbicide MCPA (Supporting Information, Figure S.9).Zooplankton had a probability of more than 50% to be in the clear effect state (i.e., clear risk) for a baseline+50 scenario (Supporting Information, Figures S.2 and S.4), and of approximately 30% for scenario 1 (baseline in 2008).Phytoplankton and the periphyton had a probability of 25%.Macrophytes and community metabolism were mostly unaffected, and fish and molluscs were predicted to have no risk.The predicted clear risk for any endpoint was highest for climate conditions in 2008, with slightly decreasing probabilities toward climate conditions in 2100.One exception was phytoplankton, for which scenarios 2 and 4 (baseline+50% in 2008 and 2050) gave a similar risk.There were few biological endpoints that all pesticides had in common; one of them was macrocrustacea (see the Supporting Information, Figures S.3, S.6, and S.9).For this endpoint, it was thus possible to compare the risks of the different pesticides (Figure 5).The fungicide was predicted to have the highest risk (22%-28%; Figure 5A), in accordance with general information about azoxystrobin having moderate toxicity to invertebrates (as well as other endpoints in general; Table 5).The herbicide MPCA was predicted to have a higher risk for macrocrustacea (16%-33%; Figure 5B) than the insecticide acetamiprid (9%-16%; Figure 5C).This outcome contrasts with the general information on MCPA having low toxicity and acetamiprid having moderate toxicity to invertebrates (Table 5).This apparent inconsistency suggests that indirect mechanisms have been at play in the mesocosm studies underlying the PERPEST model, resulting in a negative relationship between the herbicide exposure and macrocrustacean endpoint values.This phenomenon will be explored further in the next section at the level of endpoint groups, which allows the inclusion of all endpoints in the comparison across pesticides (Figure 6), as well as across scenarios (Supporting Information, Figures S.4, S.7, and S.10).

Risk to endpoint groups: Comparison across pesticide types
The combined risk or probability of effect to an endpoint group represents the joint probability of any biological endpoint in the group being affected (Figure 3).When the three pesticide types in the baseline scenario were compared (Figure 6, upper panel), the insecticide acetamiprid was found to have the lowest risk; the probability of effect ranged from 36% (invertebrates) to 0% to 6% (plants, vertebrates, and community processes).This insecticide has low-moderate toxicity to plants, and moderate toxicity to invertebrates (Table 5); the predicted risk to invertebrates can therefore in principle be explained by direct toxicity effects.
The herbicide MCPA posed risks to endpoint groups in the same range as the insecticide: 35% for plants, 38% for invertebrates, and below 10% for vertebrates and ecosystem  processes.Because MCPA in general has low toxicity to invertebrates (Table 5), the predicted high risk to invertebrates in this case might be interpreted as indirect effects occurring in the mesocosm studies on which the PERPEST model is based.A likely explanation could be that the herbicide has reduced the amount of plant food available to herbivore invertebrates, and thereby their growth, survival, and/or reproduction.Considering the individual endpoints in this case, for example, the 22% risk to phytoplankton is likely to affect the food source of the zooplankton and contribute to the predicted 25% risk to zooplankton.
The fungicide azoxystrobin posed the highest risk across all endpoint groups, and consequently also to the community level.The risk of the fungicide ranged from 18% (vertebrates) to 98% (invertebrates), with a combined 99% to 100% risk at the community level.Azoxystrobin has moderate toxicity to plants and aquatic invertebrates, and moderate-high toxicity to fish (Table 5).The predicted risks to the different groups can therefore be interpreted as a direct toxicity effect, although there could also be contributions from indirect food-web effects.
Comparison of risk to the different endpoint groups between the two pesticide application scenarios (Figure 6, upper vs. lower panel) can highlight the role of the different pesticide types, and potentially the importance of the indirect versus direct effects on the different components of the community.
In the baseline scenario (year 2050), the insecticide and the herbicide posed a similar level of risk to the invertebrate endpoint group (36% and 37% probability of effect, respectively).Under the baseline+50% scenario, the invertebrate risk due to the insecticide was only slightly increased by 9% points to 45%.For the herbicide, in contrast, the risk to invertebrates increased by 25% points to 63%.Considering the individual endpoints (Supporting Information, Figures S.3 and S.9), the baseline+50% application of the herbicide caused a strong increase in risk to two endpoints: zooplankton (+21% points) and macrocrustacea (+14% points).In comparison, the baseline +50% application of the insecticide only resulted in a small increase in risk (<5% points) for each of four endpoints (macrocrustacea, microcrustacea, rotifers, and insects).Hence, for invertebrates under these circumstances, the toxicity-related risk from insecticides would be outweighed by the more indirect food chain-related risk from herbicides.

Risk to the aquatic community
The community-level risk, whereby the risk to all endpoint groups is accumulated, reflects the combination of both direct and indirect pesticide risk from the different trophic levels (Figure 6; see also Supporting Information,Figures S.4,S.7,and S.10).This is most evident in the worst-case scenario (baseline+50%; Figure 6, lower panel).For the insecticide, the community-level risk (49%) was only slightly higher than the risk to invertebrates (45%).For the herbicide, in contrast, the community-level risk (86%) was considerably higher than the risk to either plants or invertebrates (57%-63%).In the case of the fungicide, for which the risk to the invertebrates group alone was 97%, the additional risk to other endpoint groups could not contribute much further to the overall communitylevel risk (99%).

DISCUSSION
We conducted a risk characterization for three selected pesticides using a novel probabilistic approach that links the inputs and outputs of exposure and effects models into a BN.The BN model developed for this purpose can predict the pesticide risk as the probability of effect on multiple biological endpoints and can aggregate the outcome as cumulative risk to endpoint groups and community.Furthermore, the BN predictions enabled comparison across different pesticide types, levels of biological organization, and scenarios of pesticide application and climate conditions.This modeling exercise has given us new insights into the dynamics of pesticide risk to aquatic communities through both direct and indirect (trophic) mechanisms, and also how these can unfold under different environmental scenarios.Such insights could only be obtained through the combination of process-based exposure modeling, case-based effect modeling, and probabilistic network modeling for integrating all of the information.
Our study has built on the work of Martínez-Megías et al. ( 2023), who used the RICEWQ model for pesticide exposure prediction in the same case study, but who used a probabilistic RQ approach for the risk characterization.The RQ is commonly used to measure the degree to which a PEC exceeds a regulatory threshold such as a PNEC, at least for the purpose of initial screening or prioritization of risk.In comparison with Martínez-Megías et al. (2023), and other studies applying a probabilistic RQ under various scenarios (see Mentzel et al., 2022), our research has further advanced the risk characterization by replacing the RQ with empirical cause-effect relationships between chemical exposure and biological effects for a range of endpoints.
Although probabilistic approaches for risk assessment such as BN have gained popularity and proved useful (Kaikkonen et al., 2021), there are still methodological challenges and shortcomings associated with these approaches, some of which we will address.Some precision of the exposure prediction model outputs is typically lost due to the discretization of variables for the BN model (Marcot, 2017;Nojavan et al., 2017).In our case, however, the resolution of the pesticide concentration is determined by the PERPEST model, and therefore a higher resolution of exposure concentrations nodes would not have made much difference (Figure 3).
The credibility of the developed BN is influenced by the assumptions and the data derived from predictions by the process-based exposure model and case-based effect model.The RICEWQ model used to predict the exposure concentrations in the rice paddy is readily available for simulation and can be used for higher tier exposure assessment (European Commission, 2003).A detailed description of uncertainties related to this model can be found in Miao et al. (2004).We chose this model because it enabled the simulation of agricultural conditions for rice production, such as the controlled release of water, overflow, and flooding, unlike other pesticide fate and transport models.Hence, it was considered a good tool for exposure prediction for rice cultures (European Commission, 2003).However, our modeling efforts were constrained by several conditions, especially relating to the available climate model projections, which were based on a single climate model.Using multiple climate models with different greenhouse gas emission scenarios allows the integration of more variability (Fernández et al., 2017).The climate change component of the model could have been improved by using more recent and alternative climate models, as recommended by, for example, Steffens et al. (2014) and Moe et al. (2022).These authors mentioned that using an ensemble of various global and regional climate models and various greenhouse gas emission scenarios would potentially enable more robust predictions of environmental impacts such as chemical exposure in the future.
Future development of the exposure concentration prediction could also include more realistic application scenarios to run the autoRICEWQ model.Such future research could also be extended to the discharge channel to Albufera Lake by using pesticide exposure prediction models that take more hydrological processes into account, for example, the RIVWQ model, which would better allow the consideration and integration of dilution, because the RICEWQ model only predicts the exposure in the water of the rice fields.
There is also uncertainty associated with the PERPEST model that is connected to the underlying microcosm and mesocosm effects database.The database on which the casebased effect model bases its predictions has limited data availability, in particular for fish and tadpoles.Furthermore, its incorporated data is primarily based on data sets from the temperate climate such as Europe and North America (Davis et al., 2013;Van den Brink et al., 2002), and therefore is limited in its representativity for the Mediterranean climate zone.Henceforth, this limitation could be overcome by updating the database with more regionally relevant bioassays.Some other uncertainties of the PERPEST model can be associated with the input information, such as pesticide properties for the model run.Consequently, we tried to minimize this factor by using the same information source for the selected pesticides whenever possible.For example, we collected toxicity data from the ECOTOXicology Knowledgebase (Olker et al., 2022) and, subsequently, used the same method to prepare the data used on MOSAIC (Charles et al., 2017) to predict the HC50.In essence, the effect prediction model has simple data requirements, making it easy to use (Davis et al., 2013).Some uncertainty is also linked to the way the PERPEST model output is integrated into the BN, because it is a gradient, and its concentration range thus far cannot be adjusted to fit better with the exposure distribution (Figure 3).Some other restraints are pointed out by Davis et al. (2013): the PERPEST output might be challenging to use and understand for stakeholders when used for risk management due to the lack of an "established threshold risk value."To overcome this limitation of the PERPEST model, Davis et al. (2013) suggested setting acceptable probabilities.In our study, we tried a different approach to enable easy communication of BN outputs integrated with a summarizing node for effect on endpoint group and community.These nodes show the probability of any of the affected biological endpoints to be true (=effect) or false (=no effect).
In the model we have presented, any effects of the scenarios on the biological endpoints are derived through the exposure, with no direct link from the scenarios to the effect module of the network.A more direct relationship between the scenarios and the biological effects was beyond the scope of our study, but could be explored in future studies, because the combined effects of climate conditions and chemical exposure are expected to affect biological endpoints through both altered toxicity and vulnerability.A more direct link between the scenarios and the effect components could thus increase the realism of the model.
Additional model development may result in better integration and use of the effect model outputs.An updated and extended PERPEST model database would greatly decrease uncertainties.Considering further development of the effect side, the characterization of risk to a group by the joint probability of effect to any of its component can be considered rather strict.In the future, other assumptions to aggregate the endpoint groups and effect on the community level could be explored using less stringent alternatives such as proportional influence.In addition, future modeling efforts could also retain all three effect levels (no, slight, and clear effect) also at the endpoint group and community levels.However, such improvements would require some changes to the model structure and were beyond the scope of this first effort to link effect assessment with PERPEST to environmental scenarios.

CONCLUSIONS AND FUTURE OUTLOOK
Our study has shown how a BN model can be developed and used to integrate the outcomes of two different predictive models-the pesticide exposure model RICEWQ and the biological effect model PERPEST-for risk characterization in the aquatic ecosystem of a rice paddy.This approach builds on previous probabilistic models for risk assessment of the traditionally used RQ calculation, displaying uncertainty transparently of all its model components.The present study has further developed this approach by including the risk calculation for individual biological endpoints and the cumulative risk for the endpoint groups and the community.With the model and scenarios we present, we could explore patterns in the risk of the different pesticide types (herbicide, insecticide, and fungicide) to the various individual biological endpoints.Using joint probability calculation for cumulative risk to endpoint groups, we could identify the presumably indirect risk of pesticides caused by trophic interactions, such as the high probability of effects of the herbicide to invertebrates.Finally, we could evaluate the combined risk from the different endpoint groups at the community level; for example, a high overall risk of the herbicide, possibly caused by both direct effects on plants and indirect effects on invertebrates.
Future research efforts could incorporate more scenarios, such as additional crop types, application patterns, and an ensemble of climate models, to generate more realistic pictures of pesticide risks to the aquatic ecosystem.The BN model is currently limited to single-compound assessment, and does not yet consider pesticide mixtures.Nevertheless, in its current form, this model enables accounting for uncertainty of all compartments, which allows for transparency when communicating the effect of pesticides to various biological endpoints, endpoint groups, and the community in the aquatic environment.
Supporting Information-The Supporting Information is available on the Wiley Online Library at https://doi.org/10.1002/etc.5755.

FIGURE 1 :
FIGURE 1: Location of Albufera Natural Park near Valencia (red) and location of rice field clusters (colored areas).Adapted from Instituto Geográfico Nacional (2022).

FIGURE 3 :
FIGURE 3: Example of aggregation of risk from individual biological endpoints to endpoint groups, for the insecticide acetamiprid, for climate condition in 2050 and a baseline+50% application.(A) Predicted probability of the three effect classes to five invertebrate endpoints (insects, macrocrustacean, macroinvertebrate, microcrustacean, and rotifers).(B) Effect classes aggregated into Boolean states (true/false).(C) Aggregation of individual endpoints into the endpoint group Invertebrates, by the joint probability of effect on any of the endpoints (Equation1).

FIGURE 4 :
FIGURE 4: Example of the parameterized Bayesian network for the insecticide acetamiprid.It displays the predicted effect on the biological endpoints and endpoint groups for climate conditions of 2050 and a baseline (A) or baseline+50% (B) application scenario.

FIGURE 6 :
FIGURE 6: Bar chart showing the predicted probability of pesticide effect on plants, invertebrates, vertebrates, and ecosystem processes, as well as the community level for the selected pesticides, for a baseline and baseline+50% application under climate conditions in 2050.MCPA, 2-methyl-4chlorophenoxyacetic acid.

TABLE 1 :
Overview of the exposure peak concentration means and standard deviations used as input on the Bayesian network a

TABLE 2 :
Chemical, biological, and physical properties of the selected pesticides included in the PERPEST model

TABLE 3 :
Overview of biological endpoints and endpoint groups a

TABLE 4 :
Bayesian network node description containing the node name, type, number of states, and information source

TABLE 5 :
Lewis et al. (2016)d effects of the pesticides on the different endpoint groups according toLewis et al. (2016)