Characterizing the regional concentrations and seasonality of the emerging pollutant decamethylcyclopentasiloxane (D5) using a WRF+CHIMERE modeling approach

Modeling approaches are important tools for understanding environmental dynamics from a local to a global scale. Field-based data are often difficult and expensive to obtain, so computational simulations have been used to try to replicate, among other things, the physical and chemical phenomena of the behavior of contaminants. However, models will always need a strong observational backup in order to evaluate their findings properly. At the same time, models can also optimize the sampling strategies, creating an increasing synergy to reach the best representations of the “real world” while reducing the costs.

Naturally, the need for information is even greater when the target chemicals (often called emerging compounds) have only recently started to attract the attention of the scientific community. This is the case of volatile methyl siloxanes (VMSs), molecules with alternating atoms of silicon (Si) and oxygen (O) and methyl side chains along the backbone, arranged as linear or cyclic structures (Rucker and Kümmerer, 2015). They are extensively used in industrial applications, personal care, and other consumer products such as lubricants, moisturizers, and even food additives (Graiver et al., 2003; Wang et al., 2009; Dudzina et al., 2014). They are also found in some common formulations such as dimethicone and cyclomethicone, representing linear and cyclic VMSs, respectively. These chemicals present high stability, inertness, and easy release into the environment (Wilson et al., 2005).

With half-lives in air ≥2 days, VMSs are prone to long-range atmospheric transport (LRAT), with D5 one of such chemicals reported by the European Chemical Agency (European Chemicals Agency [ECHA], 2012), as well as to bioaccumulation and biomagnification in living organisms due to their lipophilic nature (Homem et al., 2017). These issues are, however, controversial in the existing literature. For instance, the presence of VMSs in polar regions is reported (Krogseth et al., 2013a; Sanchís et al., 2015), but some authors relate this evidence to the emissions of local settlements rather than to LRAT (Warner et al., 2010). As for bioaccumulation, the conclusions differ depending on the type of environment. While a review by Wang et al. (2013) reported no biomagnification of D4 and D5 in aquatic food webs, some authors concluded that cyclic VMSs bioaccumulate in biota such as worms and flounders (Kierkegaard et al., 2011) and freshwater and marine food webs (Borgå et al., 2013; Jia et al., 2015). The potentially toxic behavior is thus also a matter for discussion, but there are still considerable gaps in the current knowledge, and the literature shows evidence of risk, especially for cyclic congeners (Quinn et al., 2007; Brooke et al., 2009a, 2009b; Dekant and Klauning, 2016, among others). Recommendations have been made urging the development of an effective methodology to evaluate the risks to humans and the environment of D5 and related chemicals (Mackay et al., 2015). From all the existing data for D4, D5, and D6, the conclusion of Redman et al. (2012) that their levels in biota were not likely to have negative effects was later corroborated for D5 by Fairbrother et al. (2015). But in a 2016 report, ECHA considered D4 persistent, bioaccumulative, and toxic and D5 very persistent and very bioaccumulative (ECHA, 2016). This led to the recent restriction of these two chemicals within the EU (effective from February 1, 2020), limiting their content in wash-off personal care products (PCPs) to 0.1% (European Commission, 2018). The extension of these restrictions to D6, leave-on products, and other consumer and professional products is currently under evaluation (ECHA, 2020).

Hence, monitoring studies in several environmental matrices are being conducted to enhance the still scarce field-based data on these contaminants and to enhance the knowledge of their harmful effects. Project SILOQUEST, one of these studies, aimed to assess the concentrations of 8 VMSs in air, vegetation, and soils in potentially relevant sites (Ratola et al., 2016). As inhalation is the principal way humans are exposed to VMSs (Rucker and Kümmerer, 2015), this project also included a modeling approach to the behavior of airborne D5 (the main congener in the common siloxane formulations), which the current study reports.

In terms of modeling, some approaches have used the Danish Eulerian Hemispheric Model (McLachlan et al., 2010), the GloboPOP model (Xu and Wania, 2013), or the BETR (Berkeley-Trent) model (MacLeod et al., 2011). These contributions deal predominantly with regional and global domains and do not rely so much on solving dynamically the climate physics but on solving the chemical equations ruling the life cycle of the chemicals. More recently, Janechek et al. (2017) used a chemistry transport model (CTM) as the Community Multiscale Air Quality (CMAQ) model to estimate the VMS atmospheric distribution, mainly for D5, on a regional domain. Hence, our idea follows the strategy of the latter contribution, associating a CTM as CHIMERE with an online approach (where meteorology/climate equations governing the transport of pollutants are coupled online with the D5 chemistry), and applying it to a local scale, such as a domain covering Portugal and the Galicia region of Spain. A regional atmospheric modeling system was thus developed to simulate the presence of D5 in this domain. The ultimate goal is to complement measurement data and help build greater knowledge of VMSs in the environment.

