New Insights into Leaf Physiological Responses to Ozone for Use in Crop Modelling.

Estimating food production under future air pollution and climate conditions in scenario analysis depends on accurately modelling ozone (O₃) effects on yield. This study tests several assumptions that form part of published approaches for modelling O₃ effects on photosynthesis and leaf duration against experimental data. In 2015 and 2016, two wheat cultivars were exposed in eight hemispherical glasshouses to O₃ ranging from 22 to 57 ppb (24 h mean), with profiles ranging from raised background to high peak treatments. The stomatal O₃ flux (Phytotoxic Ozone Dose, POD) to leaves was simulated using a multiplicative stomatal conductance model. Leaf senescence occurred earlier as average POD increased according to a linear relationship, and the two cultivars showed very different senescence responses. Negative effects of O₃ on photosynthesis were only observed alongside O₃-induced leaf senescence, suggesting that O₃ does not impair photosynthesis in un-senesced flag leaves at the realistic O₃ concentrations applied here. Accelerated senescence is therefore likely to be the dominant O₃ effect influencing yield in most agricultural environments. POD was better than 24 h mean concentration and AOT40 (accumulated O₃ exceeding 40 ppb, daylight hours) at predicting physiological response to O₃, and flux also accounted for the difference in exposure resulting from peak and high background treatments.

1. Introduction

The air pollutant ozone (O3) reduces yield in many crops including wheat, rice, and soybean [1,2]. Ozone at the ground level forms from precursor gases—chiefly NOx and volatile organic compounds (VOCs)—in chemical reactions catalysed by sunlight and heat [3]. Concentrations over much of the Earth’s land surface have approximately doubled since pre-industrial times, mainly due to anthropogenic emissions from vehicles, industry, and agriculture [4,5,6]. Annual mean surface O3 concentrations have largely stabilized in Europe since the year 2000 as a result of emission control policies [7,8], but continued increase to 2050 is likely across South and East Asia [9,10]. The pattern of O3 exposure across regions is also expected to change over coming decades. Short peak ‘episodes’ of very high concentrations are predicted to become more frequent in India and China [10,11], while in Europe and North America a decline in peak episode frequency, alongside steadily increasing annual mean O3 concentrations, was observed between 1990 and 2010 [12]. Modelling suggests this decline in peak episode frequency in Europe and North America is likely to continue to 2050 [10].

Model-based estimates of yield loss under future climate and air pollution scenarios represent a powerful way of highlighting yield benefits that could come from reduced surface O3 [13,14,15,16]. Global O3-induced wheat yield loss for the year 2000 has been estimated as ranging from 5% to 26%, with potential additional losses of 1.5% to 10% predicted for 2030 [13,14]. However, nearly all large-scale assessments of O3-induced yield loss for wheat published to date have followed an empirical approach, where O3 concentration is simulated using a chemical transport model (CTM), concentration is linked to yield loss using response functions, and response is scaled up using crop production maps and agricultural statistics. An alternative is to use a process-based approach that could potentially produce more robust estimates of future yield, through inclusion of interactive effects between O3, CO2, and climate variables [17].

This process-based crop modelling approach builds on the development of methods for modelling O3 flux into leaves that have been made over recent decades [18,19]. These methods provide an hourly estimate of O3 dose reaching sites of damage in the leaf, creating the potential for O3 effects to be integrated into crop simulation models in a dynamic way. Studies applying O3 flux modelling have generally either used a multiplicative stomatal conductance (gsto) algorithm—described in prior research [20,21,22]—or followed a semi-mechanistic approach where gsto is estimated empirically from photosynthetic rate, which in turn is modelled using the biochemical model of Farquhar et al. [23], also described in [24,25]. Since most crop models simulate growth responses at daily (and less frequently, hourly) timesteps and can respond to a changing environment [26], integration of O3 effects into crop models is feasible, if plant response to O3 can be represented in the model formulation. Attempts have been made to integrate O3 effects and crop modelling [25,27], but few estimates of O3-induced crop yield loss using a dynamic approach have been published to date [28,29]. Reasons for slow progress include the challenge of upscaling responses from the leaf to the canopy, the need for species and cultivar-specific model parameters, and the incomplete understanding of physiological mechanisms driving O3-induced yield reduction [17]. This paper goes some way to addressing these issues and identifying future research direction that would benefit from empirical investigation targeted towards developing the models that are currently in development.

It is well established that O3 exposure can reduce yield in wheat [30,31,32] and can cause foliar injury, impaired photosynthesis, altered carbon translocation, and accelerated senescence [33,34]. However, the processes linking O3 uptake to these responses are not fully understood and it is not clear which are most important in driving ultimate yield loss. Once O3 has been taken up through stomata, reactions in the plant apoplast lead to the formation of reactive oxygen species (ROS), which can then react with and damage membranes and proteins [35]. Most plants have in-built defence mechanisms and can up-regulate antioxidants to detoxify ROS, but this comes at a carbon cost, meaning O3 damage to productivity often occurs before visible symptoms appear [35,36]. An O3-induced reduction in photosynthetic rate has been widely reported [32,37,38], but quantifying the extent to which this is a direct effect of O3 on the photosynthetic mechanism, or indirect via changes to leaf pigmentation or gsto, has been a challenge for experimentalists. Disentangling direct O3 impacts on photosynthesis from the accelerated senescence response is also difficult. Some studies have observed reduced activity of the carbon-fixing enzyme ribulose-1,5-biphosphate carboxylase/oxygenase (rubisco) in response to O3 [39,40,41], leading to the hypothesis that ‘instantaneous’ effects of O3 on photosynthesis act via effects on this enzyme. The physiological mechanism underpinning the often-observed accelerated senescence response to O3 [42,43,44,45,46] is also unknown, although it has been hypothesized that it relates to long-term respiratory costs associated with detoxification and repair [25].

Several approaches for modelling O3 effects on photosynthesis and senescence have been published. In an early paper, Reich et al. proposed that ozone effects on a number of plant groups could be expressed via a linear relationship between exposure and growth [47]. Subsequent published approaches have attempted to model the separate effects of O3 on productivity and senescence and have tried to account for differential sensitivity across species. A function for modelling ‘instantaneous’ suppression of photosynthesis was proposed by Martin et al. [27], who simulated a linear reduction in carboxylation capacity of rubisco (the parameter Vcmax in the model of Farquhar et al. [23]), above a threshold hourly flux value representative of the species or cultivar-specific detoxification capacity. A similar approach of O3 acting on Vcmax was also used by Deckmyn et al. in their O3 damage module for forest trees, alongside an overnight repair mechanism and a parameter representing the respiratory cost of detoxification [48]. Ewert and Porter [25] applied a version of the Martin et al. ‘short-term’ function alongside a ‘long-term’ algorithm for modelling O3-induced senescence and assumed that ‘short-term’ photosynthetic suppression by O3 occurs throughout the leaf lifespan. Their senescence function assumes a linear reduction in mature leaf lifespan as accumulated O3 flux increases, and senescence comprises the final third of the mature leaf lifespan, during which time Vcmax is assumed to decline linearly. In this function, onset and completion of leaf senescence therefore move progressively earlier and closer together as accumulated O3 flux increases. An alternative approach for modelling O3-induced senescence is applied in the multiplicative DO3SE model (Deposition of ozone for stomatal exchange), a gsto model which estimates accumulated O3 flux—known as the Phytotoxic Ozone Dose (POD)—to vegetation [18]. In this model, leaf senescence is induced by a threshold POD, which triggers curvilinear decline in leaf gsto with a fixed shape but variable decline rate [20,49,50]. The POD ‘trigger’ can be parameterized according to the sensitivity of the cultivar or species.

The integration of O3 damage functions, such as those described above, into crop models could improve yield estimates under O3 stress. Model development must however be guided by experimental evidence that identifies damage mechanisms and their relative importance, informs parameterization, and indicates likely degree of error. Models must also be able to replicate the response to different patterns of O3 exposure. For example, acute peaks in concentration have been observed to induce greater yield loss than consistent, moderate levels with the same 24-h mean exposure [51], and modelling methods need to be able to capture these nuances.

