Attribution of Chemistry-Climate Model Initiative (CCMI) ozone radiative flux bias from satellites

The top-of-atmosphere (TOA) outgoing longwave flux over the 9.6 µm ozone band is a fundamental quantity for understanding chemistry–climate coupling. However, observed TOA fluxes are hard to estimate as they exhibit considerable variability in space and time that depend on the distributions of clouds, ozone (O₃), water vapor (H₂O), air temperature (T_a), and surface temperature (T_s). Benchmarking present-day fluxes and quantifying the relative influence of their drivers is the first step for estimating climate feedbacks from ozone radiative forcing and predicting radiative forcing evolution.

To that end, we constructed observational instantaneous radiative kernels (IRKs) under clear-sky conditions, representing the sensitivities of the TOA flux in the 9.6 µm ozone band to the vertical distribution of geophysical variables, including O₃, H₂O, T_a, and T_s based upon the Aura Tropospheric Emission Spectrometer (TES) measurements. Applying these kernels to present-day simulations from the Chemistry-Climate Model Initiative (CCMI) project as compared to a 2006 reanalysis assimilating satellite observations, we show that the models have large differences in TOA flux, attributable to different geophysical variables. In particular, model simulations continue to diverge from observations in the tropics, as reported in previous studies of the Atmospheric Chemistry Climate Model Intercomparison Project (ACCMIP) simulations. The principal culprits are tropical middle and upper tropospheric ozone followed by tropical lower tropospheric H₂O. Five models out of the eight studied here have TOA flux biases exceeding 100 mW m⁻² attributable to tropospheric ozone bias. Another set of five models have flux biases over 50 mW m⁻² due to H₂O. On the other hand, T_a radiative bias is negligible in all models (no more than 30 mW m⁻²). We found that the atmospheric component (AM3) of the Geophysical Fluid Dynamics Laboratory (GFDL) general circulation model and Canadian Middle Atmosphere Model (CMAM) have the lowest TOA flux biases globally but are a result of cancellation of opposite biases due to different processes. Overall, the multi-model ensemble mean bias is $-\mathrm{133}\pm \mathrm{98}$ mW m⁻², indicating that they are too atmospherically opaque due to trapping too much radiation in the atmosphere by overestimated tropical tropospheric O₃ and H₂O. Having too much O₃ and H₂O in the troposphere would have different impacts on the sensitivity of TOA flux to O₃ and these competing effects add more uncertainties on the ozone radiative forcing. We find that the inter-model TOA outgoing longwave radiation (OLR) difference is well anti-correlated with their ozone band flux bias. This suggests that there is significant radiative compensation in the calculation of model outgoing longwave radiation.

1 Introduction

Tropospheric ozone (O₃) is the third important anthropogenic greenhouse gas (GHG) in terms of radiative forcing (RF) as a consequence of O₃ precursor and methane (CH₄) emission increases from pre-industrial times to the present day. Tropospheric O₃ adjusted RF ranges widely from +0.2 to +0.6 W m⁻² computed from chemistry–climate model ensembles (IPCC AR5, 2013) (Bowman et al., 2013; Stevenson et al., 2013). The large uncertainty of the tropospheric O₃ RF is driven in part by the model responses to climate change. Without a good long-term record of the historical O₃ levels (Young et al., 2017; Gaudel et al., 2018), such estimates are highly dependent on the model assumptions of past O₃ levels. Differences between models in physical climate, chemical, and radiative processes conspire to complicate the assessment of the accuracy of these RF calculations. Consequently, a method to disentangle the key players caused the model differences to observations as well as the difference between the models is critical for robust estimates of chemistry–climate coupling.

About 80 % of tropospheric O₃ RF is due to O₃ longwave absorption, with the remaining 20 % from the shortwave absorption (IPCC AR5, 2013). In the longwave, 97 % of the total longwave absorption is in the 9.6 µm O₃ band (Rothman et al., 1987). The global outgoing longwave radiation (OLR) spectra were first observed from space for a few months in 1970. Radiance observations were taken during April 1970 and January 1971 by the NASA Infrared Interferometric Spectrometer (IRIS) and then from October 1997 for 9 months by the Interferometric Monitor of Greenhouse Gases (IMG) instrument, on board the Japanese Advanced Earth Observing Satellite “Midori” (ADEOS) satellite. Harries et al. (2001) showed that the changes in the greenhouse gas features between the observed spectra taken 30 years apart by these two instruments suggest increases in greenhouse gas forcing. Over the last two decades, a new generation of thermal infrared satellite instruments has provided a unique opportunity to continuously monitor the outgoing radiances covering the 9.6 µm O₃ band globally, such as NASA's Tropospheric Emission Spectrometer (TES) and Atmospheric Infrared Sounder (AIRS), ESA's Infrared Atmospheric Sounding Interferometer (IASI), and NOAA's Cross-track Infrared Sounder (CrIS). These valuable long-term global measurements can be used to derive the top-of-atmosphere (TOA) O₃ band flux and the sensitivity of the flux to the vertical distributions of O₃, defined as instantaneous radiative kernels (IRKs) (Worden et al., 2011; Doniki et al., 2015).

The TES-observed global TOA outgoing fluxes at the 9.6 µm O₃ band in clear skies (Fig. 1) show strong geographic variations as a result of the short lifetime of O₃ (Worden et al., 2011; Bowman et al., 2013). Consequently, the global O₃ GHG effect is more unevenly distributed than long-lived GHGs, such as CO₂. In addition, the variations of the TOA fluxes are not only highly dependent on the distributions of O₃ but are also dependent on water vapor (H₂O), air temperature (T_a), and surface temperature (T_s) (Kuai et al., 2017).

Figure 1 The clear-sky TES observed TOA flux at 9.6 µm O₃ band, annually averaged in 2006.

There is an additional factor where the large-scale atmospheric structure sets the overall atmospheric opacity, which describes the fraction of the light that fails to pass through the atmosphere due to the absorption or scattering. For example, O₃ changes in more opaque regions, e.g., the western Pacific, a wet region due to convection, result in a much smaller change in TOA flux than in more transparent regions, e.g., the Middle East, a dry region due to downwelling (Kuai et al., 2017). This opacity has a direct impact on radiative forcing calculations.

Chemistry–climate models diverge significantly in the simulation of these processes, which are difficult to disentangle because it is hard to quantify the response of the TOA flux due to the change in atmospheric opacity. In this study, we introduce a method to use observational-based IRKs to quantitatively estimate the contributions of the model biases in O₃, H₂O, T_a, and T_s to the TOA flux biases.