As mentioned above, the data used for this study resulted from a previously performed sampling strategy (Ratola et al., 2016). In brief, a total of 8 sampling sites were chosen (see map in Figure 1), covering different anthropogenic pressure (i.e., population size/density and emission patterns) and seasonal variety. Two sites, Porto (3) and Braga (2) are located in large urban settings prone to the impact of intense human activity. Two others, Z.I. Mota (6) and Estarreja (5), are based in industrial parks, whose activities represent possible sources of siloxanes. Two sites in Portugal were chosen to assess the variability related to urban presence in different seasons: Esmoriz (4) on the west coast in the north and Praia da Luz (8) on the south coast. Both are beach resorts, which can be hotspots for VMS release due to the extensive use of sunblocks and other PCPs. Finally, the more remote areas of Outeiro (1) and Midões (7) were chosen to evaluate background concentrations. The target domain (Figure 1) was thus chosen to represent as closely as possible the sampling points available, under a logical regional perspective, including the whole Portuguese continental territory and Galicia. It is always desirable to have as many sampling points as possible to enable better evaluation by the modeling tools and a more accurate representation of reality, but the analysis of contaminants (especially those of emerging concern, and often because of the lack of standards or rapid protocols) comes with logistic and operational costs that need to be addressed. The sites chosen in this study intend to highlight the areas of potentially different impacts of VMSs, covering the widest domain possible. Models can also be important in defining other hotspots and further optimizing the distribution of sampling points.

The sampling setup used to reach the results is also reported by Ratola et al. (2016). In brief, sorbent-impregnated polyurethane (SIP) foam disks (14 cm in diameter, 1.35 cm thick, surface area: 365 cm², mass: 4.40 g, volume: 207 cm³, and density: 0.0213 g cm^–3) protected by two stainless steel domes were placed at each site from 1 to 2 m above the ground. In each of the sampling points and campaigns, two SIPs (for duplicate analysis) were deployed, and a field blank assay was done, namely a third SIP disk was exposed but only during the deployment of the other two and taken to the lab for analysis afterward. The sorbent used was XAD-4, a polymeric material with a large uptake capacity due to its very high surface area (750 m² g^–1; Wania and Shunthirasingham, 2020) and commonly used in siloxane studies (e.g., Genualdi et al., 2011; Ahrens et al., 2014; Ratola et al., 2016; Rauert et al., 2018). The exposure time is a crucial parameter and depends greatly on the expected concentrations of the environmental contaminants, the sorption capacity of the SIP, and the meteorological conditions. Based on the group’s expertise and on the related literature (Genualdi et al., 2011; Ahrens et al., 2014), a deployment time of 3 months was chosen. The stability of D5 in the sorbent media could be a limitation, as the possibility of degradation into D3 and D4 during sampling was suggested (Krogseth et al., 2013a), but this is more of a concern when polyurethane foam disks alone are used, since XAD sorbents increase the general stability of semivolatile organic compounds (Wania and Shunthirasingham, 2020). Shorter periods could have been chosen with different configurations or sorbents (Krogseth et al., 2013b), but 3 months also coincides with the period commonly adopted in chemistry transport and climatic modeling for a year-round seasonal behavior, with the following division: December–January–February, March–April–May, June–July–August, and September–October–November. As the intention was to study the common climatic extremes in southern Europe latitudes (winter and summer), two campaigns were performed: one in winter, from the last week of November 2013 to the last week of February 2014; the other in summer, from the first week of June to the first week of September 2014.

The analytical protocol for the analysis and quantification of D5 in SIP disks is detailed elsewhere (Ramos et al., 2016). In brief, SIPs were extracted overnight with a mixture of hexane and dichloromethane (1:1) in a 200-mL Soxhlet extractor. The extract was passed through a Na₂SO₄-filled glass column and reduced to approximately 1 mL by rotary evaporation, transferred to gas chromatography (GC) vials, evaporated to near dryness by N₂, and redissolved in 150 µL of hexane before analysis. Quantification was done by GC in a Varian 4000 GC/mass spectroscopy system (Walnut Creek, CA, USA), with a CP-1177 split/splitless injector adapted with a Merlin Microseal System and a special low-bleed Agilent DB-5 ms ultra-inert column (30 m × 0.25 mm I.D., 0.25 µm film thickness) to minimize external contaminations due to possible siloxane bleeding from the equipment. Chromatographic separation was done with helium at 1 mL min^–1 as carrier gas and the following GC oven temperature program: starting at 35ºC (held for 5 min), then raised to 160ºC at 10ºC min^–1. The temperature of the injector was 200ºC using a liner without glass wool and an injection volume of 1 µL in splitless mode. The temperatures of the transfer line and ion source were 250ºC and 200ºC, respectively. Detection of the target compound was performed using the ion-trap mass spectrometer operating in electron ionization mode (70 eV) and time-scheduled selected ion storage acquisition. M4Q was used as the internal standard, and the quantification of the analytes was done with the Mass Spectrometry Workstation 6.6 software from Varian.

