Development of an in-stream environmental exposure model for assessing down-the-drain chemicals in Southern Ontario

In recent years, as questions regarding legacy pollutants such as pesticides and priority contaminants have come into greater focus, scientists and regulators have been giving increasing attention to potential environmental risks associated with chemicals found in everyday products used by consumers (Boxall et al. 2012). Pharmaceuticals and formulated personal care and cleaning products such as shampoo, laundry detergent, dish detergent, and hand soap are used in high volumes and disposed of down the drain, resulting in residual releases of their ingredients from point source discharges into the environment (Keller et al. 2014). In order to estimate environmental exposures and better understand potential environmental risks, a modeling approach often is first used to provide predicted environmental concentrations (PECs) of down-the-drain chemicals. There are a number of modeling approaches that might be utilized (Keller 2006). Multimedia environmental fate models have the advantage of being a simple approach to generating PECs in various environmental compartments (e.g., air, water, sediment, soil). However, while multimedia fate models are capable of representing a variety of scales (i.e., local, regional, continental), they have generally been limited in their ability to represent spatial and temporal variability.

For chemicals whose dominant emission to the environment is via down-the-drain releases and wastewater treatment plant (WWTP) effluent discharges, a single medium approach can be coupled with geo-location of facilities. These geographic information system (GIS)-based models provide spatially explicit in-stream chemical concentration estimates. Kapo et al. (2016) discussed a number of GIS-based models that have been developed for Europe and North America. In addition, Grill et al. (2016) examined the key advantages and disadvantages of a multimedia vs. single medium modeling approach.

For regulatory purposes, the use of simple multimedia environmental fate models has been preferred. In Europe, the EUSES (European Union System for the Evaluation of Substances) model has been used for nearly two decades (Keller 2006). In the United States, the Environmental Protection Agency (EPA) Exposure and Fate Assessment Screening Tool (E-FAST) Down-the-Drain module is a screening-level model for estimating concentrations of chemicals in surface waters that may result from the disposal of consumer products into household wastewater (USEPA 2007). In Canada, Environment Canada has used MegaFlush, a spreadsheet-based proprietary WWTP effluent dilution model, to estimate aquatic environmental concentrations (Environment Canada 2009; Mackay et al. 2015). More recently, Environment and Climate Change Canada has been developing the Consumer Release Aquatic Model (CRAM), an Excel spreadsheet tool to quantify the level of aquatic exposure to a substance resulting from the use of consumer products released down the drain to wastewater treatment systems. CRAM was designed to address transparency issues associated with MegaFlush in that it works with a probabilistic distribution of dilution data that can be released to the public without concern for indirect disclosure of confidential business information (CBI). The use of a probabilistic dilution model is an approach similar to that used by the US EPA for E-FAST. In addition, Health Canada sponsored the development of a large-scale contaminant fate model for the provinces of Quebec and Ontario (Lehner et al. 2013) which was later enhanced to create the HydroROUT model (Lehner & Grill 2013). However, none of these regulatory and academic contaminant fate models for Canada are publicly available for use. As Lehner et al. (2013) noted, the Government of Canada, as part of its Chemicals Management Plan (CMP) initiated in 2006, has committed to addressing approximately 4,300 substances, prioritized through their categorization process, by 2020. Consequently, a publicly available contaminant fate model for Canada (or a meaningful portion) would be valuable for the regulated community and other interested stakeholders.

The American Cleaning Institute (ACI) developed the iSTREEM^® model as a single medium contaminant fate model for the conterminous United States to predict the concentration of individual ingredients once they have gone down the drain, accounting for physical removal and degradation (Kapo et al. 2016). Concentrations within the wastewater influent, effluent discharge in freshwater rivers and streams, and at drinking water intakes are estimated. The earliest version of iSTREEM® was developed in the 1980s and since then has been enhanced several times. iSTREEM® is intended to be used at regional scales for evaluation of broad-scale trends in chemical exposures. An additional module of Southern Ontario, covering a significant portion of the Canadian population, was recently incorporated into the pre-existing iSTREEM® model to test the feasibility of expanding iSTREEM® to be a comprehensive North American contaminant fate model. By integrating Southern Ontario within iSTREEM®, there is the additional benefit of estimating chemical loadings in the cross-boundary US–Canadian Great Lakes watershed. This paper details the development of the lower St. Lawrence watershed environmental exposure model.

