3DVAR Aerosol Data Assimilation and Evaluation Using Surface PM_2.5, Himawari-8 AOD and CALIPSO Profile Observations in the North China

Abstract

Based on the Model for Simulating Aerosol Interactions and Chemistry (MOSAIC) aerosol scheme of the Weather Research and Forecasting model coupled with online Chemistry (WRF-Chem) and the three-dimensional variational (3DVAR) assimilation method, a 3DVAR data assimilation (DA) system for aerosol optical depth (AOD) and aerosol concentration observations was developed. A case study on assimilating the Himawari-8 satellite AOD and/or fine particulate matter (PM_2.5) was conducted to investigate the improvement of DA on the analysis accuracy and forecast skills of the spatial distribution characteristics of aerosols, especially in the vertical dimension. The aerosol extinction coefficient (AEC) profile data from The Cloud-Aerosol Lidar and Infrared Pathfinder Satellite Observation (CALIPSO), surface PM_2.5 and Himawari-8 AOD measurements were used for verification. One control experiment (without DA) and two DA experiments including a PM_2.5 DA experiment denoted by Da_PM and a combined PM_2.5 and AOD DA experiment denoted by Da_AOD_PM were conducted. Both DA experiments had positive effects on the surface PM_2.5 mass concentration forecast skills for more than 60 h. However, the Da_PM showed a slight improvement in the analysis accuracy of the AOD distribution compared with the control experiment, while the Da_AOD_PM showed a considerable improvement. The Da_AOD_PM had the best positive effect on the AOD forecast skills. The correlation coefficient (CORR), root mean square error (RMSE), and mean fraction error (MFE) of the 24 h AOD forecasts for the Da_AOD_PM were 0.73, 0.38, and 0.54, which are 0.09 (14.06%), 0.08 (17.39%), and 0.22 (28.95%) better than that of the control experiment, and 0.05 (7.35%), 0.06 (13.64%), and 0.19 (26.03%) better than that of the Da_PM, respectively. Moreover, improved performance for the Da_AOD_PM occurred when the AEC profile was used for verification, as when the AOD was used for verification. The Da_AOD_PM successfully simulated the first increasing and then decreasing trend of the aerosol extinction coefficients below 1 km, while neither the control nor the Da_PM did. This indicates that assimilating AOD can effectively improve the analyses and forecast accuracy of the aerosol structure in both the horizontal and vertical dimensions, thereby compensating for the limitations associated with assimilating traditional surface aerosol observations alone.

1. Introduction

Aerosol pollution has a wide range of effects on regional air quality, weather, climate, and human health [1,2,3,4,5,6,7]. To achieve precise air pollution prevention and control and to carry out research on the impact of aerosols on weather, climate, and human health, it is of great scientific importance to use air quality models to predict and analyze the spatiotemporal evolution of aerosols, as these have a high application value [2,8,9,10,11,12,13,14,15,16]. However, the forecast skill of the air quality model still needs to be improved owing to uncertainties in aerosol source emissions, imperfections in its physical and chemical parameterization process, and uncertainties in the initial conditions (ICs) for its chemistry and meteorology. Aerosol data assimilation (DA) technology can effectively reduce the uncertainties of the chemistry ICs for the model and then has been generally utilized to improve aerosol pollution forecasting [17,18,19,20,21,22,23,24,25,26,27,28].