In terms of method validation, a linearity range from 1 to 500 µg L^–1 was obtained and the limit of detection (calculated using the signal-to-noise ratio of 3) for D5 was 0.02 ng m³. The mean recovery (accuracy) for D5 analysis was 91 ± 5%, and intraday precision was 11% (n = 5). A strict quality assurance/quality control protocol was put in place throughout the analysis. Indeed, all laboratory personnel avoided the use of PCPs not only in the laboratory but also during the sampling campaigns, as VMSs (in particular, D5) may be present in their formulations. Moreover, during volume reduction, the extracts were never completely dried in order to reduce D5 loss by evaporation, and all glassware or related material was baked at 450ºC and rinsed with acetone and hexane before use. Apart from the field blanks, procedural blanks were also performed regularly in the laboratory to monitor external contamination. Being ubiquitous, D5 was, as usual, found in the blanks. However, the values never exceeded 20% of the actual concentrations at the sites. The results were corrected accordingly, subtracting the values found in the blank assays.

The results of the field campaigns are presented in Supplemental Material (Table S1) and show low relative errors between duplicates (from 3% to 23%, with a mean of 11%). The concentration values were converted from ng SIP^–1 to ng m^–3 considering the air sampling rates for D5, following an Excel template (http://dx.doi.org/10.13140/RG.2.1.3998.8884) created and regularly updated by Professor Tom Harner’s research group at Environment and Climate Change Canada, with input data from Xu et al. (2007) for VMSs. This tool is used and accepted worldwide in the very difficult task of estimating sample air volumes for a wide range of priority chemicals for both polyurethane foam and SIP passive air samplers (e.g., Francisco et al., 2017; Bohlin-Nizzetto et al., 2020), taking into account the molecular and fugacity properties of each chemical. It requires the input of the deployment time and the average temperature in that period. In our case, these data were obtained from the nearest meteorological stations (IPMA—Portuguese Institute for Sea and Atmosphere and www.wunderground.com). This tool is freely accessible by the scientific community and its updates are available on a regular basis.

The meteorological model Weather Research and Forecast—Advanced Research WRF (WRF-ARW) model v3.6.1 (Klemp et al., 2007; Skamarock et al., 2008) was coupled offline to the CHIMERE CTM (Menut et al., 2013) following the modeling framework defined in Ratola and Jiménez-Guerrero (2016). The approach followed here is analogous to that of Janechek et al. (2017). These authors used the WRF-ARW model to drive a different CTM (CMAQ) applied to the domain of North America (including the United States, Southern Canada, and Northern Mexico), to three cyclic VMSs (D4, D5, and D6), and some oxidation products for the year 2004. The contribution of dry deposition to the removal of VMSs from the atmosphere was considered in Janechek et al. (2017), but it is deemed small and was not considered in our study.

WRF is a fully compressible, Eulerian nonhydrostatic meteorological model. The WRF+CHIMERE modeling system was designed to cover all the field campaign domains with high resolution (Portugal). The simulation strategy is thus based on three nested domains covering most of Europe (resolution of 27 km), the Iberian Peninsula (9 km), and Portugal and Galicia (3 km, our target domain) with 33 sigma layers covering from the ground up to 10 hPa. The simulations cover the entire years 2013 and 2014. The temporal resolution of the model simulations is 1 h; however, the output data were extracted and postprocessed for exactly the same periods as the field campaigns. Figure S1 in the Supplemental Material shows the terrain height, land use index, and vegetation fraction included for this domain.

The physicochemical options for the regional modeling system are summarized in Table 1. The WRF simulations are driven every 6 h by ERA-Interim reanalysis (Dee et al., 2011), while the climatological boundary conditions for the CHIMERE CTM are based on the INteraction with Chemistry and Aerosols chemistry and aerosol model coupled to the Laboratoire de Météorologie Dynamique (LMD) General Circulation Model, LMDz (LMDz-INCA; Folberth et al., 2006). Since the contribution of the transport from outside the model domain to ground-level concentrations (those considered in this work) can be assumed to be negligible when compared to the contribution of local emission sources, the influence of using climatological boundary conditions is limited and overwhelmed by local processes (Jiménez-Guerrero et al., 2012; Jerez et al., 2013).

The overall atmospheric fate of VMS can be simplified to 3 complex processes: emissions, advection, and reaction with OH radicals (McLachlan, 2018). Since D5 is present in the atmosphere almost exclusively in the gas phase, the latter process is predominant in its degradation. Consequently, the MELCHIOR2 gas-phase mechanism (Derognat et al., 2003) is implemented within CHIMERE, modified to include the D5 reaction with OH radicals. Atkinson (1991) measured reaction rate constants of D5 with the OH radical of (1.55 ± 0.17) × 10^–12 cm³ molecule^–1 s^–1, a value also reported by Kierkegaard and McLachlan (2013) and Janechek et al. (2017). However, more recent studies estimate a degradation rate constant of (2.6 ± 0.3) × 10^–12 cm³ molecule^–1 s^–1 (Safron et al., 2015), the one considered in this study and in agreement with that found by Xiao et al. (2015), 2.46 (2.74 – 2.20, confidence interval [CI] 95%) × 10^–12 cm³ molecule^–1 s^–1. The D5 reaction with OH accounts for the main transformations of siloxanes in the atmosphere, whereas other phenomena such as photochemical reactions with NO₃ radicals and ozone are minimal and can be considered negligible for our purpose (Gaj and Pakuluk, 2015). The activation energy suggested by Safron et al. (2015) was used to correct the degradation rate constant for temperature (4.3 ± 2.8 kJ mol^–1, approximately 1.02 kcal mol^–1), which is also in the same range as the values found by Xiao et al. (2015), 0.79 (1.38 – 0.21, CI 95%) kcal mol^–1. Alton and Browne (2020) reported yet another value for the OH degradation rate, but the value of (2.1 ± 0.1) × 10^–12 cm³ molecule^–1 s^–1 is very close to the one chosen in this work. The Cl-initiated oxidation referred to by Alton and Browne (2020) was not taken into account in our study and could have some influence on the results. With most sampling points in coastal areas, there is the possibility of a higher Cl presence in the atmosphere converted from chloride originating from sea salt particles (Chang et al., 2004). Moreover, the presence of Cl can also be higher in large urban areas—Alton and Browne (2020) note that Cl atoms could account for approximately 25%–30% of D5 loss in urban areas, following measurements in Boulder (United States) and Toronto (Canada). Even if the authors report that the impact of Cl atoms on global VMS lifetime is considered reduced, this aspect should be taken into account in future studies.

A full description of homogeneous and heterogeneous chemistry in the atmosphere, including gas-phase and aerosols pollutants, is needed for an accurate characterization of the OH radical. Therefore, although D5 does not interact with aerosols, the inclusion of both anthropogenic and natural aerosols is relevant for an accurate modeling of D5. Henceforth, in addition to gas-phase chemistry, aerosol and heterogeneous chemistry is included within CHIMERE. Aerosol chemistry distinguishes between different chemical aerosol components, namely nitrate, sulfate, ammonium, elemental and organic carbon with three subcomponents (primary, secondary anthropogenic and secondary biogenic) and marine aerosols.

However, the most problematic issue when implementing the modeling strategy for this study was the assessment of the emissions. Non-D5 anthropogenic emissions for the entire period of simulations are derived from the European Monitoring and Evaluation Programme (EMEP) database (Vestreng et al., 2009). Terrestrial biogenic emissions (Guenther et al., 2006) depend on meteorological conditions and are consequently coupled to WRF outputs. With respect to D5, most emission data are nonexistent or confidential, as can be seen in an appraisal of D5 emissions by the UK Environment Agency at a local, regional, and continental level (Brooke et al., 2009a). In any case, from the information available (see Table S2 in Supplemental Material), it can be confirmed that the formulations and use of PCPs are the main emitting sources to air and wastewater, followed by polydimethylsiloxane residues and degradation. Given that the quantification of the VMS emissions was impossible for our domain of study (Portugal and Galicia) due to the lack of data, the option was to find and apply the emission rates to air available in the literature. In total, five values were found, all from northern latitudes: McLachlan et al. (2010): 135 mg capita^–1 d^–1; Navea et al. (2011): 140 mg capita^–1 d^–1; Buser et al. (2014): 190 mg capita^–1 d^–1; Montemayor et al. (2013): 260 mg capita^–1 d^–1; and Buser et al. (2013): 310 mg capita^–1 d^–1.

Taking into account the uncertainties typically observed when estimating emissions for organic contaminants (e.g., Breivik and Alcock, 2002; Breivik et al., 2012), this may not be considered a wide range, but the lowest rate is still less than half the highest. In this study, in order to cover the extremes of uncertainty, the emission factors by McLachlan et al. (2010; MC10) and Buser et al. (2013; BU13) were selected and multiplied by the population density in each cell of the target area to quantify the emissions per cell. In this way, an indication of the emission rate that would be more suitable for our domain is also possible. The hypothesis of constant population density (no seasonal variations) and emissions factors per capita (always scaled according to the population distribution in the domain) for summer and winter (as done by Navea et al., 2011) were followed in this approach to estimate the gridded emissions of D5. Janechek et al. (2017) also distributed the emissions of VMS with the gridded population of their domain (United States, Canada, and Mexico), estimating a value for D5 by combining antiperspirant sales with consumption data for each country, as was done by McLachlan et al. (2010).

The simulations for the D5 concentrations are shown in Figure 2a and indicate a predominant distribution along the more densely populated coastal areas, with a hotspot in Galicia, where concentrations surpass 15 ng m^–3 when using the MC10 emission factor or 30 ng m^–3 if the BU13 emissions are selected. Moreover, the north of Portugal has a stronger incidence of D5 than the south (below the Lisbon area), which indicates a marked north-south gradient in D5 concentrations in the target domain. The population density (higher in the north) may be one of the reasons, linked with the presence of more sources of VMS emissions, including industry. Two general patterns can be seen in Figure 2a: (1) the ground concentrations of D5 simulated under MC10 emissions are lower by a factor of about 2 when compared with the highest BU13 emission rates, which is to be expected, given the higher emissions for MC10 compared to BU13; and (2) the concentrations are generally higher in winter than in summer for both emission scenarios.

This last result may have two main causes. The first might be related to the variable OH radical concentrations in summer and winter (higher in summer), which are the dominant sink for D5 in the atmosphere (Janechek et al., 2017). During winter, the OH radical-driven atmospheric degradation is reduced due to the typically lower OH concentrations (Stone et al., 2012; Warner et al., 2020). In this sense, one important question arises, concerning the weight of the aforementioned D5 loss process by OH at the scale of the investigated domain. To address it, a sensitivity analysis was conducted for the winter and summer periods (Figure 3). In those sensitivity analyses, the D5 degradation reaction by OH was turned off in order to assess how relevant the OH chemistry is on a regional scale. In general, the OH sinks play a negligible role on the ground-level concentrations of D5, with average contributions under 0.5% in most parts of the modeling domain (Figure 3). The low relative difference between the simulations with OH degradation chemistry turned off and turned on at the peak D5 sites (urban areas) could be ascribed to the low OH concentrations caused by urban OH sinks and is consistent with low modeled surface OH (Figure S2 in the Supplemental Material). Henceforth, OH loss processes are larger for summertime over rural areas, far from the emitting sources (in those areas, the D5 concentrations are dominated by the local emissions), reaching values up to 2% for summertime in the simulation using Buser’s emissions (BU13). This agrees with the results of Janechek et al. (2017) who state that rural sites are more dependent on lifetime with respect to reaction with OH, while source locations are less sensitive. Consequently, the oxidative loss of D5 has a very limited influence in D5 concentrations, the chemical loss is slow in relation to the atmospheric transport and therefore that process will play a minor role in determining the performance of the model. These results agree with the current scientific literature. For instance, Janechek et al. (2017) assessed the concentration of the D5 oxidation products (called generically o-D5), finding that o-D5 products are significantly lower than the parent compounds (mean OH loss processes of 2.08%). Likewise, Ferracci et al. (2018) reported that OH sinks have a relatively negligible impact on atmospheric trace chemicals such as D5, which has a relatively long atmospheric half-life (ECHA, 2016). Finally, Sakurai et al. (2019) estimated that in the D5 losses by chemical reactions in the atmosphere of the Tokyo Bay are two orders of magnitude lower than advection processes, representing only a 2.3% share.

The second cause for higher winter concentrations could be related to the behavior of the planetary boundary layer (PBL). The results shown in Figure 4 depict a generally higher PBL in summer for most of the domain, allowing a stronger dilution factor of the chemicals in the lower atmosphere (McLachlan et al., 2010; Krogseth et al., 2013b; Ahrens et al., 2014). Again, these causes do not include the opposite effect of the local emission patterns and population variations, suggested as reasons for the higher D5 concentrations found in summer in several of the passive air sampling sites of this study (e.g., Porto and Z.I. Mota). Nevertheless, from the maps in Figure 2a, it can be seen that for the Alentejo and Algarve coastlines (the southernmost areas of Portugal), the simulated concentrations exhibit a slight increase in summer with respect to the winter, as is also found with the sampling site of Praia da Luz (in the Algarve). These changes are likely attributed to seasonal meteorological differences, since population density (and hence, emissions) are assumed to be constant year-round in each cell of the modeling domain. In this sense, Figure 4 depicts that, over coastal areas located in the target domain, the growth of the summer PBL may be laminated by the stronger sea breezes, as also reported by Pérez et al. (2006). This is a consequence of the intensified contrast between the cool ocean and the warm land, which sharpens the temperature gradient in Portugal, especially near the coast (Soares et al., 2014), thus hampering the growth of the mixing layer (which averages 200 m in height in summer and 400 m in winter; see Figure 4) and favoring the accumulation of D5 in the lower atmosphere. Therefore, the D5 sink by the OH reaction rate is too slow at this spatiotemporal scale to compete with vertical and horizontal transport. The use of a very high-resolution grid to capture the transport of pollutants is supported by the scientific literature for the target area (e.g. Lorente-Plazas et al., 2015a, 2015b; Nogueira et al., 2019).

In order to characterize the uncertainty that may be obtained when selecting the MC10 (135 mg capita^–1 day^–1) and the BU13 (310 mg capita^–1 day^–1) emission factors, Figure 2b depicts the concentration differences between BU13 minus MC10 for winter (left) and summer (right) periods. The main differences are found for urban areas, especially around Lisbon (Portugal) and Vigo (Galicia) in winter, with maximum differences of nearly 20 ng m^–3 in the 3-month average ground-level concentration, and in Porto in summer, where these differences can be over 15 ng m^–3. Averaged over the whole domain (excluding ocean areas) and for the entire year, the D5 concentrations using BU13 emissions are 2 times higher than those found when using the lowest MC10 emission factors (1.5 vs. 0.7 ng m^–3 as annual mean for the domain).

To build confidence in the above-mentioned results and to evaluate which emission factor provides the best results when simulating D5, the modeled concentrations of this pollutant were compared to those quantified in the SIPs. Model data were extracted from the grid cell at the sampling site location and averaged over the specific deployment time of the sampler. In evaluating the modeled atmospheric concentrations, a number of statistical parameters were considered. Spatial correlation coefficient (r), standard deviation (SD), and mean bias (MB) values are commonly used by the modeling community (e.g., Baldasano et al., 2011; Im et al., 2018; Monteiro et al., 2018; Kushta et al., 2019, among many others) and were therefore selected according to the criteria of Ratola and Jiménez-Guerrero (2016). Spatial-r is used to assess the ability of the model to capture the spatial pattern of D5 observations. The SD provides an idea of the dispersion of the data with respect to the mean (that is, the amplitude of the range of the modeled and observed data). Lastly, MB indicates the absolute difference between modeled and observed means. Other commonly used metrics used to express the results as a percentage are the normalized mean error (NME) and the normalized mean bias (NMB; Jiménez-Guerrero et al., 2008). These statistical parameters indicate the relative errors and bias, respectively. MB, NME, and NMB are defined in Equations 1–3, where Mod stands for model data, Obs for observations data, and N for number of observations.

Figure 5 (shaded) plots the field concentrations of D5, while circles in the same figure indicate the bias of the model concentrations using the MC10 (top row) and BU13 (bottom row) emission factors against observations. In general, simulations including the MC10 factor perform much better than in the case of BU13, where large overestimations in the modeled D5 concentrations are observed. Annual NMB is calculated to be –6.8% for MC10 and 73.6% for BU13, with similar differences calculated for winter and summer (Table 2). For instance, for a mean annual field concentration of 0.861 ng m^–3, the MC10 model mean is 0.714 ng m^–3 and the BU13 model mean is 1.495 ng m^–3. These results may be conditioned by the fact that BU13 uses a local emission factor for the city of Zürich, while the assumption of MC10 is based on a hemispheric global domain that includes our area of simulation.

In general, the lowest emission factor considered in this work leads to a very slight underestimation of the D5 concentrations, especially in summer (Table 2). For instance, the NMB for the simulations with MC10 emission factor is –5% in winter and –9% in summer. On the other hand, for simulations with emissions estimated from the BU13 emission factor, NMB and NME exceed +100% for summer and winter D5 concentrations. When focusing on the normalized error, for MC10 simulations, NME is higher in summer (45% vs. 20% in winter; and an annual mean of 43%), which may suggest a larger compensation in the sign of the bias (negative vs. positive errors) during the winter months. With respect to the skill of the model to reproduce the spatial pattern of D5 over the domain, correlation coefficients of .68, .69, and .75 are obtained under the MC10 emission factor for annual, winter, and summer periods, respectively. The values are considerably lower when using the BU13 approach (.58, .65, and .34), indicating that the MC10 factor is much more accurate in capturing the spatial features of D5 over the target domain, especially during summer. Lastly, the normalized bias of the SD gives an idea of the trend of the model to over- or underestimate the variability of D5 in the area of study. The use of BU13 leads to a large overestimation (NME of the SD over 80%). However, for the MC10 approach, the performance of the model is excellent, with NMB of the SD slightly overpredicted (25.0%) in winter and underpredicted (–15.7%) in summer.

Bearing in mind all the limitations mentioned previously, the NME observed are within the ranges commonly accepted (50% indicated in the European Directive 2008/50/EC) when using the MC10 emission factor. For individual stations, and since the evaluation results show higher errors for BU13, the discussion will focus only on the results obtained with MC10 (although the errors per sampling site and season can be seen for both factors in Figure 5).

As shown in Table 3, errors generally remain below 50% in both seasons, with the modeled concentrations very close to those found in SIPs. The skills of the modeling system agree with those for North America of Janechek et al. (2017) who found relative errors ranging from 27% to 56% when comparing WRF-CMAQ simulations against sampling sites in the Midwest of the United States. Highlighting the stations with the largest errors (Table 3), in summer an important overestimation was found in Porto (NME of 60%), possibly due to lower industrial activity as many companies choose to close in August and because the local inhabitants leave for their holidays elsewhere. Moreover, the effect of Cl atoms in the D5 loss in urban areas estimated by Alton and Browne (2020) may have played a part, since the Porto region is by far the largest conurbation among the sites of study. But it is in Midões, a remote area, where the model most overestimates D5 (by a factor of over 3: 0.11 ng m^–3 observed vs. 0.356 ng m^–3 modeled), which could possibly be explained by a lower accuracy between the modeled and the actual population. On the other hand, there are underestimations for the beach sites (Esmoriz with 21% and especially Praia da Luz, which has the highest field concentration of the study, with 1.78 ng m^–3 observed vs. 1.08 ng m^–3 modeled, a 39% error) and Z.I. Mota, an industrial site which includes a wastewater treatment plant that receives the input of nearby touristic coastal areas (24%). These errors could be higher if the model accounted for the Cl oxidation reaction. The overestimation at Midões persists in winter (over 189%), and there is now an overestimation in Praia da Luz (0.33 vs. 0.57 ng m^–3, 73% of error). However, this measured-to-model factor is in line with other works that use regional models to assess D5 atmospheric concentrations. For instance, Janechek et al. (2017) reported measured-to-model factors ranging from 1.37 to 2.75 in U.S. Midwest sites. Seasonal shifts in the populations and a possible change in the patterns of use of siloxane-containing products may lead to model uncertainties in simulated D5 concentrations (in particular, underestimations of D5 in summer). Indeed, Janechek et al. (2017) stated that VMS concentrations are heavily dependent on population (higher in urban and lower in rural areas), and it is expected that in countries like Portugal, the rise in temperature in the summer months prompts the locals and foreign tourists to sunbathe on the beaches or do other outdoor activities. According to the Portuguese Statistics Institute, the number of international tourists entering Portugal in 2019 was 24.6 million and the monthly distribution of hotel stays reached over 50% between June and September (INE, 2020). Consequently, a rise in the use of protective PCPs or cosmetics is expected, although no official data are available. D5 is typically the predominant VMS in the formulations of these products (Capela et al., 2016), in particular in the case of deodorants/antiperspirants, followed by cosmetics and sun care products (Wang et al., 2009; Dudzina et al., 2014; Brothers et al., 2017), which are those whose use is likely to increase in summer. Skin care and hair care products can also reach high concentrations (Capela et al., 2016), but their consumption is more equally distributed throughout the year. Future studies aimed at identifying sources of siloxanes, their respective emission patterns, and their correspondence with the concentrations in the environment will be extremely helpful in further refining (e.g., improved emission estimates, better process descriptions) this and other modeling approaches with regard to the representation of VMSs. Taking into account the good response of our approach to represent D₅, this can also be the base from which to extend it to other VMSs, namely D₄ and D₆, which are much less studied but also raising current concern regarding their use in some products.

For now, WRF+CHIMERE is yet another valid modeling tool for the representation of D5 atmospheric concentrations in other local and regional domains, since it estimates accurately not only marked seasonal variations but also the varied population densities commonly linked to different types and amounts of siloxane sources in the atmosphere. The good results can be compared to the use of the CMAQ model by Janechek et al. (2017) in the United States. Other modeling tools applied to D5 focused on different kinds of domains (typically smaller or larger) and we were able to represent certain flexibility, given the characteristics of our domain. This can be extremely important in helping not only scientists but also stakeholders and regulatory bodies to enforce the fulfillment of the VMS limitations in an EU framework, which can more easily be monitored directly (on the restricted products), or indirectly, if a reliable link can be made between the atmospheric concentrations and to the potential risk that they can represent by itself. Moreover, it is interesting that the per capita emissions that best explain the D5 observations in a European country are the same that best explain them in the United States (Janechek et al., 2017). D5 is a good marker for PCPs and has recently been used to quantify PCP emissions as a whole (e.g., Tang et al., 2015; Coggon et al., 2018; McDonald et al., 2018; Gkatzelis et al., 2021). Emission inventories in the United States suggest that PCP emissions can be as significant as mobile source emissions in major urban regions, likely contributing to urban ozone and the formation of secondary organic aerosols (McDonald et al., 2018; Qin et al., 2020). Consequently, the air quality impacts from PCP emissions could be similar in the United States and in major European cities (e.g., Lisbon). This association could also help in setting rules not only for D5 but also for other potentially harmful chemicals present in PCP formulations, as not only aquatic systems and biota (as noted in the ECHA assessments) but also regional air quality could be affected. Naturally, no improvement of the modeling efforts can be made without more intensive field campaigns, which can also benefit from models optimizing the distribution of key sampling points.

The inclusion of D5 in the full regional air quality modeling system WRF+CHIMERE represents an innovative approach to enhance the knowledge and reduce the common gaps in emerging pollutant studies. It can be said that despite the difficulties, this first modeling approach for siloxanes in Portugal yielded very promising results, since for most stations and for both seasons, the uncertainty in the predictions remains under 50%, which is a guide value reported by the European Directive 2008/50/EC for using modeling systems in the case of regulatory pollutants. The results obtained confirm that OH chemistry has a negligible effect on the simulated D5 ground-level concentrations at the scale of the investigated domain (Ferracci et al., 2018; Sakurai et al., 2019). The findings presented here indicate that emissions are one of the key entries in a model (Breivik and Alcock, 2002; Arnot et al., 2012; Breivik et al., 2012) together with other parameters, such as the reaction rates with oxidizing compounds or the physical-chemical properties of the chemicals. Since siloxanes are only now emerging as a global concern and the first limitations on their use (D4 and D5 in the EU) have just been established (February 2020), it is very difficult to identify emission rates and patterns at a local scale. For our target domain, the lowest emission factor found in the literature produces only slight overestimations of D5 concentrations when compared to field data, which reflect the accuracy and robustness of the WRF+CHIMERE approach.

Despite the limited role of the OH chemistry found in this contribution, the use of very high-resolution CTMs is needed in regional domains, in order to capture the transport of pollutants in areas with a complex orography (e.g. Lorente-Plazas et al., 2015a, 2015b; Nogueira et al., 2019). Nevertheless, the chemical mechanism included within the CTM can be simplified (if not neglected, considering D5 as a tracer in the atmosphere) because of the very limited role of the OH degradation reaction at local-to-regional scales. Therefore, further research for modeling ground-level D5 concentrations at these working scales should focus not only on chemical developments within the CTMs but also on other topics that lead to find accurate source contribution and better emission factors primarily for the VMSs of highest concern (D4, D5, and D6), which may have a strong local fingerprint. The assumption of similar geographic and temporal emission patterns and constant population throughout the year seems reasonable on a continental or global scale but may not account for local variations in every domain. For instance, it is expected that the influx of tourists in the summer can cause much more impact and urban stress in the south of Europe than in the northernmost latitudes. To remove this limitation, the inclusion of climate change scenarios is already planned for follow-up studies with a more comprehensive scope.

The modeling data are available on the Zenodo repository (DOI: http://dx.doi.org/10.5281/zenodo.4553672) and the data from the field campaign are presented in Supplemental Material (Table S1). For further information, please contact the corresponding author: Email: [email protected].

Information is provided on values for terrain height, land use index, and vegetation fraction on the target domain, the concentrations obtained in the field sampling campaigns, and D5 emissions. This material is available free of charge via the Internet.

The supplemental files for this article can be found as follows:

Figure S1. Values for terrain height, land use index, and vegetation fraction (Loveland et al., 2000; Skamarock et al., 2008) included in WRF+CHIMERE simulations for the domain. The United States Geological Survey 24-categories are 1: Urban and Built-Up Land; 2: Dryland Cropland and Pasture; 3: Irrigated Cropland and Pasture; 4: Mixed Dryland/Irrigated Cropland and Pasture; 5: Cropland/Grassland Mosaic; 6: Cropland/Woodland Mosaic; 7: Grassland; 8: Shrubland; 9: Mixed Shrubland/Grassland; 10: Savanna; 11: Deciduous Broadleaf Forest; 12: Deciduous Needleleaf Forest; 13: Evergreen Broadleaf Forest; 14: Evergreen Needleleaf Forest; 15: Mixed Forest; 16: Water Bodies; 17: Herbaceous Wetland; 18: Wooded Wetland; 19: Barren or Sparsely Vegetated; 20: Herbaceous Tundra; 21: Wooded Tundra; 22: Mixed Tundra; 23: Bare Ground Tundra; 24: Snow or Ice.

Figure S2. Simulations for the mean concentration of hydroxyl radical (cm^–3) averaged for winter (left) and summer (right) in 2014.

Table S1. D5 concentrations obtained from the passive air sampling field campaigns.

Table S2. Appraisal of D5 emissions according to its production and uses (Source: Brooke et al., 2009).

The Regional Atmospheric Modelling Group (G-MAR) at the University of Murcia and the Laboratory for Process Engineering, Environment, Biotechnology and Energy (LEPABE) members at the University of Porto are thanked for their contributions to the fruitful scientific discussions.

This work was financially supported by Base Funding—UIDB/00511/2020 of the Laboratory for Process Engineering, Environment, Biotechnology and Energy—LEPABE, financed by national funds through FCT/MCTES (PIDDAC); Project Reference EXPL/AAG-MAA/0981/2013. PJ-G acknowledges projects REPAIR-CGL2014-59677-R and ACEX-CGL2017-87921-R, funded by the Spanish Ministry of the Economy and Competitiveness/Agencia Estatal de Investigación and the European Regional Development Fund (ERDF/FEDER). NR thanks the Portuguese Foundation for Science and Technology (FCT) for the financial support of his work contract through the Scientific Employment Stimulus—Institutional Call—(CEECINST/00049/2018) and the “Programa Jiménez de la Espada” (Ref. 19641/IV/14) from Fundación Séneca—Science and Technology Agency in the Region of Murcia.

The authors declare no competing interests.

Contributed to conception and design: PJ-G, NR.

Contributed to acquisition of data: PJ-G, NR.

Contributed to analysis and interpretation of data: PJ-G, NR.

Drafted and/or revised the article: PJ-G, NR.

Approved the submitted version for publication: PJ-G, NR.

How to cite this article: Jiménez-Guerrero, P, Ratola, N. 2021. Characterizing the regional concentrations and seasonality of the emerging pollutant decamethylcyclopentasiloxane (D5) using a WRF+CHIMERE modeling approach. Elementa: Science of the Anthropocene 9(1). DOI: https://doi.org/10.1525/elementa.2020.00137

Domain Editor-in-Chief: Detlev Helmig, Boulder AIR LLC, Boulder, CO, USA

Associate Editor: Armin Wisthaler, Atmospheric Chemistry, University of Oslo, Oslo, Norway

Knowledge Domain: Atmospheric Science