In order to corroborate the model's results and understand its limitations when applied to Canadian watersheds, case studies were conducted comparing measured surface water concentrations to iSTREEM® PECs for the chemicals carbamazepine and triclosan. Samples in this dataset were collected between 2002 and 2014 for triclosan, and 2000 to 2012 for carbamazepine. Triclosan is an antimicrobial compound used in personal care products and carbamazepine is a pharmaceutical used in the treatment of epilepsy and trigeminal neuralgia. These chemicals were selected because disposal is down the drain and data on usage estimates and environmental monitoring were available. Peer-reviewed literature and government reports were mined, and data were extracted specific to Southern Ontario. Measured surface water concentrations were aggregated and compared to iSTREEM® predicted values during mean and low flow.

Ontario is the most populous province in Canada having over 13.7 million people in 2015. It was selected as the first providence in Canada to be incorporated into iSTREEM® due to availability of data, its high population density and the cross-border watershed shared with the USA. The geographic scope was further confined to the lower St. Lawrence drainage basin, where approximately 90% of the population in Ontario is connected to a wastewater treatment facility according to the Ontario Ministry of Environment (Day, W., personal communication). The St. Lawrence drainage basin spans 1.6 million km² from the Great Lakes into southeastern Quebec and northeastern United States (Environment Canada 2013). The study area border was further refined based on sub-watershed boundaries using data available through Natural Resources Canada (2011). The region selected for inclusion within iSTREEM® ranges approximately from the Detroit River in the South to Lake Temagami in the North (Figure 1).

The underlying data on hydrology, elevation, dam locations, wastewater treatment facilities, and drinking water intakes were obtained from Canadian Government sources (detailed in this section). Elevation and dam locations were used to aid in the development of the river network by establishing flow direction and locations of flow obstruction. Hydrology, wastewater, and drinking water facility attributes and coordinates provide the spatial framework for the iSTREEM® model functions.

Only freshwater streams, rivers, and inland lakes were included in this modeling exercise. These water bodies have much lower dilution than the Great Lakes, which were excluded. The Atlas of Canada 1,000,000 National Framework Data, Hydrology – Drainage Network – St. Lawrence (Natural Resources Canada 2006) scale 1:1,000,000 provided the basis of the river network. The scale of iSTREEM® US hydrology was 1:500,000, therefore the 1:1,000,000 Canadian hydrology dataset was selected since it potentially could be integrated with more ease into the US model than the National Hydro Network, which has a resolution of 1:50,000, or better. This was a consideration given the possibility of potentially integrating the US–Canadian cross-boundary watersheds in the future.

The stream segments were spliced according to natural or anthropogenic breaks such as a water treatment facility or dam, consistent with the existing US model. The segments are an important component to the model algorithm as upstream loadings to downstream segments are considered in the computation of segment-specific PECs. The network was trimmed to the selected study area in the lower region of the St. Lawrence basin by use of a GIS file publicly available from the Canadian Government (Natural Resources Canada 2011). Mean and low flow (7Q10, i.e., the lowest 7-day average flow that occurs (on average) once every 10 years) values were primarily determined by use of the Ontario Flow Assessment Tools, version III. Flow was unable to be determined for certain waterbodies using the OFAT III model; in those instances, mean flows were estimated by using gauge station data where available (Environment Canada 2012), interpolation, or by use of the McMaster University MAC-HBV model. Interpolation was conducted based on an established protocol from a previous exercise conducted in development of the US portion of iSTREEM® (Summary S1 of Kapo et al. 2016).

Flow could not be determined for all stream and river segments due to lack of available data. This accounted for approximately 6% of segments at mean flow and 12% in low flow conditions. In these cases, to be conservative, the lowest flow value the algorithm will allow, 2.83 × 10⁻⁶ m³/s (0.0001 ft³/s) was assigned.