Early aerosol DA studies are mostly performed with offline chemical transport models and the mass concentrations of gaseous pollutants (such as O₃) and particulate matter (such as PM_2.5) [29,30,31,32]. With the continuous development of atmospheric chemistry model technology and the improvement in computation power, the air quality model has developed from an offline model to a meteorological–chemical coupled model [9]. The DA for an online coupled model, such as WRF-Chem model, has also gradually developed. Meanwhile, the assimilated state variables are no longer limited to a single state variable, but can reflect comprehensive analyses of species and size distributions in an aerosol scheme. Li et al. [18] developed a three-dimensional variational DA system for the MOSAIC aerosol scheme in WRF-Chem. Then, the surface PM_2.5 and speciated concentration observations were assimilated and evaluated based on this DA system. The results showed that DA has a positive effect on both ICs and PM_2.5 forecasts for up to 24 h. With the increasing abundance of aerosol observational data, remote sensing data as well as conventional mass concentration measurements are assimilated based on different assimilation schemes and have significantly improved the forecast skills for aerosols [19,20,21,22,23,24,25]. Liang et al. [26] assimilated AEC profiles from five lidars in China and found that the root mean square error of surface PM_2.5 in the initial field for the model was reduced by 10.5 μg m⁻³ (17.6%). They also found that a larger reduction occurred when the AEC and PM_2.5 data were assimilated simultaneously. Ye et al. [33] assimilated the AEC from the CALIPSO during two pollution processes in North China and showed that DA can effectively reduce the simulation error of the three-dimensional structure of aerosols and improve the forecasting skills for pollutants. Although lidars on ground or onboard satellites can provide high-quality AEC vertical distribution, the horizontal spatial distribution of lidar data is very confined and the device maintenance is at high cost; therefore, there are many limitations to improving operational air quality forecasting by assimilating AEC data. Satellite derived AOD is the integral of the AEC of the entire atmosphere in the vertical direction. Compared with AEC data, AOD products derived from satellites, especially geostationary satellites, have high temporal and spatial resolutions and coverage; hence, they can compensate for the lack of resolution and coverage of lidar data. Therefore, research on satellite AOD assimilation has received considerable attention. Liu et al. [17] developed an AOD assimilation algorithm within the National Centers for Environmental Prediction (NCEP) Gridpoint Statistical Interpolation (GSI) 3DVAR DA system extended with the Community Radiative Transfer Model (CRTM) to compute the AOD. By implementing the GSI-3DVAR DA system, assimilation of the Moderate Resolution Imaging Spectroradiometer (MODIS) AOD substantially improved aerosol analyses and subsequent forecasts. Similar improvements were found in studies using the GSI-3DVAR DA algorithm [19,23,34,35,36,37]. However, the GSI-3DVAR method is based on the Goddard Chemistry Aerosol Radiation and Transport (GOCART) aerosol scheme, which finely describes dust and sea salt components but lacks nitrate and ammonium salt components and the corresponding secondary reaction process. In contrast, the MOSAIC scheme treats aerosol as eight individual species, these species include black carbon, organic carbon, nitrate, sulfate, chloride, ammonium, sodium, and other unspecified inorganic species, such as silica, other inert minerals, and trace metals. A sectional approach is employed to represent size distributions of each species, and thus such a comprehensive description of aerosol species and size distributions is more accurate. Additionally, the MOSAIC scheme contains more processes such as heterogeneous reactions of gaseous precursors taking place under high humidity. Thus, MOSAIC is more suitable for expressing complex aerosol components in heavily polluted areas in China [38]. A 3DVAR DA algorithm based on the MOSAIC that could assimilate surface PM_2.5 and speciated concentration observations was developed by Li et al. [18] and has been extended to assimilate aerosol optical properties such as AEC and AOD [20,25,26,27]. Sun et al. [39] established an interface for MOSAIC and conducted experiments to assimilate ground PM_2.5, PM₁₀, SO₂, NO₂, O₃, and CO. However, research on the assimilation of multi-source aerosol observations, especially for aerosol optical data for the MOSAIC scheme is still underway.

This paper presents a case study of Himawari-8 (H-8) satellite AOD assimilation using the 3DVAR DA system developed by Wang et al. [27]. The AOD and PM_2.5 were assimilated alone or simultaneously during the pollution process in Northern China in January 2019. Moreover, CALIPSO AEC profiles, as well as AOD and PM_2.5, were used to verify the improvement effect of DA on the accuracy of the aerosol vertical structure and PM_2.5 forecasting skills.

The remainder of this paper is arranged as follows. Section 2 gives a brief description of the data, WRF-Chem configurations, 3DVAR DA system, and experimental design. In Section 3, the improvement effect of DA on the accuracy of the aerosol vertical structure and PM_2.5 forecast skills is verified by PM_2.5, AOD, and CALIPSO AEC profiles. The discussion and conclusions are presented in Section 4.

2. Materials and Methods

2.1. WRF-Chem Model and Data

WRF-Chem Version 4.1 was employed to simulate aerosol pollution. The experimental domain was double-nested, the first domain (d01) has 164 × 155 grid points with a horizontal resolution of 27 km, and the second domain (d02) has 175 × 166 grid points with a horizontal resolution of 9 km covering the central and eastern regions of China (Figure 1a). There were 40 vertical layers and the resolution gradually decreased from bottom to up. The model configuration was as follows: a WRF single-moment 5-class microphysical scheme, rapid radiative transfer model for general circulation models (RRTMG) shortwave radiation scheme, RRTMG longwave radiation scheme, Noah land surface model, Yonsei University (YSU) boundary layer scheme, Grell-3D cumulus parameterization, revised MM5 Monin-Obukhov near-surface layer scheme, carbon-bond mechanism version Z (CBMZ) chemical reaction mechanism, fast-J photolysis calculation scheme, and MOSAIC_4bin aerosol scheme. There are four particle-size bins (4bin), namely, 0.0390625–0.15625, 0.15625–0.625, 0.625–2.5, and 2.5–10 µm, for each of the eight aerosol types in MOSAIC_4bin.

Three types of aerosol observational data, including PM_2.5 mass concentration, H-8 AOD, and CALIPSO AEC, were used for assimilation and evaluation, and their spatial distributions are shown in Figure 1b. PM_2.5 and AOD data were used for assimilation and all PM_2.5, AOD, and AEC data were used for evaluation.

Hourly PM data released from the China National Environmental Monitoring Center (CNEMC) (http://www.cnemc.cn, last access: 5 May 2022) were used for assimilation and evaluation. To date, more than 2000 measurement sites were established in China. There are a total of 683 air quality monitoring sites in our studied region (d02). The quality control and preprocessing treatments for PM were similar to those described by Wang et al. [25]. As part of the main processes, PM_2.5 mass concentration exceeding 600 µg m⁻³ and less than 0 µg m⁻³ were removed, and the original measurements falling into the same model grid were averaged.

AEC profile data were obtained from Level 2 (version 4-2.0) AEC retrievals from CALIOP at 532 nm (https://www-calipso.larc.nasa.gov, last accessed: 5 May 2022). The horizontal resolution of the data was 5 km. Each profile had 399 layers, with a vertical resolution of 60 m at an altitude of −0.5~20.2 km and of 180 m at an altitude of 20.2~30.1 km. For the preprocessing and quality control of the CALIPSO data, please refer to Ye et al. [33].

AOD data were obtained from Level 2 AOD retrievals from Himawari-8 at 500 nm (http://www.jma.go.jp/jma/jma-eng/satellite/index.html, last accessed: 5 May 2022). The data has a horizontal resolution of 0.05 × 0.05. Limited by sunlight availability, AOD can only be obtained when the detection region is clear sky conditions. For the preprocessing and quality control of the CALIPSO data, please refer to Ye et al. [33].

2.2. 3DVAR DA System

The 3DVAR DA method minimizes the cost function J(x), which measures the distance of the state vector from the model background and the observations. By calculating the minimal value of the function appealing to the variational method, an “optimal” analysis field in terms of the minimum analysis error variance was obtained. From the theory of optimal estimation, 3DVAR transforms the problem of constructing the analysis field with minimum error into the problem of solving the minimum of the cost function $J\left(x\right)$. The basic formula of 3DVAR is as follows:

$$J\left(x\right)=\frac{1}{2}{(x-x_b)}^T{B}^{-1}\left(x-x_b\right)+\frac{1}{2}{(Hx-y)}^T{R}^{-1}\left(Hx-y\right)$$

(1)

Incremental approach was adopted to minimize the cost function J. The incremental cost function for the 3DVAR was calculated as follows:

$$J\left(\delta x\right)=\frac{1}{2}\delta x^T{B}^{-1}\delta x+\frac{1}{2}{\left(H\delta x-d\right)}^T{R}^{-1}\left(H\delta x-d\right)$$

(2)

Here, $\delta x$ is the incremental state variable defined as $\delta x=x-x_b$, where x is the state vector and $x_b$ is the background vector. $B$ is the background error covariance matrix associated with $x_b$ and $R$ is the observation error covariance matrix, which determine the errors of observation and model information, respectively. $B$ is related to the state variables that are from the model variables, while $R$ is related to the observed variables. The $d=y-H{x}_b$ is the observation innovation vector, where $y$ is an observation vector and $H$ is the observation operator that computes the modeled observation estimates from the state variables. The satellite derived AOD and PM_2.5 observation errors are determined by measurement and representation errors [34,40]. In this study, building on the earlier work of Yumimoto et al. [21], the observation error of Himawari-8 AOD is set to 0.06. According to Wang et al. [25], Conventional observational errors were set to half the background error standard deviations (SDs) of each control variable.

H is the observation operator, which maps the model state variables at each grid cell to the observation space. In this study, the operator for PM_2.5 mass concentrations is easily performed only using a simple linear operator (combining the first three particle-size bins of each aerosol species within MOSAIC_4bin). In order to improve computing efficiency, the Mie scattering theory was utilized as the observation operator of AOD. It should be noted that the meter scattering theory is based on the assumption of spherical particles, and the AOD calculated for non-spherical or non-compact particles may be inaccurate. Appealing to the Mie scattering calculation, optical efficiencies such as extinction efficiency can be obtained after determining particle size parameter and averaged refractive index that are computed from the state variables. Then, the extinction coefficient is the summation of extinction efficiencies over four size bins. Finally, AOD is the column summation of extinction coefficient over the vertical layers, and so this is the forward AOD operator. Note that it is very inefficient to compute aerosol optical efficiencies using the complex Mie scattering calculation; nevertheless, the optical properties module in WRF-Chem employs a polynomial expansion approximation to compute them, which is more efficient. Thus, the study also utilizes this approach. The tangent linear (TL) and adjoint operators are necessary for developing the DA system, which are described in detail in our previous study [27].

B is important for determining the performance of the DA process, which transfers observation information to the model grid cells. For high-resolution numerical models, the B matrix is huge and difficult to store and find an inverse. Therefore, a simplification of the B-matrix is necessary for numerical calculations. We followed the method used by Kalnay [41] and Li et al. [18]. Specifically, the $B$ matrix is decomposed into $B={\mathrm{DCD}}^{\mathrm{T}}$, where $\mathrm{D}$ is the background error standard deviation matrix, which determines the analysis increment magnitude, a larger background error standard deviation will generate a larger increment. $\mathrm{C}$ is the background error correlation coefficient matrix, which consists of a positive definite symmetric matrix of the horizontal and vertical error correlations and determines the analysis increment scope. For the detailed calculation of B matrix, please refer to Wang et al. [27].

2.3. Experimental Design and Evaluation Method

A control experiment (control) and two sets of DA experiments were performed to evaluate the effects of the DA. The control experiment did not assimilate observation data and 60 h forecasts were produced starting at 0500 UTC on 12 January 2019. One DA experiment (Da_PM) assimilated PM_2.5 data at 0500 UTC on 12 January 2019, and used the DA analysis field as the initial chemical field. The other DA experiment (Da_AOD_PM) simultaneously assimilated PM_2.5 and AOD data. In addition to the chemical ICs, all experiments had the same configuration. The National Centers for Environmental Prediction (NCEP) reanalysis data over a 6 h interval and 0.25-degree by 0.25-degree resolution is used to generate the initial and boundary conditions of meteorological fields.

The correlation coefficient (CORR), mean fraction error (MFE), and root mean square error (RMSE) were used to evaluate the analyses and forecast accuracy of the aerosols in the experiments. A larger CORR and smaller MFE or RMSE indicate better performance.

$$CORR=\frac{{{\displaystyle \sum }}_{i=1}^N\left(M_i-\overline{M}\right)·\left(O_i-\overline{O}\right)}{\sqrt{{{\displaystyle \sum }}_{i=1}^N\left(M_i-\overline{M}{)}^2\right)}·\sqrt{{{\displaystyle \sum }}_{i=1}^N{(O_i-\overline{O})}^2}}$$

(3)

$$MFE=\frac{1}{N}{{\displaystyle \sum }}_{i=1}^N\frac{\left|M_i-O_i\right|}{\left(M_i+O_i\right)/2}$$

(4)

$$RMSE=\sqrt{\frac{1}{N}{{\displaystyle \sum }}_{i=1}^N{(M_i-O_i)}^2}$$

(5)

$$\overline{M}=\frac{1}{N}{{\displaystyle \sum }}_{i=1}^N{M}_i$$

(6)

$$\overline{O}=\frac{1}{N}{{\displaystyle \sum }}_{i=1}^N{O}_i$$

(7)

where $M_i$ and $O_i$ represent the simulated and measured sample values, respectively. $N$ is the number of valid samples. Correspondingly, $\overline{M}$ and $\overline{O}$ are their averages.

3. Results

3.1. Consistency of PM_2.5, H-8 AOD, and CALIPSO AEC

Figure 2 shows the PM_2.5 mass concentration, H-8 AOD, and CALIPSO AEC in the d02 region at 0500 UTC on 12 January 2019. As the distribution of PM_2.5 shows (Figure 2a), most of the northeastern, central, and southwestern regions of the d02 region had severe aerosol pollution, with the PM_2.5 exceeding 200 μg m^–3 at many sites in Tianjin, Hebei, Shandong, Shanxi, northern Henan, and central Liaoning. However, the northwestern and southeastern regions of the d02, such as Inner Mongolia and southern Jiangsu, had relatively cleaner air. Compared to PM_2.5, there was a strong correlation between the distribution of the AOD (Figure 2b) and PM_2.5. High AOD values (maximum of 1.8 in Hebei, Shandong, and central Liaoning) corresponded to the high PM_2.5, and small AOD values in Inner Mongolia were consistent with relatively light pollution. This consistency may be because pollution in the atmosphere was mainly concentrated in the boundary layer, such that aerosols near the surface contributed the most to the AOD. Although most of the available AOD was consistent with PM_2.5, there were many regions in d02 where there is no AOD measurements, primarily due to the influence of clouds and partly because some AOD data were excluded during the quality control process. The CALIPSO AEC profiles were detected along the line shown in Figure 1b. As shown in Figure 2c, the AEC was relatively large with a maximum value of more than 2 km^–1 between 36°N and 39°N at heights less than 1 km, corresponding to the high PM_2.5 and AOD, while the AEC is relatively small, with most values being less than 0.4 km^–1 between 31°N and 35°N as well as 41°N and 44°N, corresponding to the low AOD and PM_2.5. Overall, there were strong correlations among the three observational datasets of PM_2.5, H-8AOD, and CALIPSO AEC. Therefore, it is reasonable and effective to use these observational data for assimilation and evaluation purposes.

3.2. Comparison to PM_2.5

The scatter plots in Figure 3 show a comparison between the PM_2.5 measurements and the PM_2.5 in the control IC (and DA analysis). The PM_2.5 in the IC of the control experiment (blue points) had a large dispersion relative to the central axis, indicating that it deviated greatly from the observation. In contrast, the distance of PM_2.5 in DA analysis (red points) from the central axis was lower, indicating that DA could considerably correct the deviation of the control IC. In addition, the difference between the distribution characteristics of the scatter points for experiments Da_AOD_PM and Da_PM were relatively small, indicating that the two DA experiments had similar improvement effects on PM_2.5. The CORR, MFE, and RMSE of the control IC were 0.20, 0.62, and 107.71 μg m^–3, while that of the Da_PM were 0.97, 0.16, and 24.94 μg m^–3, and that of the Da_AOD_PM were 0.96, 0.18 and 27.95 μg m^–3, respectively. Both DA groups had a significant improvement effects on PM_2.5, and the Da_PM experiment achieved a slightly better performance than did the Da_AOD_PM experiment.

3.3. Comparison to AOD

The AOD, being the integral of the AEC of the entire atmospheric layer in the vertical direction, can reflect the aerosol extinction characteristics of the entire atmosphere. Figure 4 shows the distribution of AOD in the control IC and DA analysis, as well as the biases of the control IC and DA increments. The AOD was relatively large in the southeastern part of the control IC (Figure 4a), with a maximum value of 2.0, and was relatively small in the northwestern part. This spatial characteristic was consistent with that of the entire measurement (Figure 2b), but the control IC underestimated the AOD in the Liaoning–Tianjin–Hebei region (Figure 4d). In addition, the performance of the control IC could not be directly evaluated for the southeastern part of the d02 region because there were no AOD measurements taken there. The DA_PM substantially corrected the low bias of the control IC in the Liaoning–Tianjin–Hebei region (Figure 4c,f), indicating that the PM_2.5 DA can improve the accuracy of the AOD analyses of the IC. However, most of the AOD in the Da_PM analysis were still smaller than those observed, especially in Hebei, indicating that the improvement effect of PM_2.5 DA on AOD still had room for improvement. Compared with the Da_PM, the AOD in the Da_AOD_PM analysis (Figure 4b,e) was larger in regions such as Liaoning–Tianjin–Hebei and northern Henan, with a maximum value as large as 2.0 km⁻¹, which is more consistent with the measurements. The results show that, after the adoption of the AOD observations, DA analysis can better simulate the distribution characteristics of the AOD. In summary, the assimilation of PM_2.5 observations can improve the distribution of ground pollutants and, hence, has a positive effect on the accuracy of AOD analyses. However, it is difficult to significantly improve AOD analyses only by assimilating PM_2.5 data, while the introduction of AOD assimilation can perform better, thereby improving the description of air pollution in the entire atmospheric layer.

The scatter plots in Figure 5 show a comparison between the Himawari-8-derived AOD measurements and the simulated AOD for the control IC (and DA analysis). The IC for the control experiment (blue points) caused a significant negative bias against AOD observations, especially for AODs larger than 1.0. The Da_AOD_PM (red points in Figure 5a) considerably corrected the negative bias of the control, such that the distance of the AOD in the Da_AOD_PM analysis from the central axis was the smallest among the three experiments. In contrast, the Da_PM (red points in Figure 5b) corrected the negative bias of the control, but the improvement effect, especially for AODs larger than 1.0, was not significant. The distribution characteristic of the scatter points for the Da_PM was similar to that of the control, indicating that the DA_PM had a limited improvement effect on the AOD analyses. The CORR, RMSE, and MFE of the control IC were 0.83, 0.38, and 0.69, respectively, while those of the Da_PM were 0.83, 0.37, and 0.59, and those of the Da_AOD_PM were 0.95, 0.16, and 0.19, respectively. The Da_AOD_PM had the best performance for the improvement of the AOD analyses, especially for the AODs larger than 1.0.

3.4. Comparison to CALIPSO AEC

The AEC profile contains information on the distribution of atmospheric aerosols in the vertical dimension. Figure 6 shows the average AEC profile of the CALIPSO measurements, control IC, and DA analysis along the line shown in Figure 1b. Only the grid cells containing all the CALIPSO, PM_2.5, and AOD measurements were selected for analyses. As a result, 69 CALIPSO profiles were included in the analyses. The CALIPSO AEC profile (black line) showed an overall trend of first increasing from the ground to a height of about 200 m and then decreasing, and there were unsmooth changes below 1 km, reflecting the complexity of the aerosol vertical distribution in the boundary layer. In contrast, the profile of the control IC (blue line) decreased smoothly from the ground to a height of 1 km, which cannot reflect the fine vertical scale details, as measured. This may be because the model mainly describes the overall condition of the atmosphere; therefore, point and sub-grid information in the boundary layer cannot be characterized successfully. There were small differences between the profile of the Da_PM analysis and the control IC below 500 m, but there was no significant difference above 500 m. This is mainly because the assimilation of PM_2.5 introduces observation information near the surface, which can help to improve aerosol analyses near the ground, but the improvement decreased dramatically as the height increased. The Da_AOD_PM (red line) significantly corrected the bias of the control IC and successfully simulated the first increasing and then decreasing trend of the AEC below 1 km. The AEC profile of the Da_AOD_PM was closest to the CALIPSO profile, indicating that assimilating AOD data can substantially improve aerosol analyses in the vertical dimension. However, it should be noted that there is still a large gap between the DA_AOD_PM result and the measurement for an AEC greater than 1 km⁻¹ near the ground, and that the fine vertical scale details in the CALIPSO profiles have not been described. This may be because the adjustment of the vertical structure of aerosols by AOD assimilation is affected by the background field profile as well as the vertical distribution characteristics of the background error covariance. This indicates that AOD assimilation without additional constraints may not perform as well as AEC profile assimilation at improving the vertical structure of aerosols.

3.5. Effects of DA on the Forecast Performance for PM_2.5 and AOD

Figure 7 shows the time series of CORR, MFE, and RMSE for the PM_2.5 model forecasts with and without DA. Owing to the improvement in the accuracy of chemical IC, both the DA_PM (green line) and DA_AOD_PM (red line) experiments had better forecast skills than the control experiment (blue line) during the 60 h forecast period, with a larger CORR, smaller MFE, and smaller RMSE, indicating that the positive effect of DA can persist for more than 60 h. The difference between the results of the DA and the control was the largest at the initial time and decreased in a fluctuating manner with the forecast time, indicating that the positive effect of DA gradually decayed as the integration time increased. The difference between the results for the Da_AOD_PM and Da_PM experiments is small, indicating that assimilating PM_2.5 data alone and simultaneously assimilating PM_2.5 and satellite AOD data have equivalent improvement performances for PM_2.5 forecasting.

The scatter plots in Figure 8 show a comparison between the AOD measurements and the 24 h AOD forecasts (at 0500 UTC on 13 January 2019) of the control and DA experiments. Similar to the performance of the ICs (Figure 5), the 24 h forecasts for the control experiment (blue points) also had a negative bias against the AOD observation; the CORR, RMSE, and MFE for the control forecasts were 0.64, 0.46, and 0.76, respectively. The Da_AOD_PM results (red points in Figure 8a) were more consistent with the AOD observations than they were with the Da_PM and control. In contrast, the Da_PM (red points in Figure 8b) had better forecasting skills, but the improvement was limited. The CORR, RMSE, and MFE of the Da_PM forecasts were 0.68, 0.44, and 0.73, which were 0.04 (6.25%), 0.02 (4.35%), and 0.03 (3.95%) better than those of the control, respectively. The CORR, RMSE, and MFE for the Da_AOD_PM were 0.73, 0.38, and 0.54, which were 0.09 (14.06%), 0.08 (17.39%), and 0.22 (28.95%) better than those of the control, respectively. The Da_AOD_PM had the best performance for improving 24 h AOD forecasts.

4. Discussion

Early studies on AOD assimilation have mostly focused on the improvement effect on ground particulate matter or AOD and seldom studied the improvement in the vertical structure of aerosols. In this study, by comparing the results with and without DA using PM_2.5, AOD, and CALIPSO AEC profile measurements, we found that, compared with PM_2.5 assimilation, AOD assimilation had limited improvement on ground PM_2.5 forecasting, but the improvement effect on aerosol vertical structure was significant. However, it should be noted that, when the background field deviated greatly from the actual profile, especially when there were considerable unsmooth changes in the real profile, there was still a large gap between the DA analysis profile and the real profile (Figure 6). This may be because the adjustment of the vertical structure of the aerosol by AOD assimilation was influenced by the background field profile as well as the background error covariance. In most applied 3DVAR systems, the background error covariance is usually not flow-dependent; therefore, the consistency between the vertical structure of the background error covariance and the real background bias may vary in different meteorological situations, and, hence, the effect of AOD assimilation may be unstable. The introduction of a flow dependent background error covariance or vertical aerosol information that constrains the assimilation increment in the vertical direction is a potentially effective solution to the problem and requires further study.

5. Conclusions

In this study, we evaluated the impact of assimilating Himawari-8 AOD measurements on surface and vertical aerosol analyses over northern China using the MOSAIC aerosol scheme of WRF-Chem and the 3DVAR method. One control experiment and two DA experiments were conducted. Both DA groups showed improvements in the analysis accuracy and forecasting skills of the spatial distribution characteristics of aerosols, although the effects of the two experiments differed from each other. When PM_2.5 was used for verification, the improvements of the two DA groups were similar to each other. The PM_2.5 DA experiment improved the CORR of PM_2.5 in the analysis field by 0.77, and reduced the RMSE and the MFE by 82.77 μg m^–3 and 0.46, while the results of the Da_AOD_PM were 0.76, 79.76 μg m^–3 and 0.44, respectively. When AOD was used for verification, the PM_2.5 DA experiment showed little improvement in the analysis accuracy of AOD distribution compared with the control experiment, while Da_AOD_PM showed significant improvement. The Da_AOD_PM improved the CORR of the AOD in the analysis field by 0.12 and reduced the RMSE and the MFE by 0.22 and 0.50, respectively. Improved performance of the Da_AOD_PM occurred when the AEC profile was used for verification, as when the AOD was used for verification. Da_AOD_PM successfully simulated the first increasing and then decreasing trend of the aerosol extinction coefficients below 1 km, while neither the control experiment nor the PM_2.5 DA experiment did. Both DA groups had positive effects on the PM_2.5 mass concentration forecasting skills for more than 60 h, and the Da_AOD_PM had the best positive effect on the AOD forecasting skills. The CORR, RMSE, and MFE of the 24 h AOD forecasts for the Da_AOD_PM were 0.73, 0.38 and 0.54, which were 0.09 (14.06%), 0.08 (17.39%) and 0.22 (28.95%) better than that of the control, and 0.05 (7.35%), 0.06 (13.64%) and 0.19 (26.03%) better than that of the PM_2.5 DA experiment, respectively. This indicates that assimilating AOD can effectively improve the analyses and forecasting accuracy of the aerosol structure in both horizontal and vertical dimensions and compensate for the limitations associated with assimilating traditional aerosol data alone.