This study combines data from two independent O3 exposure experiments on European wheat that took place at the same experimental facility in 2015 and 2016 and used consistent protocols for data collection. We model the accumulated O3 flux to wheat flag leaves across different O3 treatments, using the multiplicative DO3SE model, to derive the POD0SPEC metric of O3 flux exposure (mmol m−2 PLA day−1). Exposure to O3 flux is then considered alongside leaf chlorophyll, gsto, and photosynthesis responses in order to test key assumptions underpinning published O3 effect model functions. Firstly, with regards to O3 effects on leaf senescence, we (i) examine whether inter-cultivar differences in response are captured by current senescence functions and (ii) whether leaf senescence begins at an accumulated O3 flux ‘trigger’ value. Secondly, we examine whether O3 reduces Vcmax before—and therefore independent of—onset of O3-induced leaf senescence. Thirdly, we investigate if flux is a better predictor of the physiological response to O3 than concentration-based metrics and whether flux can account for differences in the pattern of O3 exposure (i.e., peak vs. background).

Our results reveal several insights about physiological responses to O3 that can add to the evidence base for designing O3 effect model functions. The two cultivars of European wheat studied here showed substantially different sensitivity to O3 in terms of their senescence response, indicating the importance of cultivar-specific parameterization in senescence functions. The study also finds that O3 effects on the photosynthetic mechanism are not observed in young flag leaves and are only seen following the onset of leaf senescence, suggesting that O3-induced accelerated senescence is more important than direct effects on photosynthesis in determining final yield.

2. Results

2.1. Ozone Treatments in 2015 and 2016

This study combines data from two independent experiments, which took place at the Centre for Ecology and Hydrology (CEH) air pollution exposure facility in 2015 and 2016. Both experiments gathered data on the physiological response of wheat cultivars to O3 exposure using consistent experimental design and data collection techniques. The main differences in these experiments were in the ozone concentration profiles to which plants were exposed. In 2015, two wheat cultivars, ‘Mulika’ and ‘Skyfall’, were exposed to eight O3 treatments. The treatments varied in their profile, with some characterized by acute peaks in concentration and others by a consistent background level (described in full in Section 4.1). In 2016, cv. ‘Skyfall’ was exposed to five O3 treatments, all with a ‘peak’ style profile. The treatments are named in this paper according to their severity and profile, these are ‘low background’ (LB), ‘low peak’ (LP), ‘medium background’ (MB), ‘medium peak’ (MP), ‘high background’ (HB), ‘high peak’ (HP), ‘very high background, (VHB), and ‘very high peak’ (VHP).

Table 1 summarizes the O3 treatments using several exposure indices including the 24-h mean (ppb) and seasonal AOT40 (ppm h). It also quantifies exposure using O3 flux, or phytotoxic ozone dose (POD), which was modelled in this study using the multiplicative DO3SE model (described in more detail in the Methodology and Appendix A). Several different flux metrics are included in the table. Because most plants have the ability to detoxify a small quantity of O3, methods for quantifying O3 dose often use a threshold for flux accumulation with the threshold value representing the detoxification capacity. Previous experimental work in wheat has found that a threshold of six produces the closest fit between exposure and yield [50]. We therefore applied a threshold of six when calculating accumulated flux, to produce the POD6SPEC metric (where SPEC refers to the species-specific nature of the metric). However, as thresholds for physiological effects are not as well established, we have also calculated accumulated flux with no threshold, known hereafter as the POD0SPEC. The mean daily POD0SPEC was also calculated, to represent the average level of exposure on a given day in each treatment.

Table 1

Summary of ozone (O3) treatments administered in the 2015 and 2016 experiments. 24-h mean, AOT40, and mean daily peak O3 have been calculated over the full O3 exposure period, whereas the mean POD0SPEC, POD0SPEC, and POD6SPEC quantifies exposure in the flag leaf only (i.e., calculated over the period following flag leaf emergence).

Season	Ozone Treatment	24-h Mean (ppb)	AOT40 (ppm h)	Mean Daily Peak O3 (ppb)**	Mean Daily POD0SPEC (mmol m−2 PLA day−1)	POD0SPEC (mmol m−2 PLA)	POD6SPEC (mmol m−2 PLA)
2015	LB	26.94	0.002	33.21	0.43	22.87	6.64
	LP	30.39	0.02	36.44	0.46	25.19	8.17
	MB	37.42	4.19	47.74	0.57	29.91	13.03
	MP	40.39	14.51	67.59	0.69	30.99	15.95
	HB	50.06	12.49	56.73	0.71	31.10	15.8
	HP	50.14	28.56	91.90	0.82	31.79	18.48
	VHB	56.81	19.45	66.28	0.78	31.42	17.36
	VHP	55.73	40.03	116.55	1.07	32.16	20.55
2016	LP1	23.42	0.01	31.44	0.36	17.66	3.47
	LP2	22.05	0.03	30.73	0.34	17.11	2.93
	MP	30.41	6.003	55.75	0.54	27.36	11.90
	HP	39.72	21.25	81.04	0.78	31.87	17.39
	VHP	50.14	37.54	113.93	1.04	33.91	20.72

2.2. Effect of O3 on Senescence

Leaf senescence in the two cultivars across different O3 treatments was compared using chlorophyll content index (CCI) as a proxy measure for senescence. The growing seasons were divided into six segments or ‘bins’ of equal thermal time, and statistical analysis testing for the impact of O3 was carried out within each time-bin, in order to identify the point in the season when O3 effects on CCI can be seen to occur. CCI declined over the course of the growing seasons in both cultivars and in both years, and O3 accelerated this senescence (Figure 1). A substantial difference in senescence response of the two cultivars was observed. In the 2015 experiment, cv. Mulika exhibited O3-induced early senescence only in the highest treatment (VHP) (Figure 1A), whereas for cv. Skyfall in the same year, all treatments exhibited accelerated senescence relative to the lowest treatment (Figure 1B). In 2016, for Skyfall the three highest O3 treatments exhibited accelerated senescence (Figure 1C). A statistical summary of this analysis is reported in Table A3 of the Appendix A.

Figure 1

Average relative chlorophyll content index (CCI) in flag leaves for six thermal time groups, (A) cv. Mulika in 2015, (B) cv. Skyfall in 2015 and (C) cv. Skyfall in 2016. Time-bins where a statistically significant association between CCI and accumulated POD0SPEC was observed are marked with an asterisk (*). The direction of O3 effect, i.e., positive (+ve) or negative (−ve) effect on CCI, is also shown. Ozone treatments which exhibited a significant early decline (defined as a decline in CCI of 10% or more relative to the control treatments) are marked in the figure keys with asterisks (*), and those which showed no effect are marked as n.e. (no effect).

Table A3

Outcome of linear mixed models (LMM) analysis investigating which O3 treatments in 2015 and 2016 exhibited an early decline in leaf chlorophyll (CCI) relative to the control treatments.

Cultivar–Year Combination	Ozone Treatment Code	Accelerated Senescence in Treatment Relative to Control? **	Mean POD0SPEC*Thermal Time Interaction Variable p-Value
Mulika, 2015	LP	No	n.e.
	MB	No	n.e.
	MP	No	n.e.
	HB	No	n.e.
	HP	No	n.e.
	VHB	No	n.e.
	VHP	Yes	p < 0.001
Skyfall, 2015	LP	Yes	p < 0.001
	MB	Yes	p < 0.0001
	MP	Yes	p < 0.0001
	HB	Yes	p < 0.0001
	HP	Yes	p < 0.0001
	VHB	Yes	p < 0.0001
	VHP	Yes	p < 0.0001
Skyfall, 2016	LP1	No	n.e.
	MP	Yes	p < 0.0001
	HP	Yes	p < 0.05
	VHP	Yes	p < 0.0001