The presence of clouds is the primary control on atmospheric opacity. Under the cloudy sky conditions, the roles of these variables other than clouds on TOA flux are much weaker. In addition, the variation in clouds could affect model estimates not only of the ozone but also of the flux sensitivity to ozone and other variables. Both ozone and sensitivity will impact the ozone radiative flux but in opposite directions. With cloud cover, the O₃ loss will be reduced. That means too many clouds would lead to more ozone production. The presence of the cloud would also cause weaker flux sensitivity to O₃ and other variables (IRKs). Therefore, the cloud effect is a battle between the impact on ozone estimation and the radiative sensitivity to ozone (IRK). The differences in cloud variations between the models will complicate the radiative effect. Furthermore, the study of the cloud effect is also currently limited by the global observations of total cloud cover and IRK product under realistic cloud conditions. Without knowing which models have better cloud cover, we benefit from using IRK based on the observed cloud-free data by TES. Therefore, here, we first try to access the role of O₃, H₂O, T_a, and T_s in the variation of the TOA flux without cloud effect.

Worden et al. (2008) first attempted to disentangle these effects from satellites. They subsequently developed the IRK in Worden et al. (2011) for O₃, which is used in this study as a powerful tool to attribute model variability. IRKs for O₃ represent the sensitivity of TOA fluxes to the vertical distributions of the observed O₃. Aghedo et al. (2011) applied the TES IRKs to evaluate the O₃ radiative effect of chemistry–climate models' O₃ biases in the Atmospheric Chemistry Climate Model Intercomparison project (ACCMIP) (Lamarque et al., 2013). Bowman et al. (2013) found model OLR bias due to O₃ is correlated with RF in the ACCMIP models. This correlation helped to reduce the intermodel divergence in RF by about 30 % (Myhre et al., 2013). Doniki et al. (2015) updated the IRKs' calculation with a more accurate but computationally more complicated method, a five-angle Gaussian integration (GI) method, to replace the anisotropic approximation. They computed the O₃ IRKs with IASI observations and also showed that between the two methods there are about 20 % differences in IRKs and about 20 %–25 % differences in the longwave radiative effect (LWRE). They also found that the day and night difference of LWRE is mainly controlled by the T_s change instead of O₃ amount change. Kuai et al. (2017) updated the computational method for the TES O₃ IRK product with the five GI method and revealed the hydrological controls on the global distribution of the O₃ GHG effect. The study showed that H₂O, T_a, and T_s affect the O₃ IRK strength through relative humidity.

Therefore, the TOA flux in the 9.6 µm band depends on more than O₃. Consequently, in this study, we expand the TES observation-based IRKs to other quantities, including H₂O profiles, T_a profiles, and T_s. We apply these IRKs to help understand the reasons for the model divergence in the TOA flux.

The questions that have never been answered before include the following. (1) How do the model-based flux and the flux sensitivity compare to the observational-based flux and sensitivity? (2) How do they compare between the models? (3) How do the flux biases in models relate to the RF variation? Thus, benchmarking present-day O₃ band flux is the first step in answering all these questions and would help to further understand the correlations between the bias in TOA flux and the bias in O₃ RF, and eventually improve the estimation of the climate feedbacks from O₃ forcing.

To benchmark the model-simulated geophysical quantities, a recently developed multi-species multi-satellite Tropospheric Chemistry Reanalysis (TCR) product (Miyazaki et al., 2015) is used in this study to compare to the model results. This chemical reanalysis assimilates data from multiple satellites with sensitivity over complementary parts of the atmosphere, which provides better information than single-species chemical data assimilation. Satellite observations have the occasional issue of temporal discontinuity due to instrument performance and irregular spatial coverage, which can be circumvented by chemical data assimilation. Miyazaki et al. (2015) showed that statistically the model error against independent aircraft and ozonesonde observations in the assimilated species, e.g., O₃, NO₂, and CO, is significantly reduced. The multi-species assimilation improves the Northern/Southern Hemisphere OH ratio and provides the emission estimates with interannual variation. The comparison of O₃ reanalysis to the ACCMIP ensemble O₃ simulation in Miyazaki and Bowman (2017) quantified the model discrepancies in terms of seasonal amplitude, spatial variability, and interhemispheric gradient. For example, the ensemble mean is 6–11 ppb too high in the northern extratropics, while up to 18 ppb too low in the southern tropics over the Atlantic in the lower troposphere. In this study, we use the same O₃ reanalysis data (Miyazaki and Bowman, 2017) to understand the model bias in the CCMI project (Morgenstern et al., 2017), a follow-up model intercomparison study for ACCMIP. The multi-species assimilation also provides the opportunity to optimize the chemical-related species of O₃ and the emission sources of the precursors simultaneously. Further work by Miyazaki et al. (2017) showed that the surface emission of nitrogen oxides (NO_x) over a 10-year period (2005–2014) has a positive trend in regions including India, China, and the Middle East, but a negative trend over the US, southern Africa, and western Europe. The global total emission stays almost constant between 2005 (47.9 Tg N yr⁻¹) and 2014 (47.5 Tg N yr⁻¹). Therefore, the O₃ reanalysis data from TCR represent the state of the art for the current knowledge of the global distribution of tropospheric O₃ by combining the complementary information from model and satellite observations for O₃ and its precursors.

In this paper, we demonstrate a method to use the IRK products and the model biases relative to the reanalyzed tropospheric composition (O₃ and H₂O) and atmospheric state (T_s and T_a) to quantitatively attribute the radiative biases of the flux in a suite of CCMI models to these dominant components. The method and IRKs are described in Sect. 2. The models and reanalysis data are introduced in the next section. Section 4 discusses the intercomparison between models' flux biases, the bias attribution to the dominant components, and the geospatial and vertical distribution of the biases. Lastly, the conclusion and future directions are summarized in Sect. 5.

2 Instantaneous radiative kernels (IRKs) for the climate variables

The TOA flux in the 9.6 µm O₃ band (Fig. 1) is defined as

$$\begin{array}{}\text{(1)}& {\displaystyle }{F}_{\mathrm{TOA}}=\underset{v}{\int }\underset{\mathrm{0}}{\overset{\mathrm{2}\mathit{\pi }}{\int }}\underset{\mathrm{0}}{\overset{\frac{\mathit{\pi }}{\mathrm{2}}}{\int }}{L}_{\mathrm{TOA}}(v,\mathit{\theta },\mathit{\varphi },q)\sin \left(\mathit{\theta }\right)\cos \left(\mathit{\theta }\right)\mathrm{d}\mathit{\theta }\mathrm{d}\mathit{\varphi }\mathrm{d}v,\end{array}$$

where v is the frequency, integrated over the O₃ band from 980 to 1080 cm⁻¹. $L_{\mathrm{TOA}}(v,\mathit{\theta },\mathit{\varphi },q)$ is the upwelling TOA radiance at frequency v, zenith angle θ, and azimuth angle ϕ. We assume here that the radiance is symmetric in the azimuthal direction. The outgoing TOA radiances, L_TOA, are also a function of the atmospheric state, which is represented by variable “q”, e.g., H₂O, O₃, and T_a, that is in turn a function of altitude, z.