WWTP information was provided by the Ontario Ministry of Environment (Day, W., personal communication). Pertinent information for the algorithm included general facility coordinates (latitude and longitude), treatment technologies (e.g., activated sludge), average system flow, and population served. Effluent discharge locations were obtained through manual extrapolation using satellite imagery and by following a pre-established protocol to ensure consistency (Table 1). Extrapolation tools included ArcGIS version 10.0, Google Earth, Google Maps, and Bing Maps. Facilities that land applied effluent were excluded, as were those facilities that discharged directly into the Great Lakes. Lakes, streams, and rivers are lower dilution systems and therefore are more likely to be adversely impacted from trace chemicals in effluent discharge. Dilution into the Great Lakes is assumed to be much higher than in smaller waterbodies. The dataset was reviewed for accuracy and further reduced. Facilities that were more than three miles away from a river segment in the National Atlas of Canada hydrologic network were removed. WWTPs located at the end of a terminal reach were excluded because they were not applicable to the model, as there were no further streams to contribute to the chemical load. The resulting facility dataset included 188 Ontario WWTPs which were integrated into the model.

Drinking water facility locations and sources were obtained from the Ontario Ministry of the Environment Drinking Water Surveillance Program (Ontario Ministry of the Environment 2011) and a Google Earth File (.KMZ) formerly available on the Ministry's website (Ontario Ministry of the Environment 2010). The datasets were compared to ensure completeness. Data were analyzed and filtered in ArcGIS version 10.0 to include only intakes that were part of a municipal system sourced from surface water within the study area but not from a Great Lake. Any facilities on waterways that were not effluent impacted were removed from the dataset. The iSTREEM® model is focused on residential drinking water systems, therefore, non-residential drinking water systems were excluded. Many of the drinking water facilities remaining were part of larger distribution systems. These systems obtained water from other communities and were not included.

Since exact coordinates of intakes were not provided, best estimates of the intake locations were determined by internet queries to find drinking water treatment facility addresses or information pertaining to the intake itself. All information gathered was substantiated visually using satellite imagery in ArcGIS, Google Earth, or Google Maps. After the location of the water treatment facility was identified, coordinates were established in the waterbody to estimate the location of the intake pipe, and recorded. The intake location was inferred by a combination of publicly available information and satellite imagery. For remaining facilities, the township or facility operator were contacted to obtain an approximate location of the intake. Only two facilities were not included in the dataset due to the inability to establish intake coordinates. Ultimately, a total of 59 drinking water intakes were incorporated into iSTREEM® for the St. Lawrence region.

A literature search was conducted for surface water monitoring data for Southern Ontario. Triclosan monitoring data were extracted from Arlos et al. (2015), Health Canada & Environment Canada (2012), and Environment and Climate Change Canada & Health Canada (2016). As a result, 561 points from the monitoring data were identified for triclosan and were representative of the iSTREEM® geography. Water samples were collected from years ranging from 2002 to 2014. Data points reported as less than the method of detection limit (MDL; 284) were assigned the MDL, when reported. In those cases where the MDL was not reported, an average of the MDLs from the various studies was used. For carbamazepine, monitoring data from sites in Southern Ontario were extracted from Arlos et al. (2015), Metcalfe et al. (2003), Li et al. (2010), and Hao et al. (2006). The majority of those sites were within the Grand River watershed. Out of 49 samples, two had values below the MDL and another four had values below the limit of quantification (LOQ). In those cases, the MDL or LOQ value was assigned. For both the triclosan and carbamazepine datasets, the authors of the studies were contacted to obtain raw data if they were not publicly available.

The input parameters necessary to run simulations of the model require the investigator to have basic knowledge of the use and fate of the chemical of interest. The investigator must supply values for removal efficiency of the chemical for each of the six wastewater treatment technology types identified in the model and the national per capita use volume. In order to calculate in-stream attenuation beyond dilution by the receiving waterbody, the investigator must supply an in-stream decay rate.