** Control treatment was defined as the lowest treatment in terms of mean daily POD0SPEC for each season (i.e., treatment LB in 2015 and treatment LP2 in 2016).; n.e. = no significant effect at p < 0.05.

Analysis conducted within different thermal time groupings indicated that a significant negative effect of O3 on CCI was observed substantially earlier in the season for Skyfall compared to Mulika. For Skyfall in 2015, accumulated POD0SPEC was significantly negatively associated with flag leaf CCI from the third thermal time group onwards (1109–1337 °C days), after 25–36 days of O3 exposure. For Mulika in 2015, accumulated POD0SPEC was significantly negatively associated with CCI only at the fifth thermal time group (1568–1796 °C days), after 49–59 days of O3 exposure. The limited CCI data for Skyfall in 2016 supports the 2015 results, with a significant negative association between accumulated POD0SPEC and CCI observed in the 3rd and 4th thermal time-bins. A significant positive association between accumulated POD0SPEC and CCI was observed for Skyfall in 2015, in thermal time group one (649–879 °C days) spanning the first 14 days of O3 exposure. A statistical summary of the time-bin analysis is presented in Table A4 of the Appendix A.

Table A4

Statistics summarising the effect of accumulated POD0SPEC on flag leaf chlorophyll (CCI) in the six thermal time groups.

Cultivar–Year Combination	Thermal Time-Bin (°Cdays)	Number of Observations	POD0SPEC Variable p-Value	Direction of Effect
Mulika, 2015	649–879	44	n.e.
	880–1108	18	n.e.
	1109–1337	16	n.e.
	1338–1567	34	n.e.
	1568–1796	52	p < 0.01	−
	1797–2026	17	n.e.
Skyfall, 2015	649–879	36	p < 0.05	+ve
	880–1108	16	n.e.
	1109–1337	16	p < 0.05	−
	1338–1567	38	p < 0.05	−
	1568–1796	38	p < 0.001	−
	1797–2026	16	n.e.
Skyfall, 2016	649–879	0
	880–1108	0
	1109–1337	50	p < 0.001	−
	1338–1567	55	p < 0.01	−
	1568–1796	0
	1797–2026	0

n.e. = no significant effect at p < 0.05.

The timing of leaf senescence onset and completion was influenced by O3 exposure. For both cultivars in 2015, leaf senescence onset occurred earlier in O3 treatments with higher mean daily POD0SPEC, although this trend was only statistically significant for Mulika (Figure 2A). On average, O3-induced senescence onset occurred later in the season for Mulika (1725 °C days) compared to Skyfall (1216 °C days). Senescence completion also occurred earlier in O3 treatments with a higher mean daily POD0SPEC in both cultivars, according to a linear relationship (Figure 2B). Senescence completion occurred at a similar thermal time on average for both cultivars (Mulika = 1841 °C days, Skyfall = 1867 °C days). The total duration of the O3-induced senescence period was therefore longer for Skyfall than Mulika.

Figure 2

Effect of O3 on the onset and completion of leaf senescence in 2015. (A) Thermal time at senescence onset versus the mean daily POD0SPEC in each treatment. (B) Thermal time at senescence completion versus mean daily POD0SPEC in each treatment. (C) Area under the post-anthesis section of the CCI-thermal time curve versus the mean daily POD0SPEC in each treatment. Solid trend lines indicate a significant regression (p < 0.05), dashed lines indicate that the trend was not significant.

Skyfall also exhibited a linear reduction in the CCI-thermal time curve integral as the mean daily POD0SPEC increased (Figure 2C). This indicates that Skyfall exhibited reduced CCI in the flag leaf throughout the post-anthesis period in elevated O3. No significant association between mean exposure and curve integral was found for Mulika, although the highest treatment in terms of mean flux exposure (VHP) did exhibit a reduced integral compared to the other treatments.

2.3. Ozone Flux at Onset of Early Senescence

Leaf senescence data was also analysed to test whether a particular quantity of O3 flux, or ‘trigger’ flux value, could be identified as inducing the onset of leaf senescence. Accumulated O3 flux at leaf senescence onset for all 2015 treatments which exhibited significant accelerated senescence is shown in Table 2. In the highest O3 treatment (VHP), senescence onset occurred at a substantially lower accumulated POD0SPEC for the cultivar Skyfall (25.7 mmol m−2) compared with Mulika (30.1 mmol m−2), indicating differential sensitivity across the cultivars. When accumulated POD0SPEC and POD6SPEC at senescence onset are compared across the different O3 treatments for the cultivar Skyfall, senescence onset was observed to occur across a fairly wide range of accumulated flux (15.3–25.7 mmol m−2 POD0SPEC, 6.5–18.6 mmol m−2 POD6SPEC). The range of flux at senescence onset was more narrow when flux was calculated without a detoxification threshold (POD0SPEC flux range = 10.4 mmol m−2, POD6SPEC flux range = 12.1 mmol m−2) and considerably more narrow when only the five highest O3 treatments which exhibited the strongest accelerated senescence response are considered (POD0SPEC flux range of five highest treatments = 3.7 mmol m−2, POD6SPEC flux range of five highest treatments = 5.7 mmol m−2).

Table 2

Accumulated flux (PODYSPEC) at the onset of O3-induced senescence for 2015 treatments which exhibited significant accelerated senescence.

Cultivar	O3 Treatment (2015)	POD0SPEC at Senescence Onset (mmol m−2)	POD6SPEC at Senescence Onset (mmol m−2)
Skyfall	LP	17.8	6.5
	MB	15.3	7.9
	MP	22.0	12.9
	HB	24.7	14.1
	HP	25.1	16.2
	VHB	22.9	14.4
	VHP	25.7	18.6
Mulika	VHP	30.1	20.6

2.4. Response of Photosynthesis and gsto Over Time and in Elevated O3

Figure 3 presents combined datasets for four leaf-level physiological parameters capable of short-term or ‘instantaneous’ change in response to environmental stimuli: A (Light-saturated photosynthetic rate), V (maximum carboxylation capacity of rubisco), J (maximum rate of electron transport) and g (stomatal conductance). Data from flag leaf measurements have been combined across both cultivars and both experiments, and the hue of each data point corresponds to the POD0SPEC that had accumulated at the time of measurement (an equivalent figure indicating accumulated POD6SPEC at time of measurement is shown as Figure A3 in the Appendix A). The average physiological values for high and low O3-treated plants within each time-bin are also shown on the plots. The average ‘low’ value represents the mean value for the lowest 2015 treatment (LB) and lowest 2016 treatment (LP2) combined. The average ‘high’ value represents the mean value for the highest 2015 treatment (VHP) and the highest 2016 treatment (VHP) combined. A decline in the photosynthetic parameters (A, V, J) was observed across the growing season, and this decline was accelerated in high O3. g did not decline over time in low O3 but did decline throughout the season in high O3.

Figure 3

Flag leaf data for (A) A, (B) V, (C) J, and (D) g, combined across all cultivar–year combinations. The hue of each data point corresponds to the accumulated POD0SPEC at the moment of measurement. Mean values of physiological parameters in low O3-treated plants (averaged across 2015 LB and 2015 LP2 treatments) and high O3-treated plants (averaged across 2015 VHP and 2016 VHP treatments) are shown as black data points on the plots.

Figure A3

Flag leaf data for (A) Asat, (B) Vcmax, (C) Jmax, and (D) gsto, combined across all cultivar–year combinations. The hue of each data point corresponds to the accumulated POD6SPEC at the moment of measurement. Vertical lines on the plots indicate the divisions between thermal time groups. Mean values of physiological parameters in low O3-treated plants (averaged across 2015 LB and 2015 LP2 treatments) and high O3-treated plants (averaged across 2015 VHP and 2016 VHP treatments) are shown as black data points on the plots.