The IRKs (Eq. 2) represent the sensitivities of the TOA radiative flux in the 9.6 µm O₃ band to the changes in the vertical distribution of an atmospheric variable.

$$\begin{array}{}\text{(2)}& {\displaystyle \frac{\partial F_{\mathrm{TOA}}}{\partial q\left(z_l\right)}}=\underset{v}{\int }\underset{\mathrm{0}}{\overset{\mathrm{2}\mathit{\pi }}{\int }}\underset{\mathrm{0}}{\overset{\frac{\mathit{\pi }}{\mathrm{2}}}{\int }}{\displaystyle \frac{\partial L(v,\mathit{\theta },\mathit{\varphi },q)}{\partial q\left(z_l\right)}}\sin \left(\mathit{\theta }\right)\cos \left(\mathit{\theta }\right)\mathrm{d}\mathit{\theta }\mathrm{d}\mathit{\varphi }\mathrm{d}v,\end{array}$$

where z_l is altitude in discretized level l. When q represents the T_s, z_l becomes a single surface value at l=0. The partial derivative term on the right side of the equation is the spectral radiance Jacobians calculated analytically by the TES radiative transfer model.

In this study, we expanded the TES global O₃ IRKs to IRKs with respect to H₂O, T_a, and T_s. These TOA flux sensitivities still refer to the spectral window region in the 9.6 µm O₃ band for the flux. All the kernels are computed with the five-angle GI method (Doniki et al., 2015; Kuai et al., 2017). Figures 2a, c, and e show examples of IRK profiles for O₃, H₂O, and T_a for 2006. The TOA flux is most sensitive to each variable at very different vertical levels. The O₃ IRK peaks in the middle and upper troposphere (600 to 200 hPa), a higher level than the peaks in both H₂O and T_a IRKs. The middle and upper tropospheric O₃ near 500 hPa has the largest impact on the TOA flux change (close to 1 mW m⁻² ppb⁻¹ in the tropics). The H₂O IRK peaks near 700 hPa, a little higher than the T_a IRK. The T_a IRK is maximal closest to the surface, suggesting that the O₃ band flux is most sensitive to boundary layer T_a near 900 hPa. The strength of the peaks decreases with increasing latitude for all the three variables but the peak altitude does not change significantly except for the H₂O IRKs in the polar region, which peaks at a slightly higher level than in lower latitudes.

Figure 2 TES 2006 IRK for four primary components (O₃, H₂O, T_a, and T_s). Panels (a), (c), and (e) are IRKs of latitudinal band averages in the tropics (30^∘ S–30^∘ N), midlatitudes of both hemispheres (30–60^∘), and high latitudes of both hemispheres (60–90^∘). The figures below them are the pole-to-pole vertical distribution of the zonally averaged IRK. The global distribution of IRK for T_s is plotted for winter season (December to February) (g) and summer season (June to August) (h).

In addition, the T_s IRK is greater than zero, which means increases in T_s would increase the outgoing TOA flux. However, the IRKs for the GHGs, i.e., H₂O and O₃, are negative, because the increase in gas concentrations reduces the upwelling flux at TOA due to radiative absorption by the gas.

The global vertical distributions of the zonal averaged kernels for O₃, H₂O, and T_a are also shown below their profile plots in Fig. 2b, d, and f. The sensitivities of the TOA flux to these three variables are strongest in the tropics and decrease with latitude. Furthermore, the IRK for T_s is also shown in Fig. 2g and h. Unlike the other IRKs, the T_s IRK is not a function of altitude, so we show the winter (December–February) and summer (June–August) seasonal average of its global distribution. The flux sensitivities are found to be largest over the major deserts, like the Sahara, the Middle East, and Australia, corresponding to the regions with the highest values of T_s. We also notice that the values of the T_s IRKs in the Intertropical Convergence Zone (ITCZ) are much lower than those in the subtropics, which suggests that the atmosphere opacity has an impact on the strength of the T_s IRKs.

3 A method to attribute the flux biases

The flux biases between observations and models under the clear-sky conditions can be described as

$$\begin{array}{}\text{(3)}& \begin{array}{rl}& \mathit{\delta }{F}_{\mathrm{TOA}}^i\approx {\sum }_{l\in L}{\displaystyle \frac{\partial F_{\mathrm{TOA}}^i}{\partial {{\mathrm{O}}_{\mathrm{3}}}^{i,l}}}\mathit{\delta }{{\mathrm{O}}_{\mathrm{3}}}^{i,l}+{\sum }_{l\in L}{\displaystyle \frac{\partial F_{\mathrm{TOA}}^i}{\partial {{\mathrm{H}}_{\mathrm{2}}\mathrm{O}}^{i,l}}}\mathit{\delta }{{\mathrm{H}}_{\mathrm{2}}\mathrm{O}}^{i,l}\\ & +{\sum }_{l\in L}{\displaystyle \frac{\partial F_{\mathrm{TOA}}^i}{\partial T_{\mathrm{a}}^{i,l}}}\mathit{\delta }{T}_{\mathrm{a}}^{i,l}+{\sum }_{l\in L}{\displaystyle \frac{\partial F_{\mathrm{TOA}}^i}{\partial T_{\mathrm{s}}^{i,l}}}\mathit{\delta }{T}_{\mathrm{s}}^{i,l},\end{array}\end{array}$$

where $\partial F_{\mathrm{TOA}}^i$, is the total TOA flux bias in the O₃ band at the ith location. The four terms on the right-hand side of the equation are the products of the IRKs and the biases in the geophysical quantities (i.e., O₃, H₂O, T_a, and T_s). These biases are then vertically integrated on index l, over domain L, which in our case is the troposphere. The summation is the vertical integral from the surface to the tropopause.

Here, we assume that the biases due to other physical processes, e.g., surface emissivity or other atmospheric species, have much less influence on the TOA flux variation. For example, the model bias in global emissivity is not accessible but is believed to be quite small compared to O₃, H₂O, T_a, and T_s. We also assume that the nonlinearity terms are much smaller than these four first-order terms.

Following Bowman et al. (2013), the delta terms in Eq. (3) are the model biases with respect to the reanalysis data, defined as below:

$$\begin{array}{}\text{(4)}& {\displaystyle }\mathit{\delta }{q}^{i,l}=q_{\mathrm{mod}}^{i,l}-q_{\mathrm{assim}}^{i,l},\end{array}$$

where $q_{\mathrm{mod}}^{i,l}$ and $q_{\mathrm{assim}}^{i,l}$ represent the model and reanalysis O₃, H₂O, T_a, or T_s at the ith location and the lth altitude level, respectively.

The mean flux bias or the mean bias components from tropospheric uncertainties are calculated from Eqs. (3) and (4) as