iSTREEM® simulations for Southern Ontario were conducted for the chemicals triclosan and carbamazepine in order to compare the modeled output (distribution of PECs) to environmental monitoring data. iSTREEM® simulations were run at mean and low flow scenarios using parameters detailed in Table 2. Measured concentrations of triclosan in surface waters were available for the periods 2002–2005 and 2012–2014. The source of the Canadian per capita use rate of triclosan for the iSTREEM® simulations was obtained from an assessment by the Canadian Government (Health Canada & Environment Canada 2012) which estimated the value based on influent concentrations from across the country sometime between 2002 and 2011. This per capita use rate is a factor of eight lower than that used by de Zwart et al. (2006) in their model estimation of riverine triclosan concentrations in Ohio. However, the difference in use rate between Canada and the United States appears to be corroborated when the estimate by Health Canada & Environment Canada (2012) is compared to that reported for the USA by Perez et al. (2013) in a meta-analysis of triclosan environmental occurrence data covering 1999–2012.

Measured concentration data for carbamazepine were available from 2000 to 2006 and 2012. The per capita use rate of carbamazepine was calculated from data reported in Lehner et al. (2013) based on reported use of the drug, taking into account metabolic loss. Simulation parameters for treatability of carbamazepine in activated sludge WWTPs and oxidation ditches, and the in-stream decay rate were derived from Cunningham et al. (2010), and a conservative assumption of no removal by primary or secondary (RBC, trickling filter, and lagoons) was applied.

The model output contains four distinct downloadable database files (.dbf) and a Microsoft Access database (.mdb). Each file contains specific information on the WWTPs, DWIs, rivers, and user input parameters. The Microsoft Access database is a conglomerate of the output files with the addition of spatial data, such as locations of WWTPs and DWIs, hydrology and flows. This allows users to view or analyze the results in a GIS. The PECs are calculated for WWTP effluent, DWI at the point of the intake, at the upstream and downstream portion of each stream segment, and the average PEC for each segment. Besides providing PECs, attributes such as population served per DWI or WWTP, and the technology type per WWTP are included in the output.

The average river reach PECs were used for comparison against monitoring data given the measured data contained averages of replicate samples and since it was not possible to match sampling location against particular stream segment junctures.

Model output from iSTREEM® simulations for triclosan and carbamazepine included estimated in-stream concentrations (μg/L) of each chemical at mean annual and low (7Q10) flow for 11,700 km of effluent-impacted streams in Southern Ontario. Cumulative distributions of iSTREEM® modeled concentrations (based on the total number of impacted river segments) and monitoring data were compared for each flow scenario. The entire distributions for both triclosan and carbamazepine can be found in the Supplementary materials, Figures S1–S5 (available with the online version of this paper). In Figure 2, measured concentrations of carbamazepine from years 2000 to 2012 were compared to the iSTREEM® simulations at mean and low flow scenarios. Above the 14th-centile, the data points reported are measured detected values.

The monitoring data generally fell between the predicted low and predicted mean flow concentrations. The maximum concentration of carbamazepine reported in the monitoring data was 0.329 μg/L. The 95th-centile concentration predicted by iSTREEM® at mean and low flow were 0.0960 μg/L and 0.7386 μg/L, respectively.

Approximately 51% of measured triclosan concentrations were below the MDL, resulting in the appearance of ‘stacked’ vertical sections in the distribution of measured points on the graph. This high frequency of non-detect values for triclosan limits the ability to directly compare the modeled and measured distributions; however, it can be assumed that the ‘true’ distribution shape of the measured data would be shifted to the left, and thus would likely fall between the modeled mean and low flow distributions. This trend is more apparent when observing the distributions above the 80th-centile of measured values where only concentrations that had been detected in source water remain (Figure 3).