The outcome of statistical analyses carried out on each cultivar–year combination and in each thermal time group, for the parameters A, V, J and g as measured in the flag leaves, is shown in Figure 4. A statistical summary of this analysis is presented in Table A5, Table A6, Table A7 and Table A8 of the Appendix A. Grey regions on plots denote the period following the observation of a significant negative effect of accumulated POD0SPEC on flag leaf CCI. Across all cultivar–year combinations, no significant negative effects of accumulated POD0SPEC on any of the instantaneous physiological parameters was observed before negative effects of accumulated POD0SPEC on CCI were observed. A significant negative association of accumulated POD0SPEC on the parameters V and J was not observed until the 5th thermal time-bin (1568–1796 °C days).

Figure 4

Plots showing the response of A, V, J, and g to O3 flux. The y-axis represents the accumulated POD0SPEC-physiology slope in the ‘best’ linear mixed regression model for each thermal time group. Positive slope indicates a positive effect of O3 on the physiological variable; a negative slope indicates a negative effect. (A) cv. Mulika in 2015, (B) cv. Skyfall in 2015, (C) cv. Skyfall in 2016. Coloured symbols indicate a significant POD0SPEC–physiology association; black symbols indicate no statistically significant physiological response. Grey regions on plots indicate the period following an observed significant effect of O3 on flag leaf CCI.

Table A5

Statistics summarising the effect of accumulated POD0SPEC on flag leaf Asat.

Cultivar–Year Combination	Thermal Time-Bin (°Cdays)	Number of Observations	POD0SPEC Variable p-Value	Direction of Effect
Mulika, 2015	649–879	9	n.e.
	880–1108	6	n.e.
	1109–1337	0
	1338–1567	0
	1568–1796	16	p < 0.05	−
	1797–2026	0
Skyfall, 2015	649–879	8	n.e.
	880–1108	8	p < 0.01	+ve
	1109–1337	0
	1338–1567	0
	1568–1796	7	p < 0.05	−
	1797–2026	8	p < 0.001	−
Skyfall, 2016	649–879	0
	880–1108	52	p < 0.001	+ve
	1109–1337	28	p < 0.05	−
	1338–1567	118	p < 0.05	−
	1568–1796	56	p < 0.001
	1797–2026	66	p < 0.001

n.e. = no significant effect at p < 0.05.

Table A6

Statistics summarising the effect of accumulated POD0SPEC on flag leaf Vcmax.

Cultivar–Year Combination	Thermal Time-Bin (°Cdays)	Number of Observations	POD0SPEC Variable p-Value	Direction of Effect
Mulika, 2015	649–879	9	n.e.
	880–1108	6	n.e.
	1109–1337	0
	1338–1567	0
	1568–1796	16	p < 0.01	−
	1797–2026	0
Skyfall, 2015	649–879	8	n.e.
	880–1108	8	n.e.
	1109–1337	0
	1338–1567	0
	1568–1796	7	n.e.
	1797–2026	8	p < 0.001	−
Skyfall, 2016	649–879	0
	880–1108	54	p < 0.001	+ve
	1109–1337	28	n.e.
	1338–1567	118	n.e.
	1568–1796	56	p < 0.001	−
	1797–2026	66	p < 0.001	−

n.e. = no significant effect at p < 0.05.

Table A7

Statistics summarising the effect of accumulated POD0SPEC on flag leaf Jmax.

Cultivar–Year Combination	Thermal Time-Bin (°Cdays)	Number of Observations	POD0SPEC Variable p-Value	Direction of Effect
Mulika, 2015	649–879	9	p < 0.05	+ve
	880–1108	6	n.e.
	1109–1337	0
	1338–1567	0
	1568–1796	16	p < 0.01	−
	1797–2026	0
Skyfall, 2015	649–879	8	n.e.
	880–1108	8	n.e.
	1109–1337	0
	1338–1567	0
	1568–1796	7	p < 0.05	−
	1797–2026	8	p < 0.001	−
Skyfall, 2016	649–879	0
	880–1108	16	n.e.
	1109–1337	28	n.e.
	1338–1567	48	n.e.
	1568–1796	0
	1797–2026	14	n.e.

n.e. = no significant effect at p < 0.05.

Table A8

Statistics summarising the effect of accumulated POD0SPEC on flag leaf gsto.

Cultivar–Year Combination	Thermal Time-Bin (°Cdays)	Number of Observations	POD0SPEC Variable p-Value	Direction of Effect
Mulika, 2015	649–879	9	n.e.
	880–1108	6	n.e.
	1109–1337	0
	1338–1567	0
	1568–1796	15	p < 0.05	−
	1797–2026	0
Skyfall, 2015	649–879	8	n.e.
	880–1108	59	p < 0.01	+ve
	1109–1337	28
	1338–1567	116
	1568–1796	63	n.e.
	1797–2026	74	p < 0.001	−
Skyfall, 2016	649–879	8
	880–1108	7	n.e.
	1109–1337	0	n.e.
	1338–1567	0	p < 0.01	−
	1568–1796	7	p < 0.001	−
	1797–2026	8	p < 0.001	−

n.e. = no significant effect at p < 0.05.

Some evidence of heightened physiological performance in the early-season in high O3 was observed across both cultivars and both years, although the pattern was not consistent. A significant positive association between accumulated POD0SPEC and physiology in either the first or second time-bin was observed for i) J in Mulika in 2015, ii) A and g in Skyfall in 2015, and iii) A and V in Skyfall in 2016.

2.5. Comparison of O3 Exposure Metrics for Predicting Physiological Response to O3

For all measured physiological parameters (CCI, A, V, J, g) and for both cultivars, a flux-based metric of exposure was better at predicting physiological response of wheat to O3, compared to the concentration-based metrics (24-h mean and AOT40) (Table 3). For four out of the ten model sets created in this analysis, the accumulated POD0SPEC (i.e., without a threshold for accumulation) produced the best model fit. For the other six model sets, the accumulated POD0SPEC and POD6SPEC metrics were equally good at predicting physiological response. The O3 flux metric with no threshold for accumulation was therefore equal to, or better than, the O3 flux metric with a detoxification threshold at predicting the physiological response to O3.

Table 3

Summary of linear mixed model regression analysis to determine whether accumulated POD0SPEC, accumulated POD6SPEC, 24-h mean, or AOT40 represent the best predictor of physiology in the 5th thermal time-bin. The lowest AIC (Akaike Information Criterion) for each parameter and cultivar, indicating the best model, is highlighted in grey. The outcome of model selection to determine if the profile of O3 exposure (i.e., peak vs. background) was important in the flux-based models is also shown.

Parameter	Cultivar	AIC: POD0SPEC	AIC: POD6SPEC	AIC: AOT40	AIC: 24-h Mean	O3 Profile Important in POD0SPEC Model?	O3 Profile Important in POD6SPEC Model?
CCI	Mulika	−62.8	−60.5	−59.2	−58.3	No	No
	Skyfall	3.7	10.4	12.3	15.0	No	No
Asat	Mulika	−0.5	1.3	4.1	3.2	No	No
	Skyfall	−82.9	−68.0	−48.1	−54.2	No	No
Vcmax	Mulika	6.1	6.9	9.9	9.0	No	No
	Skyfall	−63.0	−62.9	−34.6	−32.3	No	No
Jmax	Mulika	−1.7	0.2	3.2	2.3	No	No
	Skyfall	6.3	7.3	9.2	9.5	Yes	No
gsto	Mulika	13.4	14.8	17.3	16.5	No	No
	Skyfall	−19.1	−7.3	2.7	−1.6	No	No

The inclusion of an explanatory variable describing the profile (i.e., peak or background) of O3 exposure in the ‘best’ model did not improve fit in nine out of the ten model sets created with the accumulated POD0SPEC metric, and in all models created with the accumulated POD6SPEC metric. Using O3 flux as the metric of exposure therefore accounts for differences in the O3 exposure resulting from peak-dominated treatments and those featuring a consistent background level, in the majority of cases.

3. Discussion