$$\begin{array}{}\text{(5)}& {\displaystyle }\mathit{\delta }{F}_q^j={\displaystyle \frac{\mathrm{1}}{N_j}}{\sum }_{i\in D_j}{\sum }_{l\in L}{w}_i{\displaystyle \frac{\partial F_{\mathrm{TOA}}^i}{\partial q^{i,l}}}(q_{\mathrm{mod}}^{i,l}-q_{\mathrm{assim}}^{i,l}),\end{array}$$

where w_i is area weighted for the latitude bands, D_j is a set of observed locations, N_j is the number of locations in the domain of D_j and tropospheric levels of L up to the tropopause. We use the chemical tropopause O₃ =150 ppb (Naik et al., 2005; Hansen et al., 2007; Bowman et al., 2013; Kuai et al., 2017). The domain of D_j can be zonal bands for the zonal mean or global area for the global mean, respectively. The global mean of the flux bias and its components will be denoted as $\stackrel{\mathrm{‾}}{\mathit{\delta }F}$ and $\stackrel{\mathrm{‾}}{\mathit{\delta }{F}_q}$, respectively.

4 Chemistry–climate models and the reanalysis data

4.1 Models and simulations

We analyze six models from the CCMI study (Table 1) (Hegglin and Lamarque, 2015; Morgenstern et al., 2017; Eyring et al., 2013). It is a combined activity of the International Global Atmospheric Chemistry (IGAC) and Stratosphere-troposphere Processes And their Role in Climate (SPARC) (Randel et al., 2004). The CCMI coordinates a number of model experiments that capture the variability and evolution of air quality, tropospheric chemistry, stratospheric O₃, and global climate. This approach builds on the legacy of previous chemistry–climate model intercomparisons, such as the Chemistry-Climate Model Validation (CCMVal, Morgenstern et al., 2010; SPARC, 2010, Eyring et al., 2010) and the ACCMIP. In this study, we use the experiment REF-C1, which is analogous to the REF-B1 experiment of CCMVal-2 (Table S30 in Morgenstern et al., 2017). REF-C1 requires using historic forcing and observed sea surface conditions. The models are free-running and simulate the recent past (1960–2010). We did not choose to use REF-C1SD (specified dynamics) because specified dynamics nudged the wind and temperature of the model to be constrained to the reanalysis data. The long-term climatological biases relative to the reanalysis between the models are minimized. Our study aims to find a correlation between the present-day radiative bias and the RF from present day to future by the model predictions. Therefore, we prefer to keep the model differences in simulating longer-term climatology between their free runs.

Table 1 The chemistry–climate models and their experiment simulations.

We note that SOCOL3 and EMAC are both based on different versions of the ECHAM5 climate model. We also added two additional model simulations with AM3 from NOAA and CESM from NCAR. These two simulations are not the specific CCMI experiment run; however, these two models have been used in many studies, and including them in this study provides more useful information on the TOA flux diversity among the most recent models.

4.2 Tropospheric Chemistry Reanalysis (TCR-1) data

We computed the biases in the geophysical variables between the model and the reanalysis data. To compute the O₃ bias in models, we used the satellite-based O₃ reanalysis from multi-constituent multi-satellite data assimilation: Tropospheric Chemistry Reanalysis version 1 (TCR-1) (Miyazaki et al., 2015; Miyazaki and Bowman, 2017) as the best synthesis of the observations. The reanalysis provides comprehensive spatiotemporal and multi-variable evaluation of model performance that complements direct comparisons against individual measurements, which may suffer from significant sampling bias (Miyazaki and Bowman, 2017).

TCR-1 assimilated multiple species data from multiple satellite products for the period from 2005 to 2017, e.g., combined TES and Microwave Limb Sounder (MLS) observations for O₃, integrated Ozone Monitoring Instrument (OMI), Scanning Imaging Absorption Spectrometer for Atmospheric Chartography (SCIAMACHY), and Global Ozone Monitoring Experiment-2 (GOME-2) for tropospheric NO₂ column, Measurements Of Pollution In The Troposphere (MOPITT) for CO, and MLS for HNO₃. TCR-1 used a global chemistry–transport (CTM) Model for Interdisciplinary Research on Climate with chemistry (MIROC-Chem; Watanabe et al., 2011) as a forecast, which includes 92 species and 262 reactions. The model has 2.8^∘ horizontal resolution with 32 vertical layers up to 4 hPa. The data assimilation was based on an ensemble Kalman filter with 32 ensemble members, which was used to simultaneously optimize concentrations and emissions of various species.

As summarized by Miyazaki and Bowman (2017), the mean bias in the reanalysis dataset against the World Ozone and Ultraviolet Data Centre (WOUDC) ozonesonde observations is from −3.9 to −2.9 ppb at the NH high latitudes (55–90^∘ N); −0.9 to −0.1 ppb at the NH midlatitudes (15–55^∘ N); and −1.0 to −0.1 ppb at the SH midlatitudes (55–15^∘ S), between 850 and 500 hPa. On average, the bias is about 0.9 ppb at the tropics and midlatitudes between 500 and 200 hPa. These biases are much smaller than biases in the model simulation without data assimilation, demonstrating that the multi-satellite data assimilation provides comprehensive constraints on the entire tropospheric profile of O₃.

For the purpose of consistency, we also use outputs of H₂O, T_a, and T_s from the reanalysis to estimate the model biases. In the reanalysis calculation, meteorological fields simulated by the atmospheric general circulation model MIROC-AGCM (Watanabe et al., 2011) were nudged toward the 6-hourly ERA-Interim meteorological reanalysis (Dee et al., 2011) for zonal wind (τ=1 d) and temperature (τ=3 d) to reproduce past meteorological fields while simulating short-term (<6 h) meteorological variations, which were used to drive the CTM, as similarly employed in CCMI C1SD simulations. Thus, the reanalysis dataset provides realistic and comprehensive estimates for both chemical and meteorological fields required for the TOA flux evaluations.

5 Results

5.1 The latitudinal distribution of the TOA flux bias

Figure 3 shows the latitudinal distribution of the zonal and annual mean of the TOA flux bias from each model relative to the reanalysis. The largest divergence between the models is located at the tropics where most models underestimate the flux, with the exception of CMAM. The low bias in the model ensemble implies the model atmosphere is more opaque than the chemical reanalysis, leading to a 133 mW m⁻² outgoing flux reduction on average. The TOA flux in an opaque atmosphere is less sensitive to the changes in tropospheric composition than a more transparent one. Under those conditions, the models would underestimate the radiative feedback from composition since the IRKs estimated under an opaque atmosphere will be weaker than those under a realistic (more transparent) atmosphere.

Figure 3 The latitudinal distribution of the zonal flux bias (model – reanalysis) with latitude weight.

Two models that have larger low biases at the equatorial region than other models are SOCOL3 and MRI-ESM1r1. Their global means of the flux bias are more than −200 mW m⁻² (Table 2). The following analysis will help to clarify the source of the bias in the models.

