Developing Ozone Risk Assessment for Larch Species

Ozone (O3) risk assessment for the protection of forests requires species-specific critical levels (CLs), based on either O3 concentrations (AOT40) or stomatal uptake (PODY) accumulation over the growing season. Larch (Larix sp.) is a genus with O3-susceptible species, widely distributed in the northern hemisphere and with global economic importance. We analyzed published and unpublished data of Japanese larch (Larix kaempferi) and its hybrid F1 (Larix gmelinii var. japonica × L. kaempferi) stomatal responses for developing a parameterization of stomatal conductance model and estimating PODY-based CLs with two Y thresholds, that is, 0 and 1 nmol m–2 s–1 projected leaf area (PLA). In parallel, we estimated AOT40-based CLs. The results show that the AOT40-based CLs for a 2% and 4% biomass loss in Japanese larch were 5.79 and 11.59 ppm h, that is, higher than those in hybrid larch F1 (2.18 and 4.36 ppm h AOT40), suggesting a higher O3 susceptibility of the hybrid. However, the use of PODY reconciled the species-specific differences, because the CLs were similar, that is, 9.40 and 12.00 mmol m–2 POD0 and 2.21 and 4.31 mmol m–2 POD1 in Japanese larch versus 10.44 and 12.38 mmol m–2 POD0 and 2.45 and 4.19 mmol m–2 POD1 in the hybrid, for 2% and 4% biomass loss, respectively. Overall, the CLs were lower than those in other forest species, which suggests a relatively high susceptibility of these larches. These results will inform environmental policy-makers and modelers about larch susceptibility to O3.


INTRODUCTION
Tropospheric ozone (O 3 ) is the most widespread phytotoxic air pollutant (Mills et al., 2018). In the period 1995-2014, control measures were effective in North America and Europe, as indicated by a decrease of O 3 concentrations, while a significant increase in O 3 concentrations occurred in East Asia (Chang et al., 2017;Mills et al., 2018). Ozone has a strong oxidative capacity and may cause severe injury to forests (Paoletti, 2007;Li et al., 2017). To assess O 3 risk to forests, different metrics have been developed (Lefohn et al., 2018). One of the most common metrics is AOT40, that is, the accumulated exposure over an hourly threshold of 40 ppb during the growing season, although there is a general consensus that the accumulated stomatal O 3 flux -or phytotoxic ozone dose (POD) -is more biologically meaningful as it estimates the amount of O 3 actually entering the plants through the stomata (Paoletti and Manning, 2007). A flux threshold Y below which O 3 uptake is not expected to be injurious to plants has been postulated. For all tree species, a uniform threshold of Y = 1 nmol m −2 s −1 projected leaf area (PLA) was recommended by the Convention on Long-Range Transboundary Air Pollution (CLRTAP, 2017) based on Büker et al. (2015). For easier calculation, a Y threshold of 0 nmol m −2 s −1 PLA was also recommended, if we assume that all O 3 molecules induce a physiological reaction after uptake (De Marco et al., 2015Anav et al., 2016), which is a plausible assumption in the light of low-dose adaptive responses (Agathokleous et al., 2019).
For the protection of susceptible vegetation from O 3 , critical levels (CLs) are recommended, defined as the "concentration, cumulative exposure or cumulative stomatal flux of atmospheric pollutants above which direct adverse effects on susceptible vegetation may occur according to present knowledge" (CLRTAP, 2017). CLs are derived for either a 2% (Norway spruce) reduction or a 4% (beech/birch, Mediterranean deciduous and evergreen species) reduction in annual new growth (based on aboveground, root, or whole-tree biomass) of young trees up to 10 years old. AOT40-based CLs for tree biomass loss (5%) are available for Fagus sylvatica and Betula pendula in a previous version of the ICP Vegetation manual (CLRTAP, 2014; AOT40-based CLs are not included in the latest version) and for some other species in the literature (e.g., 18 Japanese species including two larch species, Yamaguchi et al., 2011; Populus deltoides cv. "55/56" × P. deltoides cv. "Imperial" and Populus euramericana cv. "74/76, " Shang et al., 2017). Stomatal flux-based CLs are available for F. sylvatica, B. pendula, Picea abies, Quercus faginea, Quercus pyrenaica, Quercus robur, Quercus ilex, Ceratonia siliqua, and Pinus halepensis in the ICP Vegetation manual (CLRTAP, 2017) and for few other species in the literature (Zelkova serrata, Hoshika et al., 2012;Quercus pubescens, Hoshika et al., 2018b;Pinus pinea, Hoshika et al., 2017;hybrid poplars, Zhang et al., 2018;Feng et al., 2019b; Fagus crenata, Quercus serrata, Quercus mongolica var. crispula, and Betula platyphylla var. japonica, Yamaguchi et al., 2019). For estimating PODY (phytotoxic ozone dose above a threshold Y nmol m −2 s −1 )-based CLs, a speciesspecific parameterization of the stomatal flux or DO 3 SE model is required (Emberson et al., 2000;Büker et al., 2012). There is a need of more species-specific CLs for biomass loss in forest species, especially for forest species in Asia, where elevated O 3 pollution levels are a serious risk for forests at present Mills et al., 2018;Feng et al., 2019a).
Larch (Larix sp.) is a widely distributed genus (Pinaceae family) with global economic importance, which includes some of the few deciduous conifer species. Larch is among the dominant tree species of northern hemisphere boreal forests. Hence, its natural distribution range is very wide and spans from Siberia to Canada, passing through Europe, mountainous China, and Japan. Larches provide high-quality wood and are commercially valuable (Bardak et al., 2019). As any pioneer species, larches have a relatively high growth rate and stomatal conductance (Streit et al., 2014;Agathokleous et al., 2017;Hoshika et al., 2018c). Although their susceptibility to O 3 has been investigated in several papers (Wieser and Havranek, 1996;Matsumura, 2001;Watanabe et al., 2006;Koike et al., 2012;Agathokleous et al., 2017;Sugai et al., 2018Sugai et al., , 2019, a comprehensive risk assessment including parameterization of the stomatal conductance model and definition of CLs for biomass losses is missing. Previous studies focused on the biomass responses to O 3 of Japanese larch (Larix kaempferi) and its hybrid F 1 (Larix gmelinii var. japonica × L. kaempferi). Hybrid F 1 displays heterosis and is important for timber production and afforestation due to more desirable characteristics compared to its parents, with a significant superiority in terms of growth rates (Ryu et al., 2009;Kita et al., 2009;Agathokleous et al., 2017;Sugai et al., 2018). A question arises whether hybrid clones, selected for fast-growing capacities, are representative of natural forest responses to O 3 when used in manipulative experiments (e.g., Di Baccio et al., 2008;Hu et al., 2015;Dusart et al., 2019;Podda et al., 2019).
Our aim was to collate published and unpublished data from previous experiments for developing a parameterization of the DO 3 SE model for Japanese larch and its hybrid F 1 and estimating the CLs not to be exceeded for the protection of these larch species from O 3 . Based on published research documenting a higher O 3 susceptibility of the faster-growing hybrid F 1 than the slowergrowing Japanese larch (Agathokleous et al., 2017;Sugai et al., 2018), we hypothesized that the CLs of hybrid F 1 have a lower susceptibility than that of the wild Japanese larch.

MATERIALS AND METHODS
A literature survey was conducted in Web of Science (9 December 2019), with the keywords "ozone" and "larch" or "larix" (search method: Topic). All the identified papers (n = 33 and 36 for each combination; most were duplicates) were reviewed for relevance, including whether they reported O 3 and biomass data. Finally, data on O 3 concentrations, exposure duration, and total biomass were collected from six published experiments carried out in open-top chambers (OTCs) ( Table 1: Matsumura, 2001;Watanabe et al., 2006;Koike et al., 2012;Wang et al., 2015;Sugai et al., 2018Sugai et al., , 2019) and used to calculate AOT40 and percentage losses of biomass relative to controls in low-O 3 air. Data from combined experiments, such as O 3 with either fertilization or CO 2 , were not included. Data of Dahurian larch (L. gmelinii var. japonica) from the same experiments were not included because of scarcity, thus being insufficient for analysis.
Individual measurements of stomatal conductance across a range of environmental conditions were obtained from the authors Sugai et al. (2018Sugai et al. ( , 2019 and Agathokleous (unpublished). Measurements by Agathokleous (unpublished) were carried out in field-grown 2-year-old larch seedlings at the Sapporo experimental forest, Hokkaido University, in Japan ( Table 1). All measurements were carried out by means of Li-Cor 6400 gas analyzers (Li-Cor Inc., Lincoln, NE, United States). As soil water content measurements were missing, we used the following simplified formula for the estimation of the stomatal conductance g sto in the DO 3 SE model (CLRTAP, 2017): where g max is the maximum stomatal conductance of either Japanese larch or its hybrid F 1 , f min is the species-specific minimum stomatal conductance, and f light , f temp , and f VPD account for the effects of photosynthetic photon flux density (PPFD), air temperature (T), and vapor pressure deficit (VPD), respectively, on stomata. Parameterization was carried out using a boundary line analysis (Alonso et al., 2008;Braun et al., 2010;Hoshika et al., 2012). First, the g sto data were divided into classes with the following stepwise increases for each variable: 200 µmol photons m −2 s −1 for PPFD (when the values were less than 200 µmol photons m −2 s −1 , PPFD classes at 50 µmol photons m −2 s −1 steps were adopted), 2 • C for T, and 0.2 kPa for VPD. A function was fitted against each model variable based on 95th percentile values per class of environmental factors. Values of g max and f min were calculated as the 95th percentile and 5th percentile, respectively Bičárová et al., 2019). For details of f light , f temp , and f VPD , see CLRTAP (2017). Stomatal O 3 uptake (F st ; nmol m −2 s −1 ) was calculated as follows: where r c is the leaf surface resistance [= 1/(g sto + g ext ); s m −1 ] and g ext is the external leaf or cuticular conductance (= 0.0004 m s −1 , CLRTAP, 2017). The standard DO 3 SE model considers the leaf boundary layer resistance (r b ): where the factor 1.3 accounts for the difference in diffusivity between heat and O 3 , 150 is the empirical constant, L d is the cross-wind leaf dimension (0.008 m for conifers, CLRTAP, 2017), and u is the wind speed. The wind speed data were not available in collected literatures. However, in OTCs, since a constant ventilation from the blowers is realized, r b is less important compared with stomatal resistance (r sto ) (Unsworth et al., 1984;Uddling et al., 2004;Tuovinen et al., 2009). This is supported by the fact that the r b /r c ratio was small in the present study when assuming that r sto was r sto_min (= 1/g max ) and wind speed was constant inside a chamber (r b /r c = 0.07 and 0.06 at 1 m s −1 and 0.05 and 0.04 at 2 m s −1 of wind speed in hybrid and Japanese larch, respectively). Here, we assumed that r b was negligible for the calculation of F st . PODY (mmol m −2 ) was estimated from hourly data as follows: where Y is a species-specific threshold of stomatal O 3 uptake (nmol m −2 s −1 ) and t = 1 h is the averaging period. F st,i is the ith hourly stomatal O 3 uptake (nmol m −2 s −1 ), and n is the number of hours included in the calculation period. Y is subtracted from each F st , i when F st , i > Y. PODY was then estimated based on hourly data of air temperature, solar photosynthetic active radiation, and VPD as registered locally and accumulated over the duration of the experiments from the six papers (Table 1). Data from Matsumura (2001) and Watanabe et al. (2006) were excluded from this analysis because of missing meteorological data.
To establish PODY-based dose-response relationships, two representative values of Y (= 0 or 1 nmol m −2 s −1 ) were tested. This is because CLRTAP (2017) suggested POD1 to be suitable for biomass assessment in elevated O 3 while several studies reported a better performance of POD0 rather than POD1 for O 3 risk assessment (e.g., Sicard et al., 2016). CLs were estimated for a total biomass reduction of both 2% as suggested for deciduous species and 4% as suggested for non-Mediterranean conifer species (CLRTAP, 2017). In addition, since CLRTAP (2017) provided an AOT40-based CL corresponding to a 5% biomass reduction for forests, the CLs for the 5% biomass reduction were also shown. For PODY, CLs were calculated, referring to a "REF10" PODY calculated at a constant O 3 concentration of 10 ppb referring to a "pre-industrial" O 3 concentration, as recommended by CLRTAP (2017).
Simple linear regression analyses were used to assess the relationships between O 3 indices (AOT40, POD0, and POD1) and relative biomass. In addition, to compare the g max values between the two larches, Student's t-test was performed on values within the top five percentile in g sto data. Results were considered significant at p < 0.05. All the analyses were performed using R 3.5.1 (R Core Team, 2018).

RESULTS
The parameterization of the stomatal conductance model (Figure 1) resulted in very similar values for Japanese larch and its TABLE 2 | DO 3 SE model parameters for Japanese larch and hybrid F 1 , where g max is maximum stomatal conductance; f min is minimum stomatal conductance; f light_a is a parameter determining the shape of the hyperbolic relationship of stomatal response to light; T max , T opt , and T min are the maximum, optimal, and minimum temperatures, respectively, for calculating the function f temp that expresses the variation of g sto with temperature; VPD min and VPD max are the vapor pressure deficit for attaining minimum and maximum stomatal aperture, respectively (f VPD ).

Parameter
Japanese hybrid F 1 ( Table 2). The g max in hybrid larch was slightly higher than that in Japanese larch although g max values in the two larches were not statistically different (p = 0.48, Student's t-test for the values within the top five percentile in g sto , data not shown). On the other hand, f min was slightly higher in Japanese larch than in hybrid larch F 1 . All the dose-response relationships were significant. When AOT40 was applied, in particular, a higher slope was found for hybrid larch F 1 than for Japanese larch (Figure 2).
The CLs calculated on the basis of these dose-response relationships were 2.7 times higher in Japanese larch than in its hybrid F 1 when AOT40 was used, while PODY-based CLs were similar between the two species when using either no Y threshold or a Y threshold of 1 nmol m −2 s −1 PLA (Table 3).

DISCUSSION
The boreal area in the northern hemisphere where larches are widely distributed is at risk of changes due to the FIGURE 2 | Dose-response relationships for total biomass losses in two species of larch seedlings on the basis of AOT40, POD0, or POD1 in different experiments. Black lines denote the regressions, and gray lines denote the 95% confidence intervals. Asterisks show the level of significance: ***p ≤ 0.001; *p ≤ 0.05.
Frontiers in Forests and Global Change | www.frontiersin.org TABLE 3 | Critical levels for larch protection from ozone corresponding to a total biomass loss of 2%, 4%, or 5% and based on the dose-response relationships in Figure 2. potential O 3 impact on photosynthetic carbon assimilation , as estimated by several global atmospheric chemistry transport models and representative concentration pathways emission scenarios. For a realistic estimate of O 3 risks to forests, CLs should be developed for the major forest species or types. Even though natural areas and plantations for larch trees are very wide and larch is a major genus of the forest category defined as boreal deciduous species, PODY-based CLs were not yet available for larch and are suggested here for the first time.
Organismic "sensitivity" may be defined as "the response of an organism (i.e., biological deviation) above or below a homeostatic state (control) of a set of biological traits, after sensing some environmental stress-inducing agents" (Agathokleous and Saitanis, 2020). However, "the organismal predisposition to be inhibited or adversely affected by or die of a xenobiotic, " as expressed by "negative (inhibitory or adverse) effects induced by diseases or environmental challenges, " is termed susceptibility (Agathokleous and Saitanis, 2020). Hence, organismic susceptibility can be assessed by studying dose/exposureresponse relationships and, in particular, by comparing CLs among organisms (Agathokleous and Saitanis, 2020). Since the CLs are affected by the O 3 metric used to develop dose/exposure-response relationships, susceptibility rankings can be different depending on the O 3 metric used (Agathokleous et al., 2019).
So far, CLs have been estimated for a total biomass reduction in either deciduous broadleaf and Mediterranean conifer species (recommended biomass loss: 2%) or non-Mediterranean evergreen conifer species (recommended biomass loss: 4%) (CLRTAP, 2017). As larch is both a deciduous species and a non-Mediterranean conifer species, we decided to calculate the CLs for both the loss thresholds of 2% and 4%. We decided also to calculate the CLs for AOT40, although this metric is known for not being able to assess how much O 3 enters the leaf through the stomata (Paoletti and Manning, 2007). However, it is still the legislative standard in Europe (Directive 2008/50), is used in many other continents (e.g., Agathokleous et al., 2018;Pleijel et al., 2019) because it is simple to calculate, and helps in the comparison with other results in the literature. The AOT40-based CL suggested so far for O 3 -susceptible deciduous broadleaves (F. sylvatica and B. pendula, 5 ppm h for a 5% biomass loss; CLRTAP, 2014CLRTAP, , 2017 is similar to that of hybrid larch F 1 (5.45 ppm h for 5% biomass loss), while Japanese larch showed a markedly higher AOT40-based CL corresponding to 5% loss (i.e., 14.48 ppm h). Based on a reanalysis of only two of the papers investigated here (Matsumura, 2001;Watanabe et al., 2006), Yamaguchi et al. (2011) had already suggested high O 3 susceptibility for Japanese larch. In fact, the AOT40based CLs that they recommended were consistent with those found in our work (i.e., 8-15 ppm h). In addition, our results would suggest a higher susceptibility to O 3 of the hybrid and confirm previous studies where ecophysiological responses of the hybrid were more severely affected by O 3 exposure than those of Japanese larch (Koike et al., 2012;Sugai et al., 2019).
An accurate parameterization of stomatal conductance model is essential for the flux-based O 3 risk assessments (Emberson et al., 2000). For larch, the information of leaf-level g sto parameters was limited, although some studies tried to estimate O 3 uptake at stand level by sap-flow measurements (Nunn et al., 2007) and at forest level by eddy covariance (Finco et al., 2017). Wieser and Havranek (1995) previously reported just stomatal VPD responses to estimate stomatal O 3 uptake in European larch (Larix decidua). Our study is the first one to achieve a proper leaf-level parameterization (g max , f min , f light , f temp , and f VPD ) in larch trees to develop a fluxbased approach. The maximum value of g sto in European larch by Wieser and Havranek (1995) was 150 mmol O 3 m −2 PLA s −1 , which was comparable to the g max values in our findings. Interestingly, hybrid larch F 1 showed a slightly higher g max (140 vs. 120 mmol O 3 m −2 PLA s −1 in Japanese larch). As g max is known to play the most important role in determining PODY (Tuovinen et al., 2007), the small difference in g max between the two species translated into a higher stomatal uptake of O 3 by the hybrid at similar AOT40 levels; that is, the higher susceptibility of the hybrid under similar O 3 exposures was due to a higher stomatal uptake. It is well known that fast-growing species with high stomatal conductance are susceptible to O 3 because of an elevated stomatal uptake (Feng et al., 2018;Hoshika et al., 2018a). When the CLs are calculated on a PODY basis, in fact, the two species showed surprisingly similar CLs: 9.40 and 12.00 mmol m −2 POD0 and 2.21 and 4.31 mmol m −2 POD1 in Japanese larch versus 10.44 and 12.38 mmol m −2 POD0 and 2.45 and 4.19 mmol m −2 POD1 in the hybrid, for 2% and 4% biomass loss, respectively. These POD1-based values are below the CL recommended for non-Mediterranean trees (5.7 mmol m −2 ; CLRTAP, 2017), suggesting that these larches are more susceptible to O 3 even when evaluated on the basis of stomatal flux. Different susceptibilities to O 3 injury in the two larch species may be also due to different antioxidant capacities (Di Baccio et al., 2008). Although monoterpene emissions from leaves were preliminarily studied (Mochizuki et al., 2017), the role of antioxidants, secondary metabolites, and other leaf defensive molecules in the response of these two species to O 3 remains elusive.

CONCLUSION
Based on a reanalysis of literature results and new measurements, we conclude that Japanese larch and its hybrid F 1 should be classified as species with considerable O 3 susceptibility as compared to the CLs available so far for other forest species. We also found that AOT40 and PODY can give very different results when assessing a species' susceptibility to O 3 . While AOT40 suggested a higher susceptibility of hybrid F 1 , PODY did not highlight marked differences between the two species. Future research should clarify the O 3 susceptibility of hybrid clones versus their wild forest species and increase the number of forest species with a species-specific parameterization and PODY-based CLs, especially in the Asian continent. This kind of information is needed for improving our modeling capacities, assessing O 3 risks to local-to-global forests, and transferring this knowledge to environmental policy-makers.

DATA AVAILABILITY STATEMENT
Basic raw data are available with YH (Italy) or TK (Japan).

AUTHOR CONTRIBUTIONS
EP conceptualized the work and wrote the manuscript. EA, TS, and TK provided the data. YH analyzed the data. All authors reviewed the manuscript.

ACKNOWLEDGMENTS
This study is partly supported by the LIFE15 ENV/IT/000183 project MOTTLES and JST-2019 (Grant No. JPMJSC18HB).