The first aim of the analysis presented here was to assess whether published approaches for modelling O3-induced senescence can account for inter-cultivar variation in response. Both cultivars exhibited accelerated senescence in response to O3, but the pattern of response differed according to cultivar. In 2015, significant accelerated senescence was observed in seven O3 treatments for Skyfall, but only in the highest treatment for Mulika, suggesting higher O3 tolerance in Mulika (Figure 1). This differential tolerance is also indicated by the earlier appearance of significant O3 effects on leaf CCI across all treatments for Skyfall compared to Mulika. Senescence completion occurred progressively earlier, hence total leaf duration became progressively shorter, in both cultivars as average O3 flux (mean daily POD0SPEC) in the treatment increased (Figure 2B), according to a linear relationship. Completion of leaf senescence occurred at a similar thermal time in both cultivars (Mulika = 1841 °C days, Skyfall = 1867 °C days), meaning that the total senescence duration was longer for Skyfall. While O3-induced senescence in Mulika was characterised by a sudden drop in leaf CCI in the late-season, Skyfall exhibited a more gradual O3-induced decline in CCI.

The linear relationship between mean flux and total leaf duration observed in this study for both cultivars gives support to the senescence function of Ewert and Porter [25], which assumes a linear decline in mature leaf lifespan as O3 exposure increases. However, evidence presented here suggests that the duration of leaf senescence is likely to vary with O3 exposure. Additionally, the differential senescence duration in the two cultivars suggests that a key assumption of the Ewert and Porter function—that leaf senescence comprises the final third of the mature leaf lifespan—may not hold true for all cultivars. For example, in 2015 for Skyfall, leaf senescence in the highest O3 treatment comprised 76.7% of the total flag leaf lifespan (flag leaf emergence = 877 °C days, leaf senescence onset = 1075 °C days, senescence completion = 1725 °C days)—substantially greater than one third. The inter-cultivar variation in senescence response observed in this study would therefore only be captured by a model function that allows for the proportion of leaf lifespan comprising leaf senescence to be parameterised according to cultivar and ozone exposure. Our results suggest that to effectively model variation in the pattern of O3-induced senescence, the timing of senescence onset, and the rate (or duration) of senescence, models need to be calibrated for particular species and cultivars.

However, the analysis in this paper highlights the uncertainty associated with the approach of Danielsson et al. [20] for modelling onset of O3-induced senescence using a threshold of accumulated flux. Following this approach, senescence onset may occur at different points in time at different levels of mean exposure but should occur at approximately the same value of accumulated flux. This method was designed in the absence of a known mechanism for induction of senescence by O3 but could be interpreted mechanistically if accumulated O3 flux is assumed to be proportional to increased respiratory effort accumulated over the season, which has been proposed as a potential trigger for O3-induced senescence [25]. For Skyfall, across the five highest O3 treatments in 2015, onset occurred across a POD0SPEC range of 22.0–25.7 mmol m−2. Given the limitations associated with the method used to identify senescence onset—arbitrarily defined as a 10% reduction in leaf CCI relative to the control—as well as the inherent variation that exists between seedlings, this flux range can be considered relatively narrow. However, when all treatments that exhibited a significant O3 effect on senescence are considered for Skyfall in 2015, the range of flux at senescence onset is considerably wider (17.8–25.7 mmol m−2 POD0SPEC). These results provide an estimate of the degree of error potentially associated with applying this approach in models and suggest that the accumulation of respiratory effort does not fully explain the triggering of O3-induced senescence.

A second objective of this study was to test the hypothesis that O3 reduces photosynthetic rate in the short term by reducing carboxylation capacity of rubisco (V). The assumption that O3 reduces V is central to the ‘instantaneous’ O3 effect function of Martin et al. [27]. A version of this function is also applied by Ewert and Porter [25], where O3 reduces photosynthesis in a short-term and reversible way, in addition to and independent of the O3 senescence effect. The analysis presented here found that a significant negative effect of O3 on photosynthesis and g was only observed concurrent with O3-induced leaf senescence (Figure 4). This result was consistent across all cultivar–year combinations. For Skyfall in 2015 and 2016, significant negative effects of O3 on A were observed before a negative association between O3 and V, suggesting that reduced carboxylation capacity is not responsible for the initial reduction in photosynthetic capacity observed in these experiments. There was therefore no evidence of an ‘instantaneous’ effect of O3 on the photosynthetic mechanism, in the period preceding leaf senescence.

These results contradict several studies which observed short-term reduction in photosynthetic rate in response to O3 [39,40,41]. One possible explanation for this contradiction is that instantaneous reduction of carboxylation capacity by O3 is only relevant at acute concentrations. The reduced carboxylation efficiency reported by Farage et al. [41] was observed following 4–16 h of exposure at unrealistically high O3 concentrations of 200–400 ppb, considerably higher hourly concentrations compared to those used in the experiments in this study, which more closely mimic ambient conditions (maximum hourly O3 exposure of 117 ppb). The results presented here therefore indicate that accelerated senescence is likely to be more important than direct effects on photosynthesis for determining crop yield loss in most agricultural landscapes, where O3 concentrations are typically moderate (ranging from ~20 to ~45 ppb) with occasional peaks in concentration [5,52]. Understanding and simulating the early senescence response to O3 should therefore be the priority for O3 experimentalists and modellers.

Alternatively, our results could be explained by a differential response to O3 in younger and older leaves. Bernacchi et al. [53] and Morgan et al. [54] observed in field experiments with soybean that O3 effects on photosynthesis and g were not apparent in new fully expanded leaves and Reichenauer et al. [55] saw similar results in three wheat cultivars. Younger leaves may have a higher tolerance to O3, or alternatively the O3 effect on photosynthesis may be associated with a cumulative build-up of O3 damage in leaves or leaf age. Either way, the age-dependency of O3 effects is an important consideration in O3 effects modelling. The function described by Ewert and Porter [25] for modelling short-term effects of O3 on photosynthesis allows for leaf age to influence the rate of overnight recovery from O3 damage but not the threshold for damage. The role of leaf age in determining O3 flux thresholds would benefit from further investigation.

A surprising result from the data analysis is that O3 had a significant positive effect on some physiological parameters early in the season. CCI, A, V, and g all exhibited a positive association with O3 exposure for one or more of the cultivar–year combinations, in either the first or second thermal time group (up to 32 days following beginning of exposure). Stimulation of growth at low doses of a toxin is known as hormesis, and this phenomenon has been observed previously in a number of plants in response to exposure to low concentrations of O3 [56]. Stimulation of photosynthesis and g in wheat during the first few weeks of O3 exposure was also observed by Mulholland et al. [57] in their open-top chamber experiment, although generally there are few reported cases of this phenomenon in crop species. Ozone-induced physiological stimulation could be an adaptive response associated with plant defence responses, for example, heightened g and photosynthesis may enable the upregulation of antioxidant synthesis. Observations that yield can be stimulated at low O3 exposure concentrations have also given rise to the theory that free radicals, at low concentrations, can act as growth promotors in plants [58]. An alternative hypothesis is that the observed early-season physiological boost in this study is related to disruption of stomatal control by O3, as has been observed in some grassland species [59,60], leading to heightened g and an associated boost in other physiological parameters. Our results suggest that a bi-phase hermetic dose-response curve, such as proposed by Agathokleous et al. [56], could be considered in ozone effect functions. More experimental data is needed in order to establish if the early-season physiological boost induced by O3 in this study is consistent across other plant species and environments.

The third aim of this study was to test whether O3 flux would be a better predictor of physiological response than concentration-based exposure (AOT40 and 24-h mean). Flux was superior at predicting the response to O3 of five physiological variables (CCI, A, V, J, g) in regression models, for both cultivars (Table 3). Previous studies have reported that flux is better than AOT40 at predicting the spatial distribution of O3 injury [61] and is a better predictor of the response of wheat yield, poplar biomass, and assimilation rate in urban trees to O3 exposure [62,63,64]. However, few have compared the association between leaf-level physiology and different exposure metrics in crop species. Our results align with the general consensus in the O3 research community that O3 flux represents a more biologically relevant metric of O3 exposure than ambient concentration [61,65,66,67] and indicate that O3 flux should be the preferred metric of exposure in O3 effect model functions. More surprising is the fact that the flux metric without an accumulation threshold, POD0SPEC, was a better or equal predictor of physiological response compared to POD6SPEC, which employs an accumulation threshold of six. POD6SPEC produced the closest correlation between flux and yield of wheat in a previous analysis testing different flux accumulation thresholds [50] and has been applied in several assessments of O3 impacts in wheat [22,68]. More research is therefore needed to establish how much the capacity to detoxify O3 varies between cultivars and why the threshold flux required to induce leaf-level physiological changes appears to differ from the threshold required to reduce yield. This result also highlights the limitations of an empirical rather than a mechanistic approach to modelling ozone tolerance and detoxification capacity. The development and integration of mechanistic ozone detoxification modules in models, such as that proposed by Plöchl et al., could potentially lead to improved model accuracy [69].

The view that O3 flux should be the metric of exposure in O3 effect modelling is also supported by the fact that O3 flux accounted for the different levels of exposure in treatments dominated by peaks in concentration, versus those characterised by a consistent background level, for the majority of physiological parameters. Flux is therefore likely to perform well as a predictor of physiological response across different world regions which are currently experiencing divergent trends in the pattern of O3 exposure [10,11,12].

Limitations of this study need to be considered when interpreting and applying the results. Calculated values of O3 flux were not verified by leaf-level measurements of gas flux through stomata. However, the decision to apply the multiplicative DO3SE model in this study was based on the fact that fluxes produced by this model have previously been evaluated in several independent studies which have demonstrated the model’s predictive capability [19,21,22,70,71,72]. A further limitation is that estimates of the onset of leaf senescence were based on leaf chlorophyll content, which would have represented both the chlorophyll loss resulting from leaf injury, as well as chlorophyll loss relating to senescence. In addition, the analysis is based on only one crop species and two cultivars. As the variation in yield response to O3 exhibited by different crops, and different cultivars within the same crop species, is well established [2,73], caution must be used when extrapolating results presented here to other wheat cultivars and crops. It should however be noted that the observation in this study that no O3 effect on photosynthesis could be observed in young wheat leaves—indicating the senescence response is more important than direct effects on photosynthesis—is supported by previous work in other wheat cultivars [55] and by other experimental work in soybean [53,54]. Considerations when applying the results presented in this study, particularly when attempting to up-scale modelled responses from the leaf to canopy level, include the fact that the response observed in the wheat flag leaf may differ from the responses of lower-canopy leaves, and exposure to O3 during early seedling and leaf development may also alter the sensitivity to O3 observed in the flag leaf.

In conclusion, this study has shown that current approaches for modelling O3 effects on leaf longevity and photosynthesis in crops have some limitations and are not fully supported by the experimental data presented here. Model functions representing O3-induced senescence must allow for parameterisation of the timing of senescence onset and rate of senescence, if inter-cultivar variation in response is to be accurately simulated. Further research aimed at understanding the mechanistic ‘trigger’ of O3-induced senescence should be a priority, as this understanding may allow for the development of a more effective mechanism in models for inducing the senescence response. The results also suggest an age-dependency in the response of photosynthesis to O3—which is not currently fully considered in modelling methods—and indicate that acceleration of senescence is more important than direct effects of O3 on photosynthesis in determining final O3-induced yield loss, at O3 concentrations that crops are likely to be exposed to on a day-to-day basis. Building functions that can accurately represent the O3-induced senescence effect in crops should therefore be the priority for O3 effect modellers.

4. Materials and Methods

4.1. Experimental Site and Treatments

Both experiments that provided data for this study took place at the CEH air pollution exposure facility in Abergwyngregyn, North Wales (53.2°N, 4.0°W). The two experiments were conducted independently and differ in terms of the number of cultivars and treatments, however, the experimental design was similar across both years and the same protocols for physiological measurements were applied. Data was therefore subsequently pooled across the two years to produce a more robust dataset. In 2015, two European wheat cultivars (Triticum aestivum L., ‘Mulika’ and ‘Skyfall’) were exposed to O3 for 82 days. In 2016, only the more O3-sensitive cultivar, ‘Skyfall’, was used and was exposed for 92 days. Timelines for sowing, emergence, O3 exposure, and harvest in both experiments are presented in Figure A1 of the Appendix A. In both experiments, plants were grown in 25-litre containers (40 cm × 35 cm × 38 cm) filled with John Innes No.3 compost (Westland Horticulture Limited, Huntingdon, UK), and soil was inoculated shortly after sowing with microbial communities using a soil slurry taken from a nearby wheat field. Seeds were sown in rows 7 cm apart at a density of ~260 seedlings per square metre, which aligns with recommended seedling density for field conditions [74]. Four containers per cultivar–treatment combination were planted and placed alongside each other, producing a canopy of ~144 plants per cultivar–treatment combination. In both years, ammonium nitrate fertiliser was applied once mid-season (80 kg/ha). In 2015, fungicide (‘Unix’, cyprodinil, 1.6 kg/ha—Syngenta, Bracknell, UK) was applied once and insecticide (pyrethrum, 1 mL/L—Bayer, Monheim Am Rhein, Germany) applied three times. In 2016, fungicide was applied twice (1st application: Trifloxystrobin, 0.12 g/L—Bayer, Mondheim Am Rhein, Gernany; tebuconazole, 0.125 g/L—Bayer, Mondheim Am Rhein, Germany; 2nd application: cyprodinil, 2.25 g/L—Syngenta, Bracknell, UK) and insecticide applied once (thiachloprid, 0.15 g/L #x2014;Bayer, Mondheim Am Rhein, Germany).

Figure A1

Timeline for sowing, seedling emergence, O3 exposure, and plant harvest for (A) the 2015 experiment and (B) the 2016 experiment.

Ozone exposure took place within ‘solar domes’, hemispherical glass domes three metres in diameter and two metres in height, described previously [59,75]. Air entering domes was carbon-filtered to remove O3 and other air pollutants (i.e., NOX and SO2), and the solar dome site is an area with some of the lowest NOX and SO2 emissions in the country and is not downwind of any major sources [76], minimising any other pollutants that could have entered the domes after filtration. Following filtration, a precision-controlled quantity supplied by an O3 generator (Dryden Aqua G11, Edinburgh, UK) linked to an oxygen concentrator (Sequal 10, Pure O2, Manchester, UK), was added to incoming air. Injection concentrations were determined by a computer-controlled O3 injection system (Lab VIEW, version 8.6, National Instruments, Austin, TX, USA). Air within domes was circulated at a rate of two air changes per minute, and the O3 concentration within each dome was recorded on a 30-min cycle using two O3 analysers of matched calibration (Envirotech API 400A, St Albans, UK). Exposure profiles for each treatment are presented in Figure 5. Treatments spanned a range of seasonal mean concentrations and represented different O3 exposure patterns, representing potential future profiles of increasing background or decreasing peak O3. In 2015, four treatments consisted of a low night-time background level, with high peaks during the day, classified as ‘peak’ treatments, while the other four treatments comprised of consistent concentrations with only small peaks, classified as ‘background’ treatments. In 2016, all treatments were ‘peak’ in profile. Although O3 treatments were not replicated, numerous studies have established the statistical validity of conducting unreplicated experiments using the solar dome facility [56,77,78], and previous work has shown that no solar dome effect on air or leaf temperature is detectable [79].

Figure 5

Average hourly O3 exposure concentrations in (A) 2015 and (B) 2016. Values are shown for a one-week period, averaged over the entirety of each growing season. Each treatment has been categorised based on the 24-hour mean exposure (Low, Medium, High, Very High) and the characteristic profile of exposure (peak or background). Treatments were applied five days out of seven to mimic real-world O3 exposure.

Climatic conditions fluctuated naturally in the solar domes according to ambient conditions. Air temperature, photosynthetically active radiation (PAR), relative humidity, and wind speed were monitored in one solar dome during both experiments using an automatic weather station (Skye instruments Ltd, Llandridod Wells, UK) to obtain data for stomatal flux modelling. Plants were well-watered throughout, and soil moisture content was continuously monitored in selected plant containers to a depth of 10 cm using Theta Probes (Delta-T Devices Ltd., Cambridge, UK).

4.2. Leaf Chlorophyll and Gas Exchange Measurements

Chlorophyll content was measured non-destructively as an index (chlorophyll content index, CCI) using CCM-200 and CCM-200+ instruments (Opti-sciences, Hudson, NH, USA). A regression line fit to paired measurements was used to standardise observations made using the two instruments. In 2015, 684 measurements were made over 70 days; in 2016, 105 measurements were made over 22 days.

To assess the effect of O3 on photosynthetic capacity and g, response curves of net photosynthetic rate (A) to intercellular CO2 concentration (Ci), i.e., A–Ci curves, were constructed using a portable infrared gas analyser (Li-Cor 6400XT; LI-COR Biosciences, Lincoln, NE, USA). In 2015, measurements were made in the two lowest O3 treatments at the beginning of exposure (20–26 May). Further measurements in the two lowest treatments (LB and LP) and two high O3 treatments (VHB and VHP) were made in the mid-season (8–17 June) and late-season (16–24 July). Measurements were made in the youngest fully expanded leaf of randomly selected plants (represented by the flag leaf) from 28th May onwards. In 2016, four sets of A–Ci curve measurements were made at approximate two-weekly intervals spanning 6 June to 29 July. Measurements in 2016 were made in all treatments at each of the time intervals, except for the final measurement set in late July, when plants in HP and VHP treatments were too senesced for measurements to take place. All 2016 measurements were made in the flag leaf. For both years, four A–Ci measurements were made per treatment and per cultivar at each timepoint, and leaves were tagged following measurement so that the same leaf could be measured throughout the season.

All response curve measurements were conducted at a selected light saturation (minimum photosynthetic photon flux density = 1500 µmol m−2 s−1, LED light source), and sample chamber relative humidity was maintained between 50 and 80%. Photosynthetic rate and g were allowed to stabilise in the leaf chamber at ambient CO2 (400 µmol mol−1). The A–Ci curve was constructed by measuring A at a minimum of nine air CO2 concentrations, ranging from ca. 50 to 2000 µmol mol−1. Asat and associated g values were determined from the ambient CO2 measurements (400 µmol mol−1) from each A–Ci curve.

Additional measurements of Asat and associated g were made in 2016 over six days (16 June, 1 July, 8 July, 14 July, 20 July, 26 July). Measurements were made at ambient CO2 concentration (400 µmol mol−1) under the same light and relative humidity conditions as described above.

4.3. Derivation of Vcmax and Jmax

Maximum rate of carboxylation (V) and maximum rate of electron transport (J) were derived from A–Ci curves using the estimating utility and methodology described by Sharkey et al. [80]. Leaf temperature and atmospheric pressure were the input parameters, which were measured using the Licor 6400XT simultaneously with all photosynthesis measurements. V and J values calculated from curves were adjusted to 25 °C.

The V dataset was extended by applying the ‘one-point method’ of deriving V from A as described by De Kauwe et al. [81]. Estimation of V when only A is known using the one-point method relies on the assumption that photosynthetic rate at ambient CO2 is rubisco-limited [81]. As the measurements of A at 400 µmol mol−1 CO2 in the measured A–Ci curves typically fell within the rubisco-limited section of the curve (i.e., before the transition point), this assumption was thought to be likely to hold true for the two cultivars used in this study. The one-point method also assumes, in the absence of a known daytime respiration rate (R), that R can be estimated as 1.5% of V. V was calculated from A using the following equation: where K is the Michaelis–Menten constant, given by:

The parameters K (Michaelis–Menten constant for CO2), K (Michaelis–Menten constant for O2) and Γ* (CO2 compensation point in the absence of mitochondrial respiration) were estimated at 25 °C following the equations and constants published by Bernacchi et al. [82] describing their temperature dependence in the model species tobacco (Nicotiana tabacum, L.). Equations and constants used to derive these three parameters are listed in Table A2 of the Appendix A. O represents the intercellular concentration of O2 (210 mmol mol−1) [81].

Table A2

Equations used to derive the Michaelis–Menten constant for CO2 (KC), the Michaelis–Menten constant for O2 (KO), and the CO2 compensation point in the absence of mitochondrial respiration (Γ*).

Parameter	Unit	Equation Used in Derivation
KC	µmol mol−1	404.9∗exp(79403(T_K−298.15)298.15∗R∗T_K)
KO	mmol mol−1	278.4∗exp(36380(T_K−298.15)298.15∗R∗T_K)
Γ*.	µmol mol−1	42.75∗exp(37830(T_K−298.15)298.15∗R∗T_K)

T = leaf temperature in Kelvin; R = universal gas constant, 8.314 J mol−1 K−1.

The robustness of the one-point method was evaluated by comparing V values calculated from a subset of the measured A–Ci curves with the V values calculated from each corresponding A value (i.e., the 400 µmol mol−1 CO2 value from each A–Ci curve). V values derived using both methods were adjusted to 25 °C. A very close association was observed between V values derived using the two methods (Figure 6, adjusted r = 0.95, p < 0.001), indicating that the one-point method is robust for the cultivars used in this study. V values derived using the one-point method were therefore pooled with A–Ci-derived V values for analysis, and the potential error introduced through the use of two different derivation methods was accounted for in the statistical analysis by including—in model selection—an explanatory variable describing the derivation method (explained in more detail in Section 4.7.3).

Figure 6

Plot of V values derived from A–Ci curves versus V values calculated using the one-point method [75] from the corresponding A value extracted from each curve (A at 400 µmol mol−1 CO2). The blue line represents the linear regression model fit (p < 0.001, adjusted r = 0.95; line equation, y = 0.99x + 1.33). The red dashed line represents the line of x = y. Data for this comparison comprise a subset of the A–Ci curve dataset used in this study.

4.4. Modelling O3 Flux

Stomatal O3 flux to the flag leaf was modelled in each treatment for both years to derive a measure of exposure that accounted for the environmental influence on O3 uptake and could be tracked over time. Flux was modelled using the multiplicative g module of the DO3SE model [18], which has a published parameterisation for European wheat [49,50,83] and has been applied previously to model O3 flux to this crop [22,84]. A summary of the DO3SE algorithms and parameters used in this study are presented in Section A1 of the Appendix A.

Ozone flux for wheat is accumulated above a detoxification threshold of six in the DO3SE methodology (producing the POD6SPEC flux metric, species-specific phytotoxic O3 dose above a threshold of 6 mmol m−2 PLA s−1—previously known as the POD6, with “SPEC” referring to the species-specific version of the DO3SE model) [83], as this threshold has produced the closest correlation between POD and wheat yield in previous experiments [50]. However, as thresholds of physiological effect in wheat have been far less studied, the POD0SPEC (where no threshold for accumulation is applied, previously known as the POD0) was also calculated, in order to avoid assuming a threshold of effect. Modelled POD0SPEC and POD6SPEC for 2015 and 2016 O3 treatments are shown in Figure 7.

Figure 7

Modelled O3 flux over time in the different O3 treatments. (A) POD0SPEC in 2015, (B) POD6SPEC in 2015, (C) POD0SPEC in 2016, (D) POD6SPEC in 2016. Each O3 treatment in both years was categorised based on the 24-hour mean exposure (L = low, M = medium, H = high, VH = very high) and the characteristic profile of exposure (P = peak, B = background).

4.5. Alignment of Physiological Observations with O3 Flux, and Calculation of Mean Flux Exposure (Mean Daily POD0SPEC)

Each physiological observation (CCI, A) was aligned with the treatment-specific accumulated POD0SPEC and POD6SPEC on the day of measurement and at the exact time of measurement, wherever this data was available (referred to hereafter in this paper as ‘accumulated POD0SPEC’ and ‘accumulated POD6SPEC’). This was done to allow the impact of real-time O3 flux exposure on physiology to be analysed. The mean daily POD0SPEC (i.e., the average accumulation of flux per day, mmol m−2 PLA day−1,) was also calculated for each O3 treatment to act as a metric of mean exposure intensity. Mean daily POD0SPEC values for each O3 treatment are presented in the results section in Table 1. Mean daily POD0SPEC was calculated as the average of daily POD0SPEC accumulation from A until the modelled onset of senescence.

4.6. Data Standardisation

The two experiments had different sowing and harvest calendars. In order to compare timings across the two experiments, time was therefore standardised by calculating thermal time from plant emergence onwards (daily mean temperature sum >0 °C). Physiological data was also standardised, by conversion from raw to relative values. This was done to account for differences in instrument calibration between years, and to account for differences in beginning of season ‘baseline’ physiology between cultivars. Relative values for the physiological observations were calculated by deriving a reference value for each parameter (CCI, A) and each cultivar–year combination (i.e., Mulika in 2015, Skyfall in 2015, Skyfall in 2016). The reference value—calculated as the 90th percentile value of all observations, spanning the whole season and all treatments—was used as the baseline for calculating relative change. Skyfall CCI data from 2016 comprised too few data points for the derivation of an individual reference value, 2016 and 2015 CCI data for Skyfall was therefore combined to produce a single reference value for Skyfall, as CCI data for Skyfall was found to not significantly differ by year (p = 0.06 in regression model). A comparison of CCI data for Skyfall measured in 2015 and 2016 can be found in Section A3 of the Appendix A.

4.7. Statistical Analysis

Statistical analysis was conducted in R version 3.3.2 [85], and either involved linear regression or linear mixed models (LMMs) using the package lme4 v1.17. Model selection was by AIC (Akaike Information Criterion). The model with the lowest AIC was considered the ‘best’ model of those fitted, and models differing in <2 AIC units from the best model were defined as having little empirical support [86]. Wherever relevant, a random factor describing solar dome number was included in models to account for multiple measurements made within domes, and unique pot ID was a random factor when analysis involved multiple measurements made from the same pot. p-values were obtained for terms in the optimal models using the R package lmerTest, v2.0-33 [87]. Assumptions of normality and even spread of residuals were checked using residual plots and data were transformed where necessary. Four key analyses were conducted as part of this study and are described in more detail below.

4.7.1. Identification of O3 Treatments with Significantly Accelerated Senescence

Flag leaf CCI data was analysed in all O3 treatments to test for accelerated senescence. Each elevated O3 treatment was paired in turn with the control treatment for that experiment, and the significance of the thermal time/mean daily POD0SPEC interaction term was tested using LMMs. Control treatments were defined as the lowest in terms of mean daily POD0SPEC and comprised of treatment LB for the 2015 experiment and treatment LP2 for the 2016 experiment.

4.7.2. Analysis of O3 Effect on the Timing of Senescence Onset and Completion

The impact of O3 on leaf senescence onset and completion was examined using regression models fitted to each of the 2015 O3 treatments (separately for the two cultivars). It was not possible to conduct this analysis for 2016, as 2016 CCI measurements only spanned 22 days. Regression models comprised of relative CCI as the dependent variable and thermal time as the independent variable, and the shape of response was determined by comparing linear, quadratic, and cubic models. The best model for each O3 treatment was then used to determine i) thermal time at leaf senescence onset, ii) thermal time at senescence completion, and iii) the post-anthesis curve integral (i.e., area under the curve), as shown in Figure 8. Thermal time at senescence onset was aligned with the accumulated POD0SPEC at that time for each O3 treatment, to identify the accumulated flux ‘trigger’ values for senescence onset.

Figure 8

Summary of methods used to derive: i) Thermal time at leaf senescence onset, defined as a 10% reduction in relative CCI in the elevated treatment (grey line) relative to the control (black line); ii) thermal time at senescence completion, defined by the x-abscissa of the treatment regression line; iii) the post-anthesis integral of the regression curve, indicated on the plot as shaded regions (Post-anthesis period in 2015 = 1142 °C days onwards). Diagram not to scale.

4.7.3. Analysis of Relative Timing of O3 Effects on Different Aspects of Physiology

The effect of accumulated POD0SPEC on CCI, A, V, J, and g during successive periods of the growing season was analysed, to identify when O3 began to influence physiology. The range of thermal time spanned by flag leaf physiological measurements was divided into six thermal time-bins of equal width. The effect of accumulated POD0SPEC on each parameter, within each time-bin and for each cultivar–year combination, was analysed by comparing model fit with and without accumulated POD0SPEC as an explanatory variable. An additional explanatory variable was included in model selection for V, describing the derivation method (i.e., A–Ci curve or one-point method).

4.7.4. Comparison of Flux and Concentration-Based O3 Exposure Metrics for Predicting Physiological Response

Accumulated POD0SPEC, accumulated POD6SPEC, 24-h mean concentration, and AOT40 (accumulated O3 > 40 ppb during daylight hours) were compared in their ability to predict the response of CCI, A, V, J, and g during the 5th thermal time-bin. The 5th time-bin was selected for this analysis as most physiological parameters exhibited a response to O3 exposure during this time. For each physiological parameter, LMMs constructed with each of the metrics of O3 exposure were compared for model fit. An explanatory variable describing whether O3 had been administered as a ‘peak’ or ‘background’ profile was also included in model selection, to test whether the O3 metric that produced the best model fit also accounted for different patterns of exposure.

Table A1

Parameters applied in DO3SE to derive the species-specific accumulated O3 flux (POD0SPEC and POD6SPEC) for the 2015 and 2016 experiments.

Function	Parameter	Units	Parameter Description	Parameter Value
gmax	gmax	mmol O3 m−2 PLA s−1	maximum rate of gsto_O3	500
fmin	fmin	Fraction	Fraction of gmax at minimum gsto_O3	0.01
fphen	Mid-anthesis	DOY	Decimal growth stage 65*	174 (2015); 181 (2016)
	Astart	DOY	Beginning of flux accumulation/flag leaf emergence	146 (2015); 158 (2016)
	Aend	DOY	End of flux accumulation/flag leaf senescence	207 (2015); 214 (2016)
	fphen_a	Fraction	Proportional fall in gsto between fphen_g and fphen_h	0.3
	fphen_b	Fraction	Fraction of gmax that gsto_O3 takes at the beginning of flag leaf senescence	0.7
	fphen_e	°C days	Temperature sum at Astart	−490
	fphen_f	°C days	Temperature sum at mid-anthesis	0
	fphen_g	°C days	Temperature sum at end of maximum gsto_O3 following mid-anthesis	100
	fphen_h	°C days	Temperature sum at start of flag leaf senescence	525
	fphen_i	°C days	Temperature sum at Aend	795
flight	light_a	Constant	The rate of saturation of gsto in response to photosynthetically active radiation (PAR)	0.0105
ftemp	Tmin	°C	Temperature below Topt where gsto reaches fmin	12
	Topt	°C	Optimum temperature for gsto	26
	Tmax	°C	Temperature above Topt where gsto reached fmin	40
fVPD	VPDmax	kPa	Value where vapour pressure deficit (VPD) begins to limit gsto	1.2
	VPDmin	kPa	Value of VPD where fmin is reached	3.2
	ΣVPDcrit	kPa	Sum of hourly VPD values after sunrise above which afternoon stomatal reopening will not occur	8.0
fPAW	PAWt	%	Minimum non-limiting percentage of soil water	50
fSWP	SWPmax	MPa	Maximum soil water potential (SWP) below which gsto will start to decline	N/A
	SWPmin	MPa	SWP at which gsto reaches fmin	N/A
fO3	POD0SPEC/POD6SPEC	mmol m−2	Threshold flux at which O3-induced senescence begins	28
	exponent	constant	Rate of gsto decline with increasing flux accumulation	25