The impacts of uncertainties in emissions on aerosol data assimilation and short-term PM2.5 predictions in CMAQ v5.2.1 over East Asia

For the purpose of improving PM prediction skills in East Asia, we estimated a new background error covariance matrix (BEC) for aerosol data assimilation using surface PM2.5 observations that accounts for the uncertainties in 20 anthropogenic emissions. In contrast to the conventional method to estimate the BEC that uses perturbations in meteorological data, this method additionally considered the perturbations using two different emission inventories. The impacts of the new BEC were then tested for the prediction of surface PM2.5 over East Asia using Community Multi-scale Air Quality (CMAQ) initialized by three-dimensional variational method (3D-VAR). The surface PM2.5 data measured at 154 sites in South Korea and 1,535 sites in China were assimilated every six hours during the Korea-United States Air 25 Quality Study (KORUS-AQ) campaign period (1 May–14 June 2016). Data assimilation with our new BEC showed better agreement with the surface PM2.5 observations than that with the conventional method. Our method also showed closer agreement with the observations in 24-hour PM2.5 predictions with ~44 % fewer negative biases than the conventional method. We conclude that increased standard deviations, together with horizontal and vertical length scales in the new BEC, tend to improve the data assimilation and short-term predictions for the surface PM2.5. This paper also suggests further 30 research efforts devoted to estimating the BEC to improve PM2.5 predictions. https://doi.org/10.5194/gmd-2020-116 Preprint. Discussion started: 15 May 2020 c © Author(s) 2020. CC BY 4.0 License.


