Greenhouse gases Observing SATelite, retrieval version 02.81 (GOSAT_v2.81)
Orbiting Carbon Observatory, retrieval version 9r (OCO2_v9r)
The Orbiting Carbon Observatory - 2 (OCO-2) was launched on 04 July 2014 and retrievals are available since September 2014 (Crisp et al. 2017; Eldering et al. 2017). We have used bias-corrected retrievals, version 9r (OCO-2_v9r) (Updated Lite files from https://co2.jpl.nasa.gov/#mission=OCO-2) (O’Dell et al. 2018). OCO-2 data are available for XCO2, a priori profiles (20 levels) and column averaging kernels (20 pressure levels). Data with quality flag of 0 (best data) are used in the analysis. The data coverage of OCO-2 (approximately 10 million data points annually) is much greater than that of GOSAT (approximately 0.1 million data points annually). The single-shot accuracy of OCO-2 retrievals is estimated to be 1.3 ppm over land and 1 ppm over ocean (Kulawik et al., 2020). The better quality of the XCO2 data from OCO-2 (compared to bias uncorrected GOSAT_v2.81), although without an overlap with the MIROC-ES2L simulation, provides us an assessment of the retrieval data quality on evaluation of the CO2 flux models.
MIROC Earth System model version 2 for Long-term simulations (MIROC-ES2L) and evaluation using remote sensing observations
MIROC-ES2L is the latest generation of earth system model (ESM) that is developed at the JAMSTEC, in collaboration with Tokyo University and the National Institute for Environment Studies (NIES). The ESM has fully coupled terrestrial-ocean-atmosphere model components for simulating the evolutions of carbon-nitrogen cycles and climate for the past centuries and centuries to millennium-term future projection (Hajima et al., 2020),. The horizontal resolution of the atmosphere is T42 spectral truncation (approximately 2.8◦ intervals for latitude and longitude) and the vertical resolution is 40 layers up to 3 hPa with a hybrid σ –p coordinate; the ocean grid system consists of horizontally 360 × 256 grids with 62 vertical layers.
The terrestrial carbon cycle component covers major processes relevant to global carbon cycle, with vegetation (leaf, stem, and root), litter (leaf, stem, and root), and humus (active, intermediate, and passive) pools and with a static biome distribution. Photosynthesis or gross primary productivity (GPP) is simulated based on the Monsi–Saeki theory (Monsi and Saeki 1953), and the allocation of photosynthate between carbon pools in vegetation (e.g., leaf, stem, and root) is regulated dynamically following phenological stages. The transfer of vegetation carbon into litter–soil pools is simulated using constant turnover rates, and in deciduous forests, seasonal leaf shedding occurs at the end of the growing period. The perturbation from land-use change on terrestrial carbon cycle is simulated by the transition of 5 types of land-use tiles (primary vegetation, secondary vegetation, urban, crop, and pasture) in a grid; e.g., deforestation processes are represented by the reduction in area fraction of primary/secondary vegetation (and increase of less-vegetated land-use tile like urban), and deforested carbon is gradually lost to the atmosphere. Details on carbon cycle processes in the model can be found elsewhere (Ito and Oikawa 2002; Hajima et al. 2020).
In the Ocean ecosystem component Embedded within the ocean circulation model (OECO2), ocean biogeochemical dynamics are simulated with 13 biogeochemical tracers (three types of nutrients, four biological tracers, four carbon and/or calcium, and other two inorganic tracers). The carbon cycle processes are simulated following the biogeochemical dynamics of the tracers, assuming that all organic materials have identical elemental stoichiometric ratio. Ocean ecosystem dynamics are simulated based on the nutrient cycles of nitrate, phosphorous, and iron. The nutrient concentration, in conjunction with the controls of seawater temperature and the availability of light, regulates the primary productivity of the phytoplankton. In the model, zooplankton is assumed to be independent of abiotic conditions (e.g., seawater temperature) and dependent on biotic conditions (phytoplankton and zooplankton concentrations). Detritus contains nitrate, phosphorus, iron, and carbon, most of which is remineralized while sinking downward. A fraction of the detritus that reaches the ocean floor is removed from the system. The formulations of atmosphere–ocean gas exchange, carbon chemistry, and related parameters follow protocols from the Ocean Model Intercomparison Project (OMIP) (Orr et al. 2017).
The basic performance of the land and ocean carbon cycle are assessed by comparing the simulated spatial patterns with observation data for land (gross primary production, vegetation carbon and soil carbon for the land component) and ocean (primary production, export production, atmosphere-ocean CO2 flux, dissolved inorganic carbon) (Hajima et al. 2020). Hajima et al. (2020) also showed that atmosphere-land or ocean CO2 fluxes in cumulative value (i.e., the changes in land or ocean amount) in the period of 1850-2014 are within the estimated range suggested by Le Quéré et al. (2018), although the GPP seasonality in land tropics and CO2 flux seasonality in the Southern Ocean and the North Atlantic Ocean requires improvement. The feedback strengths of carbon cycle processes, which is associated with the magnitude of natural carbon sinks in response to atmospheric CO2 concentration and climate change, are also quantified and compared with other ESMs, showing MIROC-ES2L has the intermediate values for the feedback parameters (Arora et al. 2020).
The original historical simulation of MIROC-ES2L is obtained from the CMIP6 exercise (Eyring et al. 2016; Hajima et al. 2020), where the model is driven by anthropogenic fossil fuel CO2 emission, other GHGs and aerosol emissions/concentrations, natural forcing like solar irradiance change and volcanic events, and land-use change. The CO2 concentration is explicitly simulated by the model based on the external forcing of anthropogenic CO2 emission, simulated land and ocean CO2 sink, and simulated LUC-derived CO2 emission (the experimental configuration in which CO2 concentration is prognostically simulated by models is called "emission-driven" mode, while the configuration in which CO2 concentration is prescribed as forcing data is called “concentration-driven” mode) (Jones et al., 2016). The simulation period was 1850-2014, following a continuous 2483 spin-up years by the “concentration-driven” mode and 691 years by the “emission-driven” mode. The accuracy of future prediction by the ESMs is likely to depend on the quality of the historical simulation of the carbon-nitrogen-climate feedbacks and the model response to external forcing (e.g., terrestrial response to land-change forcing). However, the sparse in situ measurements and smaller flux footprint of the surface sites limited the evaluation of the historical simulations by the ESMs. The situation has now changed because of the global coverage of XCO2 measurements from space (GOSAT and OCO-2), which are explored in this study for an evaluation of the MIROC-ES2L historical simulation.
MIROC Atmospheric Chemistry-Transport model (MIROC4-ACTM)
The MIROC (version 4) – based atmospheric chemistry transport model (MIROC4-ACTM) simulations are conducted at horizontal resolution of T42 spectral truncations with 67 hybrid pressure layers in vertical. The greater number of vertical layers covering the altitude range from Earth’s surface to 0.0128 hPa or about 80 km, and closely spaced levels in the upper troposphere and lower stratosphere region is designed for better representation of the vertical profiles of long-lived gases in the Earth’s atmosphere. MIROC4-ACTM (Patra et al., 2018; Watanabe et al., 2011) simulated horizontal winds and temperature are nudged to the Japanese 55-year reanalysis or JRA55 (Kobayashi et al. 2015), is well validated for large-scale interhemispheric transport using sulphur hexafluoride (SF6) and the Brewer-Dobson circulation patterns using the age of air derived from CO2 and SF6 (Patra et al., 2018). Thus we expect better simulations of XCO2 by accounting for the accurate meridional and vertical profiles of CO2. However, the convective (weekly or shorter time scale) transport remained invalidated due to the lack of appropriate observational parameters, e.g., Radon-222 with radioactive decay lifetime of 3.8 days (Patra et al., 2018).
We use monthly-mean CO2 fluxes from JAMSTEC’s time-dependent inversion model for this analysis (Saeki and Patra 2017). The inverse model estimated CO2 fluxes (referred to as MIROC4_Inv) for 84 regions of the globe using the MIROC4-ACTM simulations and in situ CO2 observations at 30 sites around the globe (Figure 1). The data are provided by the National Oceanic and Atmospheric Administration (NOAA) Cooperative Air Sampling Network (Dlugokencky and Tans 2019), and Japan Meteorological Agency (JMA) (Tsuboi et al. 2013). A priori MIROC4-ACTM simulations are conducted using interannually varying emissions due to fossil-fuel and cement production (FFC) from the global carbon project (Friedlingstein et al. 2019), climatological monthly-mean terrestrial and oceanic fluxes from the Carnegie-Ames-Stanford Approach (CASA) land model and upscaling model of air-sea CO2 flux measurements, respectively (Randerson et al., 1997; Takahashi et al., 2009). Thus, the interannual variabilities in the land and ocean regional fluxes estimated by the inverse model are driven entirely by the signals in observed CO2 variabilities at 30 sites.
Preparation of model XCO2 data and processing of the gridded maps and time series
Here, in order to simultaneously evaluate both the XCO2 and CO2 fluxes simulated by MIROC-ES2L, we first made two types of MIROC4-ACTM simulations – one is driven by the CO2 fluxes of MIROC-ES2L and the other is by MIROC4_Inv flux; the simulations are referred to as ACTM_ES2LF and ACTM_InvF, respectively. Thus, the atmospheric CO2 transport, driven by JRA55 reanalysis, are common for both the simulations. This experimental setting allows us to make direct comparison of XCO2 in both simulations with that of satellite monitoring. The model simulations by MIROC4-ACTM are archived at hourly time resolution and sampled for the GOSAT and OCO-2 measurement locations by linear interpolations in horizonal grid and time (Patra et al., 2017). The model simulations are then convolved with satellite a priori (CO2prior) and column averaging kernels (Ai) after calculating to their respective pressure levels (Pi; where i is the index for retrieval pressure layers; 20 for OCO-2 and 15 for GOSAT) from the original ACTM vertical grids (67 hybrid). Thus the ACTM simulated XCO2 are
XCO2ACTM = ∑i (CO2priori . dPi) + ∑i Ai . dPi (∑i CO2ACTMi - ∑i CO2priori) (1)
where, dPi is the pressure weighting function.
Both the simulations using MIROC4-Inv and MIROC-ES2L fluxes are spin-up from 1996-2009, so that the CO2 concentration responses to the flux changes near the Earth’s surface propagates throughout the atmosphere, as high as the model top of about 90km. This is important because the XCO2 measurements capture the total columns and there are systematic CO2 gradients in the stratosphere and mesosphere that are season dependent, and lags more than 6 years for the concentration increase rate in the troposphere. Both the model concentrations are adjusted for initial values so that the global mean concentrations agree with the GOSAT observations on 01 January 2010. This is a valid approximation because we do not consider any chemical production or loss of CO2 in the models, thus the concentration increases are consistent with the flux evolution (Basu et al. 2018). The initial value adjustment helps to put all plots on an identical XCO2 scales.
All the XCO2 results are gridded to 2.5o x 2.5o latitude-longitude grid for further analysis (Patra et al., 2017). Time series analysis are conducted by further aggregating XCO2 data for the 15 partitions of the land and 11 ocean regions as marked by the maps in Figure 1. The model data differences and agreements are quantified by comparing the growth rates and seasonal cycles between the observed and simulated XCO2, following decomposition of the time series using a harmonic analysis (Thoning et al., 1989). We fitted 4 harmonics to each of the regionally aggregated time series for deriving a long-term trend and a fitted time series. The seasonal cycles are calculated as fitted – trend time series, and the growth rate is calculated as time derivative of the trend time series.
Leaf Area Index (LAI) and Gross Primary Productivity (GPP)
The LAI, one-sided green leaf area per unit horizontal ground surface (m2 m-2), is one of the main drivers of primary production in the ESMs. LAI defines the canopy structure and functions of the biosphere, influenced by the weather and climate of the region. Thus, the LAI and GPP together allow us to diagnose basic ecosystem functioning for carbon assimilation in MIROC-ES2L. We used LAI from the Moderate Resolution Imaging Spectroradiometer (MODIS) Collection 6 (Yan et al. 2016) and GPP from FLUXCOM, a CO2 flux upscaling inter-comparison project (Tramontana et al. 2016). The MODIS LAI product with original 500m resolution was aggregated globally based on the same quality control with Ichii et al. (2017) and averaged over the T42 spatial resolution for comparison with MIROC-ES2L results. The FLUXCOM GPP product (Jung et al. 2020) was generated by an ensemble of estimations by multiple machine-learning algorithms with eddy-covariance CO2 flux observation network (FLUXNET) and MODIS land products (Tramontana et al. 2016). The original 1/12 degree resolution data were converted to the T42 spatial resolution.