Table 2 The global mean of the flux bias (mW m⁻²) and the dominant components due to tropospheric O₃, H₂O, T_a, and T_s. The numbers in parentheses are the standard deviation of the zonal distribution. For the ensemble, the standard deviation is computed from the variation between the models.

^a The models that have relative small global and annual averaged TOA flux bias. ^b The extreme values for the large biases.

5.2 Flux bias attribution

The total TOA flux bias is caused by biases from atmospheric composition and temperature. In order to determine the primary drivers of these biases, we apply the IRKs to the differences between model and the chemical reanalysis as described in Eq. (3). Figure 4 shows the contribution of O₃ (blue), H₂O (green), T_s (red), and T_a (yellow) for each model to the total TOA flux bias (black). The global mean bias is summarized in Table 2.

Figure 4 The attribution of the total TOA flux bias for each model to four dominant components and their latitudinal distribution. The black curves are the same as the colored curves in Fig. 3.

In general, O₃ and H₂O are the two dominant drivers for most models where the large biases are concentrated in the tropics and subtropics. There are only three models (GEOSCCM, CMAM, and CESM) whose O₃ radiative biases ($\stackrel{\mathrm{‾}}{\mathit{\delta }{F}_{{\mathrm{O}}_{\mathrm{3}}}}$ in Table 2) are less than 50 mW m⁻² and are almost negligible zonally. While the flux bias is better represented in these models, it does not follow that they represent tropospheric O₃ more accurately, as will be shown in the following section. The other five models (AM3, SOCOL3, EMAC-L47MA, EMAC-L90MA, and MRI-ESM1r1) have significant negative peaks at low latitudes (Figs. 3 and 4), actually resulting from their strong O₃ contributed biases (from 80 to 180 mW m⁻²; numbers are footnoted as b in Table 2).

The TOA flux bias from H₂O is the second largest component for most models. Similar to O₃, most models show the fluxes are biased low in the tropics due to the H₂O uncertainties with the exception of CMAM, which has the strongest global mean bias (127.9 mW m⁻²). Note that, in the reanalysis, no data assimilation (or nudging) was applied for specific humidity. Watanabe et al. (2011) demonstrated a dry bias in the lower troposphere and a wet bias in the middle and upper troposphere in MIROC-AGCM, primarily attributable to temperature biases. Nevertheless, the reported H₂O biases can be greatly reduced in the reanalysis because of the nudging applied for temperature.

The flux bias due to T_a is found to be negligible in all models, which indicates that the model T_a estimates provide reasonable radiative fluxes. T_s radiative bias is also meridionally weak relative to the flux bias in O₃ and H₂O (Fig. 4). With the exception of CMAM, the T_s ensemble global mean bias is less than 35 mW m⁻² (see Table 2). Figure 4 suggests the strong bias from T_s in CMAM (−100.2 mW m⁻²) comes from the two subtropical regions.

Interestingly, the positive flux bias due to H₂O (127.9 mW m⁻²) is compensated by the negative flux bias due to T_s (−100.2 mW m⁻²) in CMAM, leading to the lowest global mean in $\stackrel{\mathrm{‾}}{\mathit{\delta }F}$ (42.9 mW m⁻², calculated with Eq. 5). This compensation is also true for AM3 but between a positive H₂O radiative bias (87.7 mW m⁻²) and negative O₃ radiative bias (−140 mW m⁻²). This analysis reveals that these two models are both right but for wrong – and opposite – reasons.

However, all the other models have a strong negative global mean bias and are mostly driven by the two major components (O₃ and H₂O). SOCOL3 and MRI-ESM1r1 are the two models that have the strongest low bias up to −200 mW m⁻², which is mainly due to their strong O₃ radiative bias (−180 mW m⁻²). Their O₃ estimates are both biased high in the tropics and subtropics. We will show later that such bias is particularly strong in the upper troposphere.

5.3 Vertically resolved radiative bias of the O₃, H₂O, and T

The zonal flux biases among the models are both significant and mainly in the tropics. However, those biases are the vertically integrated product of the model profile bias and the IRKs both with their own vertical structures. The vertically resolved radiative bias can provide more insight into the processes leading to the biases. To further investigate, we examined the vertically resolved flux bias for O₃, H₂O, and T_a (Figs. 5–7) and the global distribution for T_s (Fig. 8). These are computed from Eq. (3) before the vertical summation. These figures show that the maximum contribution to the flux bias is a balance between the peak of the IRKs (Fig. 1) and the peak of the geophysical quantities' bias (Figs. 9–12). The positive tropical O₃ radiative bias for GEOSCCM, CMAM, and CESM is commonly centered in the midtroposphere, corresponding to the peak of the IRKs (Fig. 5). On the other hand, the primary O₃ flux bias contribution in the tropics for SOCOL3, EMAC-L47MA, EMAC-L90MA, and MRI-ESM1r1 is in the upper troposphere around 200 hPa even though the IRKs are roughly half the peak sensitivity. These strong negative biases exceed 15 mW m⁻².

Figure 5 Vertical resolved O₃ radiative bias. The black curves are the zero lines.

Figure 6 Vertical resolved H₂O radiative bias. The black curves are the zero lines.

Figure 7 Vertical resolved T_a radiative bias. The black curves are the zero lines.

Figure 8 Global distribution of the T_s radiative bias.

The strong tropical H₂O radiative bias collapses to shallower tropical regions below 400 hPa and is maximized near 800 hPa for most models exceeding 50 mW m⁻² (Fig. 6). CMAM has the unique and strongest net positive bias of above 50 mW m⁻² centered lower and close to 900 hPa. While most models' flux bias is centered near 800 hPa, particularly GEOSCCM, AM3, and CESM show a more vertically uniform – and opposing – flux bias.

Figure 7 indicates that tropical T_a radiative bias is largely negligible for vertical layers above 600 hPa as a consequence of the rapid decrease in sensitivity of the T_a IRKs. The maximum bias is in the lower troposphere between 900 hPa and the surface. CMAM and CESM both show the strongest positive bias exceeding 10 mW m⁻² over most of the tropics. However, CESM has a compensating negative bias from 700 to 800 hPa that leads to a mean global bias of only 6.4 mW m⁻² (Table 2), whereas CMAM has a positive bias throughout, leading to an atmospheric T_a radiative bias of 22.5 mW m⁻², the largest of the models studied here.

Surprisingly, the model ensemble T_s turns out to be the second largest contributor to the total bias (Table 2), instead of H₂O, driven primarily by three models: CMAM, CESM, and GEOSCCM, as shown in Fig. 8. CMAM shows a negative bias that covers all of Africa, exceeding 500 mW m⁻², and Asia centered over India. Consequently, CMAM has the largest total bias ($-\mathrm{100.2}\pm \mathrm{93.3}$ mW m⁻²). CESM and GEOSCCM T_s radiative biases, on the other hand, are centered at high latitudes in the western hemisphere over the eastern US and Canada, exceeding 300 mW m⁻².

The vertically and spatially concentrated radiative biases provide clues as to what processes are the most important for the total flux bias. These processes drive the distribution of the constituents, which we will discuss in detail in the next sections.

5.4 The spatial source of TOA flux bias

The source of the attributed flux biases can be traced back to their spatial origins, which can provide more insight into the underlying processes and the differences between the models.

5.4.1  O₃ bias

Figure 9 shows a vertically resolved zonal averaged distribution of O₃ biases between the model and the chemical reanalysis similar to that in Fig. 5. Three models (GEOSCCM, CMAM, and CESM) have the weakest globally averaged O₃ radiative bias reported in Table 2 (−33.5, −7.3, and 3.3 mW m⁻²). These three models also have the lowest O₃ bias in tropical troposphere on average (−1.1 ppb, −1.3 ppb, and 3.0 ppb, reported in Table 3) and as a consequence have weaker radiative bias in the region with the strongest O₃ IRK globally (1.0, 1.4, and 2.1 mW m⁻²; reported in Table 4). On the other hand, the global O₃ bias is greater than 7 ppb for all the other models and results in a large O₃ radiative bias in the tropics, especially SOCOL3 (13.4 ppb in the tropical O₃ bias and −9.2 mW m⁻²; see Tables 3 and 4) and MRI-ESM1r1 (13.7 ppb and −10 mW m⁻²).

Figure 9 The zonal averaged vertical–latitudinal distribution of O₃ model biases to the TCR-1 O₃ assimilation data. The black curves are the zero lines.

Table 3 Models' bias in the tropical troposphere between 25^∘ S and 25^∘ N, and below 200 hPa.

Table 4 Models' flux bias (mW m⁻²) in the tropic troposphere between 25^∘ S and 25^∘ N, and below 200 hPa.

GEOSCCM, CMAM, and CESM commonly have a vertically compensated pattern in the tropics that is biased high in the upper troposphere while biased low in the middle and lower troposphere (Fig. 9). Their O₃ low biases in the middle troposphere are approximately 5 to 10 ppb, where the peak of the IRK centered, but the high biases in the upper troposphere are about 5 ppb. Such a high–low pattern leads to compensation during the vertical integration through the troposphere into the radiative effect at the top of the atmosphere. The corresponding vertical resolved O₃ radiative bias for these three models in Fig. 5 shows the consistent tropical vertical distribution but in an opposite sign since the O₃ IRK is negative. In contrast, the other five models have vertically systematic biases high in the tropical O₃, and the biases increase from the middle troposphere to the upper troposphere. Especially SOCOL3 and MRI-ESM1r1 strongly overestimate O₃ by more than 20 ppb in a wide region of tropical upper troposphere. The O₃ radiative biases in this region remain significantly high, stronger than −15 mW m⁻², causing these two models to have the highest O₃ radiative biases in the global and annual mean (both about −183 mW m⁻²).

The systematic bias in the entire tropical tropospheric O₃ and strong overestimation of upper troposphere in SOCOL3 and MRI-ESM1r1 could be caused by several factors. For example, the transport from the lower stratosphere could be too high. Alternatively, precursor emissions of tropospheric O₃ could also be too high. The analysis with the spatially explicit biases provides important clues to implicate the specific processes that individual modeling groups can investigate.

Figure 10 The zonal averaged vertical–latitudinal distribution of H₂O biases (model to the ERAI reanalysis data). The black curves are the zero lines.

The GEOSCCM has been used to study the tropospheric O₃ response to variations in the El Niño–Southern Oscillation (ENSO), where Oman et al. (2011, 2013) compared the model to satellite observations. These regular comparisons may have led to the improved simulation of tropospheric O₃ profiles and consequently lower vertical O₃ bias. The GEOSCCM model in the CCMI study uses the tropospheric–stratospheric chemical package developed within the Global Modeling Initiative (GMI) program (Duncan et al., 2007), which has more realistic ozone chemistry, an internally generated quasi-biennial oscillation, an improved air–sea roughness parameterization and other improvements (Oman and Douglass, 2014).

Nielsen et al. (2017) showed that GEOSCCM has successfully reproduced the changes in the quasi-global (60^∘ S–60^∘ N) annual mean trend in total O₃ column since 1960s to the present day. For the present-day atmosphere, simulated tropospheric partial column O₃ from GESCCM Ref-C1 for CCMI was compared to satellite observations of OMI and MLS (Ziemke et al., 2011). The differences are mostly a few Brewer–Dobson units (DU) except in the Northern Hemisphere subtropics and middle latitudes in autumn and winter with the 4–6 DU biases which are under investigation.

Figure 11 The zonal averaged vertical–latitudinal distribution of T_a biases from models to the reanalysis data. The black curves are the zero lines.

The finding that SOCOL3 and MRI-EMS1r1 both have strong overestimates in the tropical upper troposphere is also understandable. SOCOL3 is the third generation of the coupled chemistry–climate model (CCM) SOCOL (modeling tools for studies of SOlar Climate Ozone Links). Several steps have been taken to improve the SOCOL model simulation of O₃. Stenke et al. (2013) first attempted to reduce the O₃ bias in their middle atmosphere by updating their middle-atmosphere general circulation with an advanced advection scheme. Revell et al. (2015) revealed that ozone precursor emissions are the biggest players that control the global mean change in tropospheric ozone. In a parallel study, Revell et al. (2018) developed an updated version of “SOCOL3.0”, “SOCOL3.1”, to reduce the tropospheric ozone bias. By improving the treatment of ozone sink processes, the tropospheric column ozone bias in “SOCOLv3.1” is reduced up to 8 DU, mostly due to the inclusion of N₂O₅ hydrolysis on tropospheric aerosols. We expect that the future similar analysis with the SOCOL3.1 could show a reduced flux bias for this model.

Meanwhile, the strong tropical upper tropospheric O₃ biases in MRI-EMS1r1 are believed to be related to the weak tropical convective updraft and the large lightning NO_x emissions in the model. The model with weak updraft fails to bring enough low O₃ air from the surface to the upper troposphere in the tropics or overestimates the upper tropospheric mixing of stratospheric ozone-rich air. In addition, the global lightning NO_x (LNO_x) emission used in MRI-EMS1r1 is 10 TgN yr⁻¹. The best estimate of annual mean LNO_x based on satellite data assimilation is 6.3 TgN yr⁻¹ (Miyazaki et al., 2014). The LNO_x in GEOSCCM is approximately 5 TgN yr⁻¹ (Martini et al., 2011), which shows less tropical upper tropospheric O₃ bias compared to MRI-EMS1r1. Thus, the overestimation of the O₃ precursor in the upper troposphere is another reason for too much O₃. Figure A1 shows the improvement in the radiative biases due to less O₃ bias in the experiment by half the LNO_x emissions in MRI-EMS1r1 (see the Appendix).

Figure 12 The global distribution of T_a biases from models to the reanalysis data.

In summary, the potential reasons for the prevalence of O₃ radiative bias in the tropical middle and upper troposphere in the models could be due to following facts: (1) the tropical O₃ IRK is strongest in this region (Fig. 2); (2) the largest O₃ bias in the models also centered in the same place (e.g., SOCOL3 and MRI-EMS1r1; Fig. 9); (3) the simulations with the systematic bias throughout the tropical troposphere, when vertically integrated, accumulated into a larger column bias when compared to the models with vertically random biases.

5.4.2 H₂O bias

H₂O turns out to be the primary contributor for three models (GEOSCCM, CMAM, and CESM) since their O₃ radiative bias is small. It is also the second dominant driver after O₃ in the other five models. Different from O₃, H₂O IRKs in Fig. 2 show the strongest sensitivity to the tropical lower troposphere centered at 800 hPa, where H₂O is most concentrated globally. We found the model biases in H₂O are strongest in the tropical lower troposphere. It explains why the strongest radiative bias from H₂O is also located in the tropical region near 800 hPa in all models as shown in Fig. 6. Figure 10 and Table 3 further help to indicate that H₂O is biased low only in two models, AM3 (−586.5 ppm) and CMAM (−506.6 ppm). We note that H₂O IRKs are also negative as O₃. Therefore, these two models have the unique overestimates in H₂O radiative bias at low latitudes (see Fig. 4), while all the other models are predominantly biased high in tropical H₂O concentrations, which result in the negative radiative biases.

Figure 13 The correlation of the ozone band TOA flux biases to the model-calculated broadband OLR (a) and the correlation of the attributed radiative components to the broadband OLR (b–e).

5.4.3 T_a bias

We found the T_a radiative biases in these model ensembles are all negligible. There are two reasons. One is that the T_a biases are small overall (less than 2 K) even at the tropical lower troposphere (below 1 K on average in Table 3). The other reason is that the compensation in the vertical integration helps to reduce the radiative bias at the top of atmosphere.

Figure 11 shows that the model biases in T_a range within ±2 K for all the models because the current chemistry–climate models have been well developed to simulate the global atmospheric T_a fields relative to reanalysis. The region with strongest sensitivity, identified by the T_a IRKs (Fig. 2), is the tropical lower troposphere (the region within $\pm {\mathrm{30}}^{\circ }$ and below 800 hPa). The T_a biases in the tropics shift between positive and negative vertically in most models except CMAM, which is systematically biased (Fig. 11). The oscillated T_a biases suggest that simulated air temperatures stay around the reanalyzed profiles. These models better represent the air temperature than the trace gases like H₂O and O₃. The oscillation around the reanalyzed profile leads to vertical compensation in the air T_a radiative bias. Therefore, the flux bias from T_a is a small component compared to the radiative bias from O₃ and H₂O. Figure 4 suggests the only model that has a small tropical peak in the T_a radiative component is CMAM, which has the strongest T_a radiative bias (22.5±40.5 mW m⁻² in Table 2) among all the models. Figure 11 shows that this model has a deep region with a strong bias of about 2 K in tropical areas and also persistently overestimated T_a vertically. Figure 7 further suggests the strong radiative bias in CMAM mainly comes from the tropical lower troposphere (>10 mW m⁻² below 800 hPa). The other models have vertical compensation in the tropics (less than 0.5 K bias on average; see Table 3), and therefore they are less biased in TOA flux (less than 1 mW m⁻² in Table 4). The two EMAC models both have strong biases at the southern high latitudes but still have a weak radiative effect in this region (see Fig. 7) due to much weaker T_a IRK at high latitudes.

5.4.4 T_s bias

The global distribution of the T_s biases indicates that the biases in sea surface temperature (SST) are smaller than the biases in land T_s for all the models (Fig. 12) because the CCMI experiment (REF-C1) selected in this study used the observed SST. CMAM is the model that has the strongest T_s radiative bias ($-\mathrm{100.2}\pm \mathrm{93.3}$ mW m⁻² in Table 2), which peaks in both subtropical regions (Fig. 4). These large negative biases are due to the large underestimates of the T_s over the major deserts, e.g., the Sahara, the Middle East, and Australia (Fig. 12). In other words, the real deserts' surface is hotter than the model's prediction. At the same time, the T_s IRKs at the subtropical desert region are also strongest globally, since the TOA flux is more sensitive to T_a when the atmosphere is transparent, which is due to the downdraft of the Hadley cell controlling the region (Kuai et al., 2017). The downwelling airflow results in less precipitation and less clouds, as well as higher T_s during summer over the desert surface. These factors cause the CMAM to have the largest T_s radiative bias compared to all the other models.

In contrast, two EMAC models and MRI-ESM1r1 also have strong high bias in T_s in Siberia (>4 K in Fig. 12), but the radiative bias is much less significant compared with the Middle East in CMAM in Fig. 8. The T_s IRKs are weaker during the winter season in high latitudes than the low latitudes if the T_s is low (Fig. 2). However, the IRKs at the subtropical desert region stay strong during winter. Therefore, the annual mean of the T_s radiative bias is much weaker in Siberia in two EMAC models and MRI-ESM1r1 than the Middle East region in CMAM. Consequently, the global annual means of the T_a radiative biases for two EMAC models and MRI-ESM1r1 are small, although the large biases in T_s are found in their Siberian region.

6 Correlation to the broadband OLR

The analysis up to this point has been limited to the 9.6 µm band. We posed the question as to whether biases in this band could provide any insight into biases in the entire OLR band. To that end, we found an anti-correlation ($R=-\mathrm{0.6}$) between the global mean of the O₃ band flux biases and the clear-sky broadband OLR calculated internally by the models, as shown in Fig. 13a. The CMAM OLR is inconsistent with the ensemble (more than 2 W m⁻² higher than all the other models), and therefore it is excluded in the correlation. Interestingly, a very similar regression line and anti-correlation coefficient ($R=-\mathrm{0.6}$) are found between the O₃ radiative bias and the broadband OLR (Fig. 13b). The similar regression line indicates that the O₃ radiative bias dominates the 9.6 µm TOA flux distribution, which is confirmed by the attribution analysis that O₃ radiative bias is the largest term in five of eight models. The anti-correlation suggests a radiative compensation between the 9.6 µm band and the other parts of the OLR, assuming a constant globally integrated OLR at TOA. More interestingly, a strong correlation (R=0.9) is found between T_s radiative bias and broadband OLR (Fig. 13c) because the T_s affects the entire baseline of the outgoing radiance and its radiative effect plays the same role in the O₃ band as in the entire OLR. However, there is no significant correlation found between T_a radiative bias and OLR, likely because there is no coherent bias in T_a radiative effect between the O₃ band and in the entire OLR. There is neither a correlation between the H₂O radiative bias and broadband OLR. H₂O absorption is ubiquitous in the OLR. Consequently, biases in the 9.6 µm band do not drive the magnitude of the overall H₂O absorption in spite of the H₂O biases.