The concentrations predicted by the iSTREEM® simulation at mean flow were consistently lower than the measured surface water concentrations. The concentrations predicted by the iSTREEM® simulation at low flow, were generally higher than the measured concentrations. The highest measured concentration of triclosan in this dataset was 0.691 μg/L, measured in 2002. The 95th-centile concentrations predicted by iSTREEM® at mean flow and low flow were 0.0133 μg/L and 0.0939 μg/L, respectively. Above the 99th-centile, the iSTREEM® model (particularly at low flow) highly overestimates the predicted concentrations. This is expected, particularly in low flow scenarios, due to the model's overestimation of the high-end values, as detailed in Kapo et al. (2016).

The case studies suggest that the model is highly conservative above the 99th-centile. Qualitatively, the measured data distributions for both chemicals did not consistently follow the trend of mean or low flow, but the measured distributions for both chemicals generally fell between their respective modeled mean and low flow distributions. Methodologies, seasonality, and year of sampling did not appear to have an obvious pattern among the monitoring studies used as sources. Additionally, since both chemicals were measured over a 12-year span, it is possible there was a change in the volume used by consumers over that duration or improvement in removal technologies. Either change would affect how well the measured data distributions compared to the PECs distributions. To investigate whether measured concentrations of triclosan changed over time, the concentrations were grouped by year, from 2002 to 2005 and from 2012 to 2014. Descriptive statistics were run for both datasets to determine their similarities with non-detects included (Table 3).

Differences in measured surface water concentrations for triclosan over time were apparent visually (Figure 4), with an overall decreasing trend in measured concentrations of the time period. A two-sample t-test assuming unequal variances confirmed this (P < 0.01; α = 0.05). This change in measured surface water concentrations of triclosan concentrations could be due to wastewater treatment technology enhancements or changes in consumer use.

There appears to be a temporal difference in environmental exposure for carbamazepine based on measured data from 2000 to 2006 compared to 2012 (Figure 5). However, no statistical difference was found when comparing the sampling years of 2000–2006 to 2012 by the Wilcoxon rank sum test with continuity correction. Since there were only 49 samples, it is possible a larger dataset may show a significant difference.

In both case studies, iSTREEM® low flow predictions were generally conservative when compared to measured concentrations sampled in more recent years in surface waters. The iSTREEM® model was designed to function as a screening tool to evaluate broad-scale exposure trends (distributions). This study demonstrates that despite the high variability and uncertainty in the monitoring data, the model adequately meets these objectives.

The creation of the St. Lawrence environmental exposure model was the first step in the development of a complete US–Canadian model. Future work to develop the model entails expanding east to encompass the St. Lawrence watershed in Québec. This expansion would include Montréal, Ottawa, Québec City, and the associated watersheds, thus significantly increasing the population coverage of the model. Additional work will need to occur to integrate the St. Lawrence drainage basin into the existing iSTREEM® US model. Since iSTREEM® US has been recently upgraded using the hydrology dataset of the National Hydrologic Network Plus to achieve a resolution of 1:50,000, further advancement of the Canadian model could also include upgrading the spatial data to the National Hydro Network.

A more robust environmental monitoring dataset would aid in comparison studies, as limited data were available for the case studies presented in this paper. Future collaborations with investigators in the field would serve to assist in model corroboration by providing access to monitoring data results.

iSTREEM® is a practical screening tool for environmental exposure assessment that could be used for regulatory assessments, such as those occurring under the Canadian Chemical Management Plan. For research purposes, it allows investigators to narrow project scope (e.g., field work, further modeling) to particular geographic locations using the PECs derived by the model. The St. Lawrence drainage basin iSTREEM® model is the first step in integrating a cross-border US–Canadian model to better estimate residual chemical concentrations from upstream loading. Based on the presented case studies for triclosan and carbamazepine, the PECs and monitoring data were found to generally yield similar results. However, limited environmental monitoring data were available for carbamazepine in the selected geography. As new monitoring data become available further studies comparing the PECs to measured environmental concentrations could be conducted to help deepen our understanding as to the model's limitations and strengths.

The authors would like to acknowledge Katherine Kapo for her thorough review of this work. We would also like to extend our gratitude to Xinhao Wang, Raghu Vamshi, and Chris Holmes for their technical expertise in the development of this model and to the Canadian agencies and investigators who provided us access to the raw data in their studies.