Introduction
Particulate matter (PM) affects human health and climate change (Brook et al., 2010;Chung et al., 2005Chung et al., , 2012Cohen et al., 2017). The impact of elevated levels of PM on human health, in particular, has become a central issue in Asia. One study reported that PM 2.5 causes 63 premature deaths per 100,000 population in Asia (Apte et al., 2015). According to the World 35 Health Organization (WHO), the mortality rate related to air quality in North Korea was 238 for every 100,000 people, the highest of 194 countries investigated (WHO, 2017). With regard to climate change, recent studies have revealed a significant direct and indirect aerosol effects and high PM levels on the radiation budget in Asia (Chung et al., 2010;Dong et al., 2016;Jung et al., 2019).
To investigate elevated PM conditions and their influences on health and climate in Asia, researchers frequently apply 40 chemical transport models (CTMs). CTM predictions of PM 2.5 and related parameters can show appreciable errors, as evidenced by comparisons between CTM outputs and surface observations (Choi et al., 2019;Denby et al., 2008;Song et al., 2012), aircraft measurements (Kim et al., 2013(Kim et al., , 2017bSouri et al., 2018), and space-based observations from low earth orbit (LEO) and geostationary (GEO) satellite sensors (Han et al., 2013(Han et al., , 2015(Han et al., , 2019Jung et al., 2019;Lee et al., 2016;Park et al., 2014, Park et al., 2011. The inaccuracy of CTM simulations has been associated with uncertainties in emissions of primary 45 air pollutants and meteorological fields as well as omissions of photochemical reactions occurring in chemical mechanisms (Han et al., 2013(Han et al., , 2015Kim et al., 2017a;Song et al., 2012).
To enhance the performance of CTM simulations, a number of researchers have constrained the initial conditions (ICs) in km × 15 km and 35 vertical sigma levels, respectively, with the model top at 50 hPa. We selected Statewide Air Pollution Research Center version 2007 (SAPRC07) gas-phase chemistry and AERO6 aerosol modules for this study. Other physical and chemical configurations for the WRF and CMAQ model simulations are described in Lee et al. (2016). Figure 1 shows the model domain, which covered areas of China, the Korean Peninsula, and Japan. The observation sites in China and Korea are highlighted in red and blue, respectively. All surface PM 2.5 observations from these sites were used in 100 the 3D-VAR DA. In particular, we focused on two regions: South Korea and East China designated by blue solid line for validation of the model performance. East China included the Beijing, Tianjin, and Hebei areas (known as the Jing-Jin-Ji area), where severe haze events frequently take place and dramatically affected the air quality of South Korea (Choi et al., 2019;Huang et al., 2018;Lee et al., 2017;Liu et al., 2017).
It should be noted that to construct a BEC that account for uncertainties in emissions, we used two anthropogenic emission 105 inventories. One was the Comprehensive Regional Emissions for Atmospheric Transport Experiments (CREATE) version 3.0 for the year 2015. This inventory specifically considers the economic and energy growth factors in China and Korea (Woo et al., 2012). The emissions were processed by the Sparse Matrix Operator Kernel Emissions (SMOKE) system developed by the United States Environmental Protection Agency (U.S. EPA). Details of the processes were described in Woo et al. (2012). The CREATE emissions with a resolution of 0.1° × 0.1° were converted into a 15 km × 15 km horizontal 110 resolution via a flux conserving interpolation technique.
The other inventory consisted of projected emissions from the Mosaic Asian anthropogenic emission inventory (MIX) for the year 2010 . This inventory based on a mosaic approach combines available emission inventories from various countries. The components of the MIX for South Korea, China, and Japan were the Clean Air Policy Support System (CAPSS) , the Multi-resolution Emission Inventory for China (MEIC), developed by Tsinghua University 115 (http://www.meicmodel.org), and the Regional Emission Inventory in Asia (REAS) version 2.1 (Kurokawa et al., 2013), respectively.
To account for changes in emissions between 2010 and 2016 over South Korea and China, we modified the MIX inventory.
Projection ratios for SO 2 , NO x , CO, VOC, and PM between 2010 and 2016 in China were based on a study by Zheng et al. (2018). The ratios for South Korea were based on changes in emission between 2010 and 2015 (the latest available emission 120 data for South Korea when this study was started). Information was obtained from the National Air Pollutant Emission report 2015, available at the National Air Pollutant Emission Service website (http://airemiss.nier.go.kr). Lacking information, we did not change emission data from other countries (e.g., Japan) inside the model domain. Using a flux conserving method, we converted the MIX emissions with 0.25° × 0.25° resolution into emission with a 15 km × 15 km resolution. Table 1 shows the annual emission fluxes of major air pollutants from the CREATE and MIX inventories over South Korea and 125 China. The MIX emissions show larger SO 2 (+3%), NOx (+10%), and PM 2.5 (+3%) fluxes, and fewer CO (-8%) and nonmethane volatile organic compound (NMVOC, -9%) emissions over China than the CREATE emissions. The differences in South Korea are relatively small, except for CO in the MIX emission inventory. The comparison between the emission https://doi.org/10.5194/gmd-2020-116 Preprint. Discussion started: 15 May 2020 c Author(s) 2020. CC BY 4.0 License.
inventories is demonstrated with selected species of SO 2 , NOx, and PM 2.5 over the domain (Fig. 2). In general, the CREATE emission rates are lower than the MIX emission rates, particularly over the East China. 130 We added identical biogenic and biomass burning emission data to both anthropogenic emissions data and applied the Model of Emission of Gases and Aerosols from Nature (MEGAN) v3.0 (Guenther et al., 2006) to model biogenic emissions.
To prepared input data for the MEGAN model, we used i) meteorology from the WRF model simulations, the configuration of which is described above; ii) leaf area index (LAI) data reprocessed with a 30 arc spatial resolution from eight day averaged Moderate Resolution Imaging Spectroradiometer (MODIS) data in 2016 (Yuan et al., 2011); and iii) the green 135 vegetable fraction (GVF) from the Visible Infrared Imaging Radiometer Suite (VIIRS) sensor.
The biomass burning emissions used in this study came from the Fire INventory from NCAR (FINN) v1.5 (Wiedinmyer et al., 2011). For the sake of reader convenience, we refer to the sum of biogenic, biomass burning, and two anthropogenic emissions datasets as the MIX and CREATE emission inventories.

GSI 3DVAR System 140
Figure 3 displays a schematic of the 24-hour PM 2.5 prediction system in this study. We employed Grid-point Statistical Interpolation (GSI) version 3.6 (Shao et al., 2016)

supported by the Developmental Testbed Center NCEP Environmental
Modeling Center for the 3D-VAR system. To update the aerosol initial conditions, we used surface PM 2.5 observations in South Korea and China measured at (00, 06, 12, and 18) UTC daily. The 3D-VAR algorithm optimizes assimilated fields (also known as analysis fields) with observations by iterative processes to minimize the cost function (J(x)) defined below: 145 species for PM CMAQ,k (i.e. n = 7 in Eq. (2)). The major species in PM CMAQ,i and PM CMAQ,j are the Aitken-and accumulationmode sulfate, nitrate, primary and secondary organic aerosols (OAs), elementary carbon (EC), and trace elements and those in PM CMAQ,k are mostly sea-salt and dust aerosols. We calculated the mass contributions of each aerosol species to PM 2.5 before the assimilation and then used them to allocate horizontal and vertical increments from the 3D-VAR assimilation to the CMAQ aerosol species. 165 The role of the BEC (i.e., B in Eq. (1)) is to determine how much observational information spreads to model grid points horizontally and vertically. The BEC comprises four successive matrices, all of which generally modelled independently (Descombes et al., 2015) as follow: In Eq. (3), U P represents the physical transform that defines the control variables related to the observations and their linear 170 relationships. In this study, we selected PM 2.5 as the control variable and used Eq. (2) to calculate it. Thus, U P is an identity matrix, S is a diagonal matrix representing model errors, and its diagonal components are the standard deviation of the control variables. U v is a matrix for vertical transformation, estimated by vertical correlations of the control variable (PM 2.5 in this study). Similarly, U h denotes the matrix of the horizontal transform defined by the horizontal correlation of PM 2.5 . All  Each simulation used either the CREATE or the MIX emissions along with the two meteorological fields denoted by "Met. 1" 180 and "Met. 2" in Fig. 4 and initialized at 06 and 18 UTC daily, respectively. WRF simulations for MET.1 and MET.2 were preformed using initial and boundary conditions obtained from the FNL analysis data. As mentioned in Sect. 1, this method is referred to as "time-lagged forecasts" (Hoffman and Kalnay, 1983), used in the NMC method.
After the simulation, we processed the results by the interface module "CMAQ2GENBE," developed in this study. The module extracted aerosol variables from the CMAQ model simulations, used Eq. (3) to calculate PM 2.5 fields, and converted 185 the PM 2.5 data into inputs for GEN-BEv2.0. In the GEN-BEv2.0, we combined the inputs into three cases: NMC CREATE, NMC MIX, and NEW NMC. NMC CREATE (the green box in Fig. 4) used two simulations (simulation 1 and simulation 2 in the WRF-CMAQ box) that used the same CREATE emissions, but different time-lagged meteorological fields. The NMC MIX case (the blue box in Fig. 4) was the same as the NMC CREATE case, except it used MIX emission. Using all of the simulation results to estimate the BEC parameters, the case of NEW NMC (the red box in Fig. 4) combined the uncertainties 190 of both meteorology and emissions. This approach is similar to that of a previous study by Kumar et al. (2019) for satelliteretrieved AOD DA, who added a perturbation derived from a comparison of various global and regional inventories to anthropogenic emissions. Section 3.3 compares the statistical parameters among the BECs.
Finally, we applied 3D-VAR to the analysis fields using PM 2.5 observations and background fields via the GSI program. We selected the background fields from the first simulation (Met. 1 + CREATE emissions) at assimilation times 00, 06, 12, and 195 18 UTC daily and used the "BKOBS2GSI" module to convert the background fields and observations into binary format files available for the GSI. We analyzed only two DA cases with NEW NMC BEC and NMC CREATE BEC, to clarify the comparisons clear in Sect. 3.

PM 2.5 Observations
During the simulation period May 1 to 14 June, 2016, we acquired hourly surface PM 2.5 from the Korean public air quality Since the observation data were collected without any information on measurement errors and biases, we applied an 205 additional quality control process from previous studies (Jiang et al., 2013;Peng et al., 2018;Schwartz et al., 2012). First, we assumed that observations showing concentrations smaller than 6 μg m -3 or larger than 800 μg m -3 were unrealistic in our model domain, so we discarded them. Than we applied the "buddy test," which is typically used to eliminate data anomalies in meteorological observations as follows: where O h is the 1-hour averaged PM 2.5 measured at time h (hour), O h-1 and O h+1 are observations measured before and after 1 hour in the same locations, and a and b are the empirical parameters, 50 μg m -3 and 0.15, respectively. We applied the second criterion only when O h-1 and O h+1 were measured. As mentioned in Sect. 2.2, in the observation error covariance matrix (R) in Eq. (1), the diagonal terms are defined as (50 μg m -3 + 0.0075 × PM 2.5 ) 2 (Schwartz et al., 2012), and all the off-diagonal components are zero under the assumption that all errors in PM 2.5 measurements are mutually independent. We used quality-215 controlled surface PM 2.5 observations from 154 sites in South Korea and 1,524 sites in China for the DA every six hours, took measurements at all the Korean sites and Eastern Chinese sites (542 sites), and used them for validation studies.

Results and Discussion
This section discusses the impact of various BEC matrices on the 3D-VAR DA and the capability of the BEC estimation method to predict surface PM 2.5 over northeast Asia. Section 3.1 compares the parameters of the BEC calculated from the 220 conventional NMC method and those from the NEW NMC method. To estimate the influence of the two, we validate the assimilated aerosol initial conditions against surface measurement data. Finally, we discuss the 24-hour prediction skills of surface PM 2.5 and changes in the PM 2.5 vertical profiles for each DA case.
https://doi.org/10.5194/gmd-2020-116 Preprint. Discussion started: 15 May 2020 c Author(s) 2020. CC BY 4.0 License. Comparing the two NMC methods, we found remarkable differences below 850 hPa, which is about 1 km above the surface. coarse-mode dust species had a similar vertical patterns compared to those of our NEW NMC BEC. The characteristics of the vertical and horizontal length scales, however, have not been fully explained in this study, thus requiring future investigation.

BEC Parameters from the NEW NMC and Conventional NMC Methods
In summary, the NEW NMC BEC shows larger standard deviations and larger vertical and horizontal length scales than the 260 other two conventional NMC BECs, particularly at low altitude. Therefore, the analysis fields determined by the PM 2.5 3D-VAR process are closer to the observations and the increments expand farther both vertically and horizontally. The following sections further discuss the influences of these specifications in the BECs on analysis fields and prediction skills.

Impact of the New Background Error Covariance matrix on Surface PM 2.5
To assess the impact of 3D-VAR with various BECs on surface PM 2.5 analysis fields, we used the NEW NMC BEC and the 265 NMC CREATE BEC. Hereafter, CTL (control), NMC, and NEW NMC runs will represent simulations without DA, simulation with DA using the NMC CREATE BEC, and simulation with DA using NEW NMC BEC, respectively. To focus on the impact of the changes in BECs, we conducted the CTL and other DA runs with the CREATE inventory, which was used for estimating NMC CREATE BEC, and identical meteorology (i.e., Met. 1), shown in Fig. 4. Figure 6 and 7 show the performances of the CTL, the conventional NMC, and the NEW NMC runs for metrics: index of 270 agreement (IOA) (which represents both errors and biases of analysis (Willmott, 1981)); mean bias (MB); and root mean square error (RMSE). We sampled the CMAQ PM 2.5 data from the simulations at 00, 06, 12, and 18 UTC daily between 1 May and 14 Jun 2016, and matched those to hourly model results at each location. Figures 6 and 7 focus on South Korea and East Asia, respectively, and total map of the metrics are presented in the supporting information (Fig. S1). The CTL run showed relatively higher IOAs at the sites in South Korea (a mean value of 0.59) than those in the eastern part of China (a 275 mean value of 0.47). The negative biases of the CTL run over South Korea (-7.38 μg m -3 ) are relatively larger than those over East China (-4.14 μg m -3 ).
The conventional NMC and NEW NMC produced significantly lower errors and biases for surface PM 2.5 . For example, compared to the averages from the results of the CTL runs, IOA, RMSE, and MB in the NEW NMC were improved by 56%, 59%, and 85% over South Korea, and 62%, 40%, and 22% over East China, respectively. To highlight such improvements, 280 "error reductions" (ER) as the ratio of the differences between the RMSEs calculated from the analysis (the conventional NMC or NEW NMC) and those from the CTL run to the RMSEs in the CTL run. The ERs for South Korea and East China presented in Fig. 8 indicate that the NEW NMC run produces more accurate prediction, particularly over South Korea ( Fig   S2 for total map of ERs). Table 2 summarizes the statistical performance metrics over South Korea and East China.
Collectively, we found that the predictions of the NEW NMC were more accurate than those of the conventional NMC. This 285 enhanced performance of the NEW NMC also improve its short-term prediction skills of surface PM 2.5 over South Korea and East China, the details of which are discussed in the following section.

Impacts of New Background Error Covariance Matrix on 24-hour PM 2.5 predictions
The performance of the 24-hour PM 2.5 predictions with and without DA was compared with surface PM 2.5 observations.  Table 3 also summarizes these statistical performance metrics at 6-hour time intervals.
Both the conventional NMC and NEW NMC runs showed improved model performances. We found the most influential time window to be from (+1H to + 6H). In the case of South Korea, the new NMC run, compared to the CTL run, improved 295 IOA, R, RMSE, and MB by 24%, 53%, 24%, and 71 %, respectively. In the case of East China, the NEW NMC run improved IOA, R, RMSE, and MB by 26%, 107%, 20%, and 2%, respectively. Similar to the previous discussion in Sect 3.2, the NEW NMC run performed better than the conventional NMC run, showing an increased IOA (+7 %) and R (+12 %), and a decreased RMSE (-7%) and MB (-44%) over South Korea.
Interestingly we also found that the improvements of MBs in the NEW NMC run lasted longer than 24 hours while other 300 improvements, such as those in the IOA, R, and the RMSE diminished with time. Improvements in the MBs for (+7H -+12H), (+13H -+18H), and (+19H -+24H) were 92%, 73%, and 89% over South Korea, and 83%, 93%, and 93% over East China, respectively. Improvements for the NMC run with the same time intervals were only 16%, 0%, and 0% over South Korea, and 6%, 9%, and 9% over East China.

Conclusion 320
In this study, we developed a new method of estimating the BEC for PM 2.5 DA, accounting for the uncertainty in emissions, which has not been addressed by the conventional DA of surface PM 2.5 observations. To account for such emission uncertainties in the BEC, the current study utilized not only time-lagged meteorological fields but also two versions of independent anthropogenic emission inventories. In an assessment of the impact of the new BEC on PM 2.5 analysis fields and the short-term prediction of surface PM 2.5 in the 3D-VAR framework, we designed and carried out three experiments: CTM 325 model simulations without DA, with DA using the conventional NMC BEC, and with DA using the NEW NMC BECs, over the period 1 May to 14 June 2016 in East Asia.
We found that the new approach exhibited a tendency to generate substantially increased standard deviations, vertical length scales, and horizontal length scales in the BEC. Such increases occurred particularly near the surface (below 700 hPa).
Subsequently, the use of the NEW NMC BEC positively impacted both the performance of the DA and the 24-hour 330 predictions of PM 2.5 and significantly reduced negative biases of the PM 2.5 predictions were (by as much as ~90%). Thus we conclude that the improvements were the result of the improved standard deviation and increased horizontal and vertical length scales in the NEW NMC BEC.
Based on the findings from this study, as well as the recent efforts in numerical weather prediction (NWP), we suggest two directions of research that will contribute to the construction of a more robust and sophisticated BEC matrix. First, the 335 statistical parameters in the BEC matrix can be further improved by optimization that accounts for the spatial and temporal characteristics of CTM background errors. The results of the CTL simulations (without DA) and the DA simulations showed that the performance of PM 2.5 model simulations over northeast China were poorer than that of simulations over South Korea.
The background errors over China, specifically the standard deviation in the BEC matrix, should be larger than those over South Korea. However, we did not address this issue in depth. Recent NWP DA studies have shown that the parameters in 340 the BEC matrix can be tuned and optimized region-by-region while accounting for the characteristics of background errors (Choi et al., 2017;Song et al., 2018). Additional effort should be devoted to estimating a more sophisticated BEC with regional optimization. The impact of an optimized BEC on PM 2.5 predictions can then be tested.
Another topic of research would be the application of a hybrid method recently adapted in an NWP system for estimating the BEC matrix (Kwon et al., 2018;Massart, 2018;Song et al., 2018) to the chemical DA. This method includes a weighted sum 345 of a static BEC estimated by the conventional NMC method and a dynamic BEC estimated by daily ensemble predictions. It has been reported that the hybrid BEC method improves the NMC BEC by mimicking daily changes in background errors.
In the construction of a hybrid BEC matrix in the chemical DA, perturbations in each ensemble run may be made created considering not only emission uncertainties but also other major uncertainties in CTM simulations, such as errors of aerosol parameterization, initial and boundary conditions, and chemistry mechanisms. All of topics of research will be tested in the 350 series of future studies.

Competing interests
The authors declare that they have no conflict of interest. 365 Collins, W. D., Rasch, P. J., Eaton, B. E., Khattatov, B. V., Lamarque, J.-F. and Zender, C. S.: Simulating aerosols using a chemical transport model with assimilation of satellite aerosol retrievals: Methodology for INDOEX, J. Geophys. Res.