Impacts of Surface Ozone Pollution on Global Crop Yields: Comparing Different Ozone Exposure Metrics and Incorporating Co-effects of CO2

Surface ozone (O3) pollution poses significant threats to crop production and food security worldwide, but an assessment of present-day and future crop yield losses due to exposure to O3 still abides with great uncertainties, mostly due: (1) to the large spatiotemporal variability and uncertain future projections of O3 concentration itself; (2) different methodological approaches to quantify O3 exposure and impacts; (3) difficulty in accounting for co-varying factors such as CO2 concentration and climatic conditions. In this paper, we explore these issues using a common framework: a consistent set of simulated present-day O3 fields from one chemical transport model, coupled with a terrestrial ecosystem-crop model to derive various O3 exposure metrics and impacts on relative crop yields worldwide, and examine the potential effects of elevated CO2 on O3-induced crop yield losses. Throughout, we review and explain the differences in formulation and parameterization in the various approaches, including the concentration-based metrics, flux-based metrics, and mechanistic biophysical crop modeling. We find that while the spatial pattern of yield losses for a given crop is generally consistent across metrics, the magnitudes can differ substantially. Pooling the concentration-based and flux-based metrics together, we estimate the present-day globally aggregated yield losses to be: 3.6 ± 1.1% for maize, 2.6 ± 0.8% for rice, 6.7 ± 4.1% for soybean, and 7.2 ± 7.3% for wheat; these estimates are generally consistent with previous studies but on the lower end of the uncertainty range covered. We attribute the large combined uncertainty mostly to the differences among methodological approaches, and secondarily to differences in O3 and meteorological inputs. Based on a biophysical crop model that mechanistically simulates photosynthetic and yield responses of crops to stomatal O3 uptake, we further estimate that increasing CO2 concentration from 390 to 600 ppm reduces the globally aggregated O3-induced yield loss by 21–52% for maize and by 27–38% for soybean, reflecting a CO2-induced reduction in stomatal conductance that in turn alleviates stomatal O3 uptake and thus crop damage. Rising CO2 may therefore render the currently used exposure-yield relationships less applicable in a future atmosphere, and we suggest approaches to address such issues.

Surface ozone (O 3 ) pollution poses significant threats to crop production and food security worldwide, but an assessment of present-day and future crop yield losses due to exposure to O 3 still abides with great uncertainties, mostly due: (1) to the large spatiotemporal variability and uncertain future projections of O 3 concentration itself; (2) different methodological approaches to quantify O 3 exposure and impacts; (3) difficulty in accounting for co-varying factors such as CO 2 concentration and climatic conditions. In this paper, we explore these issues using a common framework: a consistent set of simulated present-day O 3 fields from one chemical transport model, coupled with a terrestrial ecosystem-crop model to derive various O 3 exposure metrics and impacts on relative crop yields worldwide, and examine the potential effects of elevated CO 2 on O 3 -induced crop yield losses. Throughout, we review and explain the differences in formulation and parameterization in the various approaches, including the concentration-based metrics, flux-based metrics, and mechanistic biophysical crop modeling. We find that while the spatial pattern of yield losses for a given crop is generally consistent across metrics, the magnitudes can differ substantially. Pooling the concentration-based and flux-based metrics together, we estimate the present-day globally aggregated yield losses to be: 3.6 ± 1.1% for maize, 2.6 ± 0.8% for rice, 6.7 ± 4.1% for soybean, and 7.2 ± 7.3% for wheat; these estimates are generally consistent with previous studies but on the lower end of the uncertainty range covered. We attribute the large combined uncertainty mostly to the differences among methodological approaches, and secondarily to differences in O 3 and meteorological inputs. Based on a biophysical crop model that mechanistically simulates photosynthetic and yield responses of crops to stomatal O 3 uptake, we further estimate that increasing CO 2 concentration from 390 to 600 ppm reduces the globally aggregated O 3 -induced yield loss by 21-52% for maize and by 27-38% for soybean, reflecting a CO 2 -induced reduction in stomatal conductance that in turn alleviates stomatal O 3