The anti-correlation between the biases in the 9.6 µm band and in the entire OLR band would suggest some bias drivers in the 9.6 µm band must play different roles at the other part of the OLR band. The further investigation of these processes would help to explain the radiative effect of different biases on the OLR estimations from models (Huang et al., 2008, 2014).

7 Conclusions

We have demonstrated a new method to quantitatively attribute the biases in O₃ band TOA flux from chemistry–climate model ensembles to O₃, H₂O, T_a, and T_s radiative components without cloud effect using observationally constrained IRKs in the clear sky. The study also provides the first vertically and globally resolved view of the radiative bias for each component. An IRK depicts the sensitivity of TOA fluxes to the vertical distribution of the geophysical quantities, such as O₃, H₂O, T_a, and T_s. While the products of 9.6 µm O₃ band IRK for O₃, H₂O, T_a, and T_s have been developed with the satellite observations by Aura TES, the record could be extended by MetOp-IASI and SNPP-CrIS Fourier Transform Spectrometer (FTS) measurements. We compute the model biases against reanalysis data for four key variables: O₃, H₂O, T_a, and T_s. Especially for O₃ biases, the newly developed TCR-1 O₃ assimilation data (Miyazaki et al., 2015; Miyazaki and Bowman, 2017) are, for the first time, used as the state-of-the-art benchmark for tropospheric O₃ in models. These specific bias comparisons for the CCMI study cause the modelers to investigate the reasons for these biases and motivate them to improve their simulations. For example, MRI-ESM1r1 shows the reduced LNO_x emission help to improve their tropical upper tropospheric O₃ and its radiative bias.

O₃ abundance is found to be the dominant driver for the ensemble flux bias. Tropical tropospheric O₃ is too high for most models and accounts for about 70 % of the flux bias (Table 2). The second driver in the model ensemble becomes the T_s instead of H₂O because the T_s radiative components are commonly biased low in the model ensemble, while the H₂O radiative biases between models are biased randomly in both directions with large diversity. For individual models, however, H₂O is the second most important driver, a larger component than T_s, for many cases such as AM3, SOCOL3, and MRI-ESM1r1.

In addition to determining that the tropospheric O₃ and H₂O are overestimated, and the surface is too cold, the study also tells us the geolocations, in latitudes and altitudes, of the deviations in these geophysical quantities that propagate into the flux bias.

The largest spread of the flux bias between the models is found in the tropics. The principal contributors governing each model are different and controlled by different processes over different regions. The flux biases in five of the eight models (AM3, SOCOL3, EMAC-L47MA, EMAC-L90MA, and MRI-ESM1r1) are primarily driven by too much O₃ in the tropical middle and upper troposphere. H₂O is a big driver in five models (AM3, SOCOL3, GEOSCCM, CMAM, and CESM). T_s is an important contributor in CMAM in addition to its H₂O.

Although AM3 and CMAM overall have relative lower TOA flux biases globally, we found they are actually right for the wrong reasons. In AM3, the dominant positive H₂O radiative bias (87.7 mW m⁻² in Table 2) happens to be canceled by the dominant negative O₃ component (−140 mW m⁻²), while in CMAM, the large positive H₂O component (127.9 mW m⁻²) is mostly being compensated by T_s radiative bias (−100.2 mW m⁻²). The two relatively young models among the model ensembles, SOCOL3 and MRI-ESM1r1, have a large potential to be improved for their fluxes by reducing their strong negative radiative biases from both tropical upper tropospheric O₃ and tropical lower tropospheric H₂O.

On average, the model ensemble underestimates the flux by about 133 mW m⁻² due to overestimated tropical tropospheric O₃ and H₂O. The underestimate of the TOA flux implies the model atmosphere is too opaque. In a more opaque atmosphere, the change in flux will be weaker for the same change in tropospheric O₃ because the sensitivity (i.e., IRKs) is weaker. With such feedback, the O₃ RF, the changes in O₃ GHG effect from pre-industrial times to the present day would likely be underestimated. The opacity of the atmosphere is controlled by climate processes, such as the hydrological cycle, that is shown can indirectly affect the O₃ GHG effect and RF, as discussed in Kuai et al. (2017).

The spatially explicit and process-focused differences could be used as a basis for emergent constraints (Bowman et al., 2013). New techniques such as hierarchical emergent constraints (HECs) can harness this spatial information so that specific processes affecting O₃ RF can be identified (Bowman et al., 2018). Moreover, if this correlation exists between the TOA flux bias and the O₃ RF, then a similar issue could be found in the RF of other GHGs, such as CO₂ and CH₄. That is a subject for future research.

Finally, although the chemical reanalysis dataset provides comprehensive information on model radiative biases, we need to understand its performance. For instance, further improvements are still needed for lower tropospheric O₃ (Miyazaki and Bowman, 2017). Ingesting more datasets and applying a bias correction procedure would be useful to improve reanalysis accuracy. The lower tropospheric O₃ analysis would benefit from the recently developed satellite retrievals with high sensitivity to the lower troposphere (Fu et al., 2018) and the optimization of additional precursor emissions.

Appendix A

Here, we compare the MRI-EMS1r1 experiment run RefC1_50 % LNO_x with its RefC1 run. The new run's emission decreases by about 50 % compared with the original run. The global lightning NO_x emission annual mean in 2006 simulated in the experiment run is reduced from ∼10.79 TgN yr⁻¹ in RefC1 to ∼5.21 TgN yr⁻¹. The 10-year average changes from ∼10.44 to ∼5.18 TgN yr⁻¹.

We found the total flux bias is much reduced due to the improved O₃ radiative bias (Fig. A1a, b). As we expected, the vertical resolved O₃ radiative bias shows that the overestimation of the tropical upper tropospheric O₃ radiative bias is much weaker in the new run (Fig. A1c, d). This improvement is due to the lower O₃ biases in this region caused by reduced LNO_x emission (Fig. A1e, f).

Figure A1 The comparison of the MRI-ESM1r1 experiment of half LNO_x in total flux bias (a, b), O₃ radiative bias (c, d), and O₃ bias (e, f). (a, c, e) New run with half LNO_x. (b, d, f) RefC1 run.

We also see some changes in the latitudinal distributions of the H₂O radiative bias. This is because the reduction of the upper tropospheric O₃ will cause the model responses in the O₃ heating rate, which would have radiative effect on the temperature, atmospheric stabilities, and convective activity (Nowack et al., 2015). All these factors would impact water vapor and cloud formation.