INTRODUCTION
Future food production is projected to increase by at least 50% by year 2050 in order to fulfill the rapidly growing global food demand, which is a combined result of population growth and a worldwide shift toward a more meat-intensive diet from staple food (Alexandratos and Bruinsma, 2012). However, food production is threatened by various environmental problems, including air pollution. Tropospheric ozone (O 3 ), which is one of the major air pollutants worldwide, has adverse impacts on natural vegetation and crops (e.g., Betzelberger et al., 2012;Feng et al., 2015Feng et al., , 2019Mills et al., 2018a). It is projected that O 3 pollution will become more severe in the future due to a warmer climate and, depending on the region, more anthropogenic emissions of O 3 precursors including nitrogen oxides (NO x = NO 2 + NO), methane (CH 4 ) and carbon monoxide (CO). For instance, under the Intergovernmental Panel on Climate Change (IPCC) Representative Concentration Pathway 8.5 (RCP8.5), the global tropospheric O 3 burden is estimated to increase from 337 Tg in year 2000 to 395 Tg in year 2100 (Young et al., 2013).
To assess the impacts of O 3 pollution on crop yields, several metrics have been developed based on results from both field and open-top chamber experiments. The first batch of these metrics are various concentration-based metrics, including daylighthour growing-season average O 3 concentration (M7 or M12), accumulated daylight O 3 concentration above certain thresholds (AOT40 and SUM06), and continuously weighted growingseason average O 3 concentration (Adams et al., 1989;Lesser et al., 1990;Lee and Hogsett, 1996;Fuhrer et al., 1997;Wang and Mauzerall, 2004;Mills et al., 2007). AOT40 and SUM06 in particular consider the baseline ability of plants to detoxify low levels of O 3 , whereby crops suffer yield losses only if the hourly O 3 concentration exceeds 40 (AOT40) or 60 (SUM06) ppb. W126 considers the non-linear response of yield loss to O 3 exposure whereby higher concentrations will result in progressively more severe yield losses. These metrics have been used extensively to evaluate global yield losses due to present-day O 3 concentrations, which are estimated to be 2-16% for the four major staple crops (wheat, rice, maize, soybean) (Ainsworth, 2017). These metrics are easy to calculate, but they typically do not consider the combined effects of O 3 with other factors that regulate plant O 3 uptake, and overlook the active responses of plant ecophysiology to ambient climatic and environmental changes and thus likely inadequate for examining yield losses in a future climate and atmospheric environment (Musselman et al., 2006). The covarying factors that may modulate O 3 -crop relationships include ambient CO 2 concentration, temperature, drought conditions, soil water stress, etc., all of which may undergo tremendous changes under future climate change conditions (Emberson et al., 2000;Tai et al., 2014). O 3 -induced plant damage is mainly caused by the stomatal uptake of O 3 into the leaf interior instead of direct plant surface deposition (Amiro et al., 1984;Nouchi, 2002). Stomatal uptake of O 3 can be strongly influenced by the extent of stomatal opening, which is usually quantified as stomatal conductance (g s ) and strongly dependent on environmental conditions such as soil water availability (e.g., Tingey and Hogsett, 1985), light intensity (e.g., Paoletti and Grulke, 2010), and CO 2 concentration (e.g., Onandia et al., 2011). Atmospheric CO 2 concentration has been rising, and plants generally respond by reducing stomatal conductance, thereby reducing water loss from transpiration and improving water use efficiency (WUE) (e.g., Franks et al., 2013). The closure of the stomata indirectly alleviates the uptake of O 3 and thus the damage inflicted on them by O 3 . Some studies showed that O 3 -induced yield loss is alleviated under elevated CO 2 concentration (Fiscus et al., 1997;Tao et al., 2017). Meanwhile, higher CO 2 itself enhance photosynthesis rates and crop yields via the fertilization effect (e.g., Rosenzweig et al., 2014), which may further complicates crop-O 3 -CO 2 relationships. Likewise, temperature extremes and the often associated drought conditions not only damage crop yields but also modulate crop physiology and thus their uptake of O 3 . Moreover, temperature and O 3 concentration are often highly positively correlated with each other (Jacob and Winner, 2009), which can confound any observed yield-temperature and yield-O 3 relationships as found in observational and experimental studies (Tai et al., 2014;Tai and Martin, 2017). For all these reasons, the concentration-based exposure-yield loss metrics described above may not hold in the future as they do not account for how crop physiological responses to the changing atmospheric environment may interfere with stomatal O 3 uptake and damage, or how O 3 may vary with other environmental conditions that also affect crops.
In the recent two decades a newer approach has been developed and popularized to assess the statistical relationships between stomatal O 3 uptake and crop yields (Emberson et al., 2000;Pleijel et al., 2007). In this, so called, flux-based approach, the stomatal control of O 3 uptake and its environmental dependence (especially in relation to soil water availability) are explicitly taken into account. This approach usually does not require mechanistic simulation of ecophysiological processes (e.g., photosynthesis, carbon allocation) and responses of crops. Typically, stomatal conductance is formulated as a semiempirical function of a baseline, standard-conditioned stomatal conductance modified multiplicatively by a series of "activity factors" accounting for variations in environmental conditions; such kind of a multiplicative function is usually termed "Jarvis-type" formulation due to its origin in the study of Jarvis (1976). Generally, the calculated instantaneous O 3 flux, above a certain threshold (e.g., Y nmol O 3 m −2 s −1 ) to account for natural detoxification, is then cumulated over the growing season to calculate a phytotoxic O 3 dose (POD Y , mmol O 3 m −2 ); the statistical relationships between POD Y and crop yields have been determined by several experiments thus far, and are used to determine O 3 impacts on crop yields. Such flux-based metrics have been shown to perform better in statistical terms than concentration-based metrics, and better predict the distribution of O 3 damage when combined with O 3 fields simulated by regional atmospheric chemistry models (Pleijel et al., 2007;Millls et al., 2011). An important reason for their better performance is that they explicitly account for stomatal control of O 3 uptake, the variation of which can largely influence crop damage even with the same ambient O 3 concentrations (Clifton et al., 2020).
Finally, O 3 impacts on crops can also be estimated using mechanistic ("process-based") models that explicitly consider the effects of O 3 flux (either instantaneous or cumulative) on photosynthesis rate and/or stomatal conductance, which ultimately affects crop yields through modifying carbon assimilation and other crop biogeochemical processes. Such an approach has been applied in tandem with land surface and terrestrial ecosystem models that include a crop component, e.g., the Community Land Model (CLM) (Oleson et al., 2013), Joint-UK Land Environment Simulator (JULES) (Osborne et al., 2015), and MOSES-TRIFFID (Sitch et al., 2007). Compared to the metric-based methods, using a mechanistic crop model can more fully simulate different agricultural practices (e.g., irrigation) and crop physiology fully coupled with photosynthesis, which responses to varying environmental conditions including longterm changes in climate and atmospheric chemical composition. However, applications of these process-based models to a large or global scale are challenging due to limited O 3 flux data for crops to more comprehensively parameterize the essential model components; most large-scale O 3 -crop impact studies to date have still relied on concentration-based or flux-based metrics.
The above three approaches to estimate O 3 impacts on crops all have different input requirements and have traditionally been used with different modeling frameworks, which makes cross-comparison among them challenging. In this paper, we describe and compare present-day O 3 -induced yield losses using these different methods including concentration-based metrics (M7/M12, AOT40, W126), flux-based metrics (POD Y ) and a mechanistic crop model, essentially keeping input requirements and modeling frameworks as consistent as possible so as to reduce inter-model discrepancies that are not due to fundamental methodological differences (e.g., input ozone fields and meteorology).

Model Description
We use the Terrestrial Ecosystem Model In R (TEMIR) to compute biogeophysical responses of terrestrial ecosystems, including crops, to changes in the atmospheric and terrestrial environment. Equipped with the same driving meteorology as the GEOS-Chem chemical transport model (CTM), TEMIR is designed to be highly compatible with GEOS-Chem so that these two models can be coupled asynchronously to examine atmospheric chemistry-biosphere interactions. TEMIR can simulate up to 24 plant function types (PFTs), with an additional bare land type, which are the same as the Community Land Model version 4.5 (CLM4.5). For these 24 PFTs, default TEMIR can calculate photosynthesis and stomatal conductance based on the Farquhar-Ball-Berry (FBB) model as described by Oleson et al. (2013), which could capture the environmental conditions including radiation, temperature, and soil water availability. Eight of the 24 PFTs are human-managed crops, including maize, soybean, winter cereal and temperate cereal, each being either rainfed or irrigated. As opposed to natural vegetation, humanmanaged crops have their unique phenology that controls the allocation of the assimilated carbon. They have a separate carbon stock that keeps track of the carbon allocated to the pod, which is turned into a crop yield at the end of the growing season. In addition to the default stomatal conductance scheme that is fully coupled with photosynthesis, we also implement a Jarvistype stomatal conductance scheme that is consistent with the Deposition of O 3 for Stomatal Exchange (DO 3 SE) model (https:// www.sei.org/projects-and-tools/tools/do3se-deposition-ozonestomatal-exchange/), which has been developed to estimate the risks of O 3 damage to European vegetation according to Long-Range Transboundary Air Pollution (LRTAP) methodologies for effect-based risk assessment. Both FBB-based and DO 3 SE-based stomatal conductances are used to calculate O 3 fluxes and POD Y . To stay consistent with Mills et al. (2018a), we choose Y = 3 nmol O 3 m −2 s −1 in this study.

Crop Land Coverage and Production
Distribution of croplands is based on the year-2000 global land surface from the Community Land Model version 5 (CLM5), replacing the cropland distribution used in CLM4.5 whereby crops are missing in the tropical regions. The new input is then regridded to a resolution of 2 • latitude by 2.5 • longitude consistent with GEOS-Chem, which provides the percentage coverage of major crops including maize, soybean, wheat and rice. In TEMIR, crops are simulated when the land coverage of the crops is >1%. Figure 1 shows global distribution of average production of the four staple crops; data are from the Food and Agriculture Organization (FAO), mapped onto the 2 • × 2.5 • grid using a data fusion technique following Tai et al. (2014).

Crop Calendar
Crop calendar is based on the global dataset from Sacks et al. (2010), which provides planting and harvesting dates of major crops including maize, soybean, wheat and rice with a resolution of 0.5 • longitude by 0.5 • latitude. For crops that have multiple growing season in a year (e.g., maize and rice in the tropical regions), only the primary growing season is considered in this study. With these inputs, planting dates can be prescribed in our simulations and harvesting dates can be either prescribed or estimated using growing degree day (GDD) (section Exposure-Yield Relationships). Supplementary Figure 1 shows the harvesting dates used in our study.

Ozone and Meteorological Data
Decidedly, the same O 3 concentration and meteorological fields are used throughout to ensure inter-model differences are methodological but not input-driven. The present-day O 3 concentration field is simulated by the GEOS-Chem global CTM (http://acmg.seas.harvard.edu/geos/), which simulates fully coupled O 3 -NO x -VOC-aerosol chemistry together with the emissions, transport and deposition of important atmospheric chemical species. For our simulations of O 3 , GEOS-Chem is driven by assimilated meteorological data at 2 • × 2.5 • resolution from the reanalysis product Modern-Era Retrospective analysis for Research and Applications, Version 2 (MERRA-2), which is also used as the input meteorology of TEMIR to ensure internal consistency in the TEMIR-GEOS-Chem framework. The surface O 3 concentrations simulated by GEOS-Chem are then used to directly calculate the various concentration-based exposure metrics, and also together with MERRA-2 meteorology as inputs for TEMIR to simulate the corresponding POD 3 , photosynthetic damage and crop yield losses.

Ozone Metrics and Damage Parameterization
Concentration-Based Exposure Metrics: AOT40, M7/M12, and W126 In this study, we consider three different ozone exposure metrics widely used to quantify ozone exposure and impacts on crops. They include two cumulative exposure metrics (AOT40 and W126) and one mean exposure metrics (M12 or M7). They are defined as: is the hourly mean O 3 concentration during the 12 h of local daylight (08:00-19:59); n is the number of hours in the 3-months growing season defined as the 3 months prior to the start of the harvesting period in every country according to the crop calendar; and w i is a weighting function assigning greater weight to higher levels of hourly O 3 (Lefohn et al., 1988), defined by: AOT40 and W126 are in units of ppm-h, and M12 in units of ppb. M7 is defined similarly to M12 but uses a 7-h (09:00-15:59) daytime window instead of 12-hour for summation. To remain consistent with the flux-based metrics, both cumulative exposure metrics and mean exposure metrics are cumulated and averaged, respectively, consistently over the same growing period, from 104 days before the prescribed harvesting date to 14 days before the prescribed harvesting date (Supplementary Figure 1). Supplementary Figures 2-4 show the global distribution of these metrics as calculated from the O 3 concentration field simulated by GEOS-Chem for year 2010.
Flux-Based Exposure Metrics: POD 3 -DO 3 SE and POD 3 -FBB In order to calculate the flux-based metrics, POD Y , we simulate stomatal conductance (g s ) using both the default TEMIR method of directly coupling g s with net photosynthesis rate (A n ), and the semiempirical, Jarvis-type formulations as in DO 3 SE. For both methods, leaf-level POD Y (mmol m −2 ) is calculated as: where the threshold in this study is Y = 3 nmol O 3 m −2 s −1 , t is the time step ( t = 3,600 s in TEMIR), and the O 3 flux, F O3 (nmol O 3 m −2 s −1 ), minus the threshold is cumulated over a 90-days growing period, from 104 days before the prescribed harvesting date to 14 days before the prescribed harvesting date (Supplementary Figure 1).
The O 3 flux is calculated as: is the ozone concentration observed or of the lowest atmospheric model layer. R a (s m −1 ) is the aerodynamic resistance, and R b (s m −1 ) is the quasilaminar (boundary-layer) resistance. We use the same set of equations to be consistent with DO 3 SE, which are detailed in Supplemental Information. Both R a and R b are bulk resistances. In contrast, the stomatal resistance, r s (mmol −1 O 3 m 2 s) = 1/g s , is at the leaf level, and its unit is converted to s m −1 using a conversion factor α ≈ 42,000 mmol m −3 . The formulation of stomatal conductance for O 3 , g s (mmol O 3 m −2 s −1 ), is what differs between FBB and DO 3 SE. In the default TEMIR FBB calculation, g s is coupled with net photosynthesis rate (A n , µmol CO 2 m −2 s −1 ) via: where g s is controlled by the leaf surface CO 2 partial pressure c s (Pa), leaf surface water vapor pressure e s (Pa) and Frontiers in Sustainable Food Systems | www.frontiersin.org temperature-dependent saturation vapor pressure e sat (Pa); m = 9 and b = 10,000 µmol m −2 s −1 for C 3 plants, and m = 4 and b = 40,000 µmol m −2 s −1 for C 3 plants; β t (ranging from 0 to 1) is the soil water stress function that accounts for water limitation, and this factor is also applied internally to A n during its computation; the factor 1.61 × 10 3 is the ratio of leaf conductance for water The conductance g s in FBB is calculated separately for sunlit and shaded leaves (i.e., g ssun and g ssha , respectively) and the canopy-averaged leaf-level g s is: where LAI sun and LAI sha are the leaf area index (LAI) of sunlit and shaded leaves, respectively. In contrast to the photosynthesis-based formulation, the stomatal conductance g s (mmol O 3 m −2 s −1 ) in the default DO 3 SE model is a semi-empirical multiplicative function: The definitions of different terms in Equation (12) that account for variations in various environmental factors, including plant-type-specific phenology (f phen ), light (f light ), temperature (f T ), vapor pressure deficit (VPD, f D ), and soil water potential (SWP, f SW ), are listed in Supplementary Table 1. Figures 2, 3 show the simulated POD 3 with the DO 3 SE multiplicative formulation and default FBB photosynthesisbased formulation, respectively. In general, they show very similar spatial patterns for all four crops, with the DO 3 SE method generally yielding higher POD 3 than FBB, reflecting systematically higher values of g s simulated by DO 3 SE. Table 1 lists the statistical relationships of crop yields with the concentration-based exposure metrics (AOT40, M7/M12, W126) and flux-based metrics POD 3 , summarized from various sources. These metrics have been used extensively to estimate global yield losses (e.g., Mills et al., 2018a). The statistical relationships between the flux-based metrics and crop yields are usually called "dose-response" functions; those between the concentrationbased metrics and crop yields are likewise called "concentrationresponse" functions. Relative yield (RY) is defined to be the ratio of O 3 -affected yield to unaffected yield at zero O 3 exposure. In this study, we find RY for each 2 • × 2.5 • grid cell i; we then calculate the globally or regionally aggregated yield loss (L, %)  from the grid cell-by-grid cell relative yield (RY i ) loss following:

Exposure-Yield Relationships
where P 0i is the theoretical production in each grid cell i without O 3 damage, P i is the actual production as shown in Figure 1, and the summation is over the whole world for globally aggregated estimates, or over the whole region for regionally aggregated estimates. This approach of aggregating regional and global yield (more accurately, production) losses is consistent with previous studies (Ainsworth, 2017;Mills et al., 2018b) with which we will compare our results. We note however that these different metrics were derived from different experiments that might not always have consistent study designs or environmental controls even for the same crop, which may further add to the uncertainties represented merely by the use of different metric definitions.

Mechanistic Crop Model and Flux-Based Ozone Damage Parameterization
Finally, we also attempt to estimate crop yield losses based on a biophysical crop model built on the TEMIR model, whereby the instantaneous stomatal O 3 uptake directly affects the simulated rate of photosynthesis and/or stomatal conductance, and subsequently impacts on crop yield via modifying carbon fixation and allocation. We may call it the "plant response-based" approach, and this is as opposed to using the cumulative O 3 uptake (i.e., POD Y ) to directly estimate relative yield (i.e., using dose-response functions). To simulate the impacts of O 3 on crop growth and yield, here we adopt the O 3 damage parameterization developed by Sitch et al. (2007) (Sitch scheme). As with the computation of POD Y , the parameterization of the Sitch scheme makes use of the instantaneous stomatal O 3 flux (F O3 ), but instead of computing a cumulative POD Y from F O3 , it considers a photosynthetic damage factor f (i.e., fractional reduction) that is a function of F O3 : where F crit represents a critical threshold accounting for O 3 tolerance, below which instantaneous O 3 exposure does not affect photosynthesis, and F crit = 5 nmol m −2 s −1 for all crops; the O 3 sensitivity parameter b (nmol −1 m 2 s) is specific to the crop with different photosynthetic pathways. Soybean is treated as a C3 grass and maize is treated as a C4 grass in this study (Supplementary Table 3). Fractional reduction factor f is multiplied directly to the net photosynthesis rate A n to represent O 3 damage, which then indirectly affects g s via the coupling between A n and g s ; f, A n, and g s are solved together iteratively. Crop yield (Y, gC m −2 ) can be inferred from the accumulated assimilated carbon (i.e., net primary production or NPP) during Frontiers in Sustainable Food Systems | www.frontiersin.org the growing season, which is given as: where R is the plant respiration rate that depends on the total biomass and A n ; a repr is the biomass partitioning fraction toward plant reproductive organs that depends on growing degree days (GDD) which is calculated by the accumulation of temperature. Similarly, leaf area index (LAI) can be calculated using the above formula by considering the carbon stock of leaf for the whole plant: where SLA is the specific leaf area (m −2 gC) of the plants and L sen is the loss of leaf carbon during the later growing stage. The parameterization of L sen in our model follows the CLM4.5 crop model (Oleson et al., 2013), in which the leaf loss rate is proportional to the leaf carbon stock. Using these two equations, we can infer that O 3 can reduce the NPP and thereby crop yield through reducing the leaf-level photosynthetic rate as well as the LAI. The reduction of A n due to the O 3 damage during the early growing stage can result in loss in leaves as plant prioritizes carbon allocation toward leaf. Whereas, in the later growing stage when plants prioritize the allocation toward pods, it can result in less carbon assimilation for the grain. These two pathways also apply for crops that are grown under elevated CO 2 levels, which per se generally increase crop yield by enhancing A n . In addition, higher CO 2 levels can lead the reduction of g s , which in turn indirectly lowers the uptake of O 3 through the stomata and thus reduce O 3 damage. These interactive pathways involving CO 2 and O 3 within crops can only be resolved using a mechanistic model like this one described above. Important parameterization schemes and parameter values are listed in Supplementary Table 2.
We also test another plant response-based approach using cumulative (instead of instantaneous) uptake of O 3 (CUO), which is similar to POD, to alter photosynthesis rate and stomatal conductance. The CUO-A n and CUO-g s relationships follow the linear statistical relationships consolidated by Lombardozzi et al. (2013) for different plant types. However, the specific statistical relationships for crops only have modest sensitivity of A n or g s to CUO beyond a given threshold, and when implemented in TEMIR-crop they result only in minimal damage on crop yield. While theoretically using CUO to construct plant responsebased metrics may be preferable to using instantaneous uptake, we resort to use the Sitch scheme only for comparing with the exposure metrics and do not calibrate the Lombardozzi scheme further, as this part of the study aims only to demonstrate how exposure metrics are unable to capture the effects of rising CO 2 .    Figures 5-9 show the maps of relative yield (RY) based on the three concentration-based and two flux-based O 3 exposure metrics. While as in the case for the metrics themselves RY for a given crop show generally consistent spatial distribution across different metrics, the magnitudes of RY appear to differ quite substantially. Below we diagnose these differences and compare them with results from previous studies both globally and regionally. Figure 10 shows the aggregated global yield losses from our study for the four crops and different exposure metrics, as well as summarized results from Mills et al. (2018b), which used the flux-based approach similar to POD 3 -DO 3 SE to estimate global losses, and Ainsworth, 2017, which presented from AOT40derived and M7/M12-derived estimates from van Dingenen et al. (2009) and Avnery et al. (2011). Even with a consistent set of O 3 concentrations as input into the calculation, our estimates of global yield losses from different metrics differ quite substantially for a given crop. W126 almost always gives the lowest losses, almost appearing to be outlying compared to the other metrics, except for wheat. The concentration-based metrics differ greatly among themselves. In contrast, the two flux-based metrics are generally close to each other, and lie close to the middle of the range covered by all metrics. POD 3 -DO 3 SE always gives slightly higher loss estimates than POD 3 -FBB, mostly reflecting the generally lower stomatal conductance calculated by photosynthetic coupling (as in FBB) than a Jarvis-type multiplicative approach (as in DO 3 SE). Generally either AOT40 or M7/M12 would give the largest losses for a given crop. When compared with estimates from previous studies, our estimates of global yield losses using AOT40 and M7/M12 are generally quite consistent with that summarized by Ainsworth (2017), except for wheat, whereby our AOT40 estimate is higher but M7/M12 estimate is lower than theirs. Our flux-based estimates are, however, always quite substantially lower than that found by Mills et al. (2018b). Their estimates are almost always on the high end of the range given by all metrics, except for wheat. The larger differences between our flux-based estimates and theirs may reflect a multitude of factors, but since the methodological approaches are largely similar and we use the same dose-response relationships, we ascribe the differences primarily to the differences in O 3 and meteorological fields used as inputs for the calculation. Figure 11 shows the regionally aggregated yield losses for four major crop-producing regions (our estimates only). As in the case for globally aggregated losses, AOT40 and/or M7/M12 generally lead to the largest yield losses for most regions and crops, W126 leads to the smallest (except for wheat), and the fluxbased metrics lead to losses in the middle of the range. What is particularly noteworthy is that for a given crop, the differences in yield losses among the metrics themselves in the same region can be often larger than the inter-regionally differences. This suggests again that the greatest uncertainty behind O 3 -induced yield losses appears to come from methodological differences, which can be even larger than that from the differences in input data. Comparing across different crops, soybean and wheat appear to suffer greater losses from O 3 damage than maize or rice, which is consistent with previous studies (e.g., Ainsworth, 2017;Tai and Martin, 2017;Mills et al., 2018b). However, when we carefully compare across metrics, it appears that very high yield losses for soybean and wheat are mostly driven by a single metric (M7/M12 for soybean, AOT40 for wheat); the average losses over all metrics and the losses derived from the flux-based metrics for soybean and wheat, though still higher than, are generally closer to that for maize or rice.

Using a Mechanistic Crop Model to Evaluate O 3 -Induced Yield Losses and Effects of CO 2 on O 3 -Induced Yield Losses
Finally, we present findings using a mechanistic crop model, which can help examine the direct impacts of O 3 on plant physiology whereas using the metrics cannot. Only the results for maize and soybean are presented due to the lack of validation data for wheat or rice. Under the influence of O 3 , global mean photosynthesis rate under present-day O 3 concentrations and ambient CO 2 concentration relative to the scenario without O 3 is reduced by 18-34% for maize and 33-56% for soybean (Supplementary Figure 5), depending on the sensitivity of crops to O 3 . Reduction of crop yield couples tightly with the reduction of photosynthesis rate, and yield decreases  by 21-38% for maize and 36-59% for soybean (Figure 12). Another impact of O 3 on plant physiology is the reduction of LAI (Supplementary Figure 6), whereby our results suggest that maximum growing season LAI is reduced by about 1-2 under low plant-O 3 sensitivity. Compared to other major maize-planting region, maize in Europe has higher maximum growing LAI and reduction in maximum LAI caused by O 3 , which can be attributed to the overestimation of the length of the vegetative stage. The Sitch scheme used to estimate O 3 damage comes with a range of damage with low and high sensitivity to O 3 . Based on the relative yield given by the model, we find the globally aggregated yield loss for maize to be 9.5 and 41% with low and high O 3 sensitivity, respectively, at preset-day CO 2 level. For soybean, the global yield loss is estimated to be 26 and 63%, respectively. These losses are larger than that given by the exposure metrics presented above in Section Comparing Different O 3 Exposure Metrics and Yield Losses in which maize loss is around 0-5% and soybean loss is around 2-16%, with the low-sensitivity estimates being closer and more reasonable as compared with the metric-based damage. We further compare our simulation results for soybean with the analysis conducted by Morgen et al. (2003), whereby crop yield decreases only by around 6 ± 5% when O 3 level is between 30 and 59 ppb and 24 ± 4% when O 3 level is between 60 and 79 ppb. For maize, based on the study by Peng et al. (2019), observed yield changes were +4.1, −2.7, and −16.6% for maize growth under daylight O 3 levels of 39, 60, and 76 ppb, respectively. In addition, the effect of O 3 on peak LAI is not significant unless plants are exposed to >100 ppb of O 3 for a long time (Dermody et al., 2006;Betzelberger et al., 2012). This suggests that even with the low-sensitivity parameters, the Sitch scheme may be overestimating the sensitivities of soybean and maize to O 3 .
At 600-ppm CO 2 level, O 3 -induced yield losses are found to be always smaller. The globally aggregated yield loss for maize is 3.5 and 5.6% for low and high sensitivity, respectively, and for soybean is 19 and 29% for low and high sensitivity, respectively. Comparing between the O 3 -induced yield losses at 390 vs. 600 ppm CO 2 , the higher CO 2 concentration can reduce the maize yield loss by 21 and 52% for low and high sensitivity, respectively, and reduce the soybean yield loss by 27 and 38% for low and high sensitivity, respectively. The results from the SoyFACE experiments suggest that at elevated O 3 (+25% ambient level), soybean relative yield losses are 10 ± 11% and 5 ± 4% under 370 and 550 ppm of CO 2 , respectively; i.e., higher CO 2 concentration reduces O 3 -induced yield loss by roughly half (Morgen et al., 2005(Morgen et al., , 2006. However, the uncertainty of yield losses is still large, which could be attributable to the insufficient numbers and durations of experiments conducted (3 years at elevated CO 2 ; 5 years at ambient CO 2 level).
Impacts of higher CO 2 concentration level (600 ppm) on plant physiology is reflected on the change of stomatal conductance; global leaf-level stomatal conductance is reduced by 14% for maize and 31-33% for soybean relative to 390 ppm CO 2 (Supplementary Figure 7), which indicates that O 3induced yield loss is less severe under higher CO 2 concentration (Supplementary Figures 8-11). Compared to the studies by Ainsworth et al. (2002) and Leakey et al. (2004), the reductions in stomatal conductance under elevated CO 2 levels are consistent for both maize (−23% at 550 ppm relative to 350 ppm CO 2 A B C D in their studies) and soybean (−36% at 600 ppm CO 2 in their studies), respectively. Supplementary Figures 8-11 show the regression between RY and AOT40 at different crop planting regions. AOT40 is used here only as an indicator of O 3 exposure for better comparison with the exposure metrics. It shows that for most of the regions, RY is higher when the CO 2 level is higher. It also suggests that the respond of RY to different AOT40 levels varies in different regions, which can be explained by the fact that AOT40 does not take account into the O 3 effect during the early growing stage.

CONCLUSIONS AND DISCUSSION
The exposure to O 3 , a phytotoxic chemical, brings tremendous harm to crops, but their impacts on global and regional crop yield even for the present day are far from certainty. The large uncertainty may first stem from: (1) large spatiotemporal variability of O 3 concentration itself that is not adequately quantified, due to incomplete observations and inter-model differences in O 3 -relevant process representation in different models; (2) different methodological approaches to quantify O 3 exposure and impacts, including the concentration-based approaches, flux-based approaches, and mechanistic crop models; (3) co-varying factors such as CO 2 concentration and other environmental conditions (heat stress, water stress, etc.) that are not explicitly accounted for in many of the commonly used approaches. In this study, we therefore attempt to explore these issues using a common methodological framework: a consistent set of simulated present-day O 3 fields from one chemical transport model, coupled with a terrestrial ecosystem model (including crops) to derive various O 3 exposure metrics and impacts on relative crop yields worldwide. Throughout, we FIGURE 13 | Summary of TEMIR-simulated relative yield at four different regions under 390 and 600 ppm CO 2 concentration. The error bar represents one standard deviation of the data. Data that is above the 95th percentile or below the 5th percentile is removed.
also review and explain in details the differences in formulation and parameterization in different approaches. Finally, we also explore the effects of elevated CO 2 on O 3 damage to crops. We find that while the spatial pattern of yield losses for a given crop is generally consistent across metrics, the magnitudes can differ substantially. Sometimes, the inter-metric differences can be greater than the inter-regional differences. Overall, pulling two concentration-based (AOT40, M7/M12) and two flux-based (POD 3 -DO 3 SE, POD 3 -FBB) together (excluding W126 due to its outlying tendency and lack of comparison with other studies), we estimate the globally aggregated yield losses for the four staple crops to be (± standard deviation): 3.6 ± 1.1% for maize, 2.6 ± 0.8% for rice, 6.7 ± 4.1% for soybean, and 7.2 ± 7.3% for wheat. The average yield losses we find are generally consistent with previous studies (e.g., Ainsworth, 2017;Mills et al., 2018b), but on the lower end of the full uncertainty range previously covered (especially for wheat). We attribute the relative large combined uncertainty we find mostly to the differences among methodological approaches (even within our own study with consistent inputs), and secondarily to differences in inputs (as compared with other studies).
The exposure metrics, either concentration-based or fluxed based, which are almost always derived from wellcontrolled, single-factor experiments, are generally ill-suited to account for the effects of co-varying and possibly confounding environmental factors, including CO 2 concentration, temperature, precipitation, drought frequencies or intensity, etc. We further explore how rising ambient CO 2 concentration in the future may affect the impacts of O 3 on crop yields using a mechanistic biophysical crop model that explicitly simulates the effects of O 3 uptake on photosynthesis rates and thus downstream on crop yields (i.e., plant response-based metrics. Crop ecophysiology embedded within the crop model is able to account for co-varying factors such as CO 2 , temperature, soil water availability. We find that, first, our mechanistic model as of the current version gives very different, much larger relative yield losses than the exposure metrics; this is expected as the ozone-induced damage parameterization that we have adopted in the simulations is generalized from studies with C3 and C4 grasses, with a lack of data for more refined calibration in different regions for specific crops. We should therefore only focus on the relative differences between two different CO 2 levels to infer their impacts. We find that increasing CO 2 concentration from 390 ppm (representative of 2010) to 600 ppm could reduce stomatal conductance by 14% for maize and 31-33% for soybean, consistent with the finding from SoyFACE experiments. The globally aggregated O 3 -induced yield loss is reduced by 21-52% for maize and by 27-38% for soybean. Although the reduction in O 3 -induced yield loss for soybean is smaller than the fieldobserved value (54%), it still reflects a CO 2 -induced reduction in stomatal conductance that in turn quite substantially alleviates stomatal O 3 uptake and thus damage to yields. Rising CO 2 may therefore render the currently used exposure-yield relationships less applicable in a future atmosphere.
Many modeling studies have used mechanistic crop models to predict future crop yields under rising CO 2 concentrations and changes in other co-varying factors such as climatic heat stress Rosenzweig et al., 2014); elevated CO 2 is generally found to substantially enhance crop yields via the fertilization effect. However, modeling results for the combined effects of CO 2 and O 3 pollution on different crops are thus far limited [e.g., Tao et al. (2017) and Guarin et al. (2019) for wheat only; Lombardozzi et al. (2018) for eight crop types], and observations or experiments that may validate the wide applicability of mechanistic crop models for the coeffects of CO 2 and O 3 are scarce and limited to few crops [e.g., Dermody et al. (2006) for soybean; Rudorff et al. (1996), Cardoso-Vilhena et al. (2004);and Feng et al. (2008) for wheat]. Meanwhile, although the exposure metrics (especially flux-based metrics such as POD 3 ) have found success in assessing current and future impacts of O 3 pollution on global crop yields and food security, they are inherently incapable of capturing the co-effects of CO 2 , which calls into questions their applicability in a future atmosphere enriched with CO 2 . Ultimately, greater effort should be put into conducting experiments for a basket of globally important crops that explicitly consider co-varying levels of O 3 and CO 2 as well as co-varying hydrometeorological conditions, and the experimental results can then be used to derive revised exposure-yield relationships that depend on changing CO 2 concentration and climatic stress. Before that is possible, a shortcut we suggest is to rely on the mechanistic understanding of CO 2 -O 3 co-effects embedded in either an FBB-like photosynthesis-stomatal conductance model or a biophysical crop model. For instance, since the FBB model explicitly considers the reduction in stomatal conductance and thus O 3 uptake at elevated CO 2 , a framework similar to POD 3 -FBB can be used directly to evaluate the CO 2 -O 3 co-impacts on crop yields based on existing POD 3 -yield relationships. This approach, however, can only account for the effect of CO 2 via modifying stomatal conductance, but cannot represent the CO 2 fertilization effect on crop yields via enhanced photosynthesis. Only a fully mechanistic biophysical crop model (e.g., TEMIR-crop, CLM-crop, JULES-crop, etc.) can incorporate the CO 2 fertilization effect. Thus, an alternative approach is to use such a crop model to derive scaling factors that account for the variations in simulated O 3 -yield relationships at different CO 2 levels; these factors can then be used to scale existing exposure-yield relationships. These approaches may be valuable in assessing the global impacts of O 3 pollution on crop production and food security over the next coming decades as ambient CO 2 concentration is expected to continue to rise, particularly before more of the long-term CO 2 −O 3 covarying experiments for a larger variety of crops can be done.

DATA AVAILABILITY STATEMENT
Data used for generating the figures in this paper can be found online (http://www.cuhk.edu.hk/sci/essc/tgabi/data.html). Full datasets can be made available on request to the corresponding author.

AUTHOR CONTRIBUTIONS
AT coordinated the study, performed data analysis, and wrote the manuscript. MS performed model development and simulations of ozone exposure metrics and yield losses. JP performed model development and simulations of biophysical crop model. DY simulated the ozone concentration fields. ZF performed data analysis and assisted in writing. All authors contributed to the article and approved the submitted version.