Air stagnation in the YRD are calculated based on meteorological data from 2001 to 2017. The data include precipitation, 500 hPa winds, 10 m winds, 2 m temperature, relative humidity, surface pressure and solar radiation. As for precipitation, the TRMM 3B42 data (https://pmm.nasa.gov/data-access/downloads/trmm) are used. TRMM 3B42 products are an estimate of precipitation rate based on the combined instrument rain calibration algorithm, with a temporal resolution of 3 h and a spatial resolution of 0.25° (Huffman et al., 2001; Mao and Wu, 2012). These data have been widely used in China in the field of meteorology, hydrology and water resources management (Sun et al., 2012; Zhang et al., 2013). For solar radiation, the monthly data from 2000 to 2017 over China provided by Feng and Wang (2021) are adopted. This dataset has the spatial resolution of 0.1°. It was generated by the geographically weighted regression method to merge the sunshine-duration-derived solar radiation data with the satellite-derived cloud fraction data (MODAL2 M CLD) and aerosol optical depth data (CERES SYN AOD). The data for other meteorological factors are obtained from the daily ERA datasets provided by European Center for Medium-Range Weather Forecasts (ECMWF) (Dee et al., 2011). The datasets have the spatial resolution of 0.25°. They are widely used and have a good application in China as well (Bao and Zhang, 2013; Wang et al., 2015).

To illustrate the overview of summertime O₃ in the YRD, the air quality monitoring data in the typical YRD cities from 2013 to 2017 are used. The data are acquired from the national air quality real-time publishing platform (http://106.37.208.233:20035), which provides hourly concentrations of six air pollutants (PM_2.5, PM₁₀, SO₂, NO₂, O₃, and CO) over China. These data are strictly in accordance with the national monitoring regulations. The hourly values for each city are calculated by averaging the concentrations at all national monitoring sites in that city. Identification and handling of invalid and lacking data are also manually performed during data processing, following the methods adopted in some previous studies (Xie et al., 2016b; Shu et al., 2016; Gao et al., 2020, 2021; Zhan et al., 2020, 2021). The YRD region consists of 26 typical cities. As shown in Figure 1B, the cities include Shanghai (SH); Nanjing (NJ), Wuxi (WX), Changzhou (CZ), Suzhou (SZ), Nantong (NT), Yancheng (YC), Yangzhou (YZ), Zhenjiang (ZJ), Taizhoushi (TZS) in Jiangsu province; Hangzhou (HZ), Ningbo (NB), Jiaxing (JX), Huzhou (HZ1), Shaoxing (SX), Jinhua (JH), Zhoushan (ZS), Taizhou (TZ) in Zhejiang province; and Hefei (HF), Wuhu (WH), Maanshan (MAS), Tongling (TL), Anqing (AQ), Chuzhou (CZ1), Chizhou (CZ2) and Xuancheng (XC) in Anhui province.

Air stagnation is usually identified by using predefined thresholds of daily upper level winds, near-surface winds and precipitation (Wang and Angell, 1999; Horton et al., 2012; Horton et al., 2014; Huang et al., 2017; Garrido-Perez et al., 2018, 2019). In this study, based on the previous researches, the used meteorological variables include wind speed at 500 hPa (as upper level wind), 10 m wind speed (as near-surface wind) and precipitation. Thus, a given day can be considered as an air stagnant day when the daily mean 500 hPa wind speed is weaker than 13 m s⁻¹, the daily mean 10 m wind speed is weaker than 3.2 m s⁻¹, and the daily precipitation is less than 1 mm (Wang and Angell, 1999; Horton et al., 2014; Huang et al., 2017). We use ASDs to represent the total number of air stagnation days throughout the summer (June–August).

To study the relation between O₃ pollution and air stagnation (Ozone Pollution and its Relation to Air Stagnation Section), S _n is further defined to indicate the air stagnant intensity in a day, given by: S n = ∑ i = 1 3 A i − B i A i (1) where A _i denotes the threshold of the ith meteorological variable and B _i is the exact daily value of the ith variable. In this study, the value of A ₁ for wind speed at 500 hPa is set as 13 m s⁻¹, A ₂ for 10 m wind speed is 3.2 m s⁻¹, and A ₃ for precipitation is 1 mm. The used meteorological data are described in Meteorological Data Section. Large S _n represents strong air stagnant intensity. Because the spatial resolutions of the data are 0.25°, the values of ASDs and S _n are calculated at each 0.25° grid over the YRD (116–123°E, 27–35°N).

The empirical orthogonal function (EOF) analysis has been widely used to split the temporal variance of spatially distributed data into orthogonal spatial patterns called EOFs (Gianelli et al., 2007; Hannachi et al., 2007; Fu et al., 2015). Each EOF has a corresponding eigenvalue, which determines the relative variance contribution to the total variance in the field. Furthermore, the eigenvector, which describes the spatial pattern of the EOFs, is associated with a time series that represents the temporal evolution of that spatial pattern. Most of the variance contribution can be contained into the first few EOFs, which can convey enough information to understand the underlying process. In this study, the NCAR Command Language (NCL) EOF coding package is used. For detailed calculation formulas of EOF as well as the manuals can refer to official website of NCL (https://www.ncl.ucar.edu/Applications/eof.shtml).

To investigate the reasons for the variations of ASDs, the EOFs of ASDs are calculated. Since the time series represent the temporal evolution of its associated spatial pattern, we can effectively obtain the dominant weather systems that affects the spatial pattern by establishing the connection between time series and existing atmospheric circulation indexes (Hannachi et al., 2007; Fu et al., 2015). The atmospheric circulation indexes used in this paper come from 130 climate system monitoring indexes (http://cmdp.ncc-cma.net/Monitoring/cn_index_130.php) issued by the National Climate Center of China. Based on the meteorological reanalysis data, these historical indexes were calculated and provided monthly. Furthermore, we calculate the linear regression coefficients between the time series and the meteorological variables (500 hPa wind speed, 10 m wind speed and precipitation) defined ASD to investigate the specific mechanism. The regression coefficients are equal to the increment of the variables when the time series change by unit 1.

Figure 2 gives the spatiotemporal distribution of ASDs in summer over the YRD from 2001 to 2017. The regional average value of summertime ASDs is about 30 days (32.3% of the entire summer days), similarly as found by Huang et al. (2017). There are significant spatial and temporal variations for ASDs. For the inter-annual variation, the regional average value of ASDs ranges from 9 to 54 days (9.2–58.4% of the entire summer days), with the maximum ASDs value in 2010 and the minimum value in 2015. This may be associated with the anomaly of the intensity and position of the western Pacific subtropical high system (discussed in Meteorological Dynamic Mechanism for Air Stagnation Days in Summer of the YRD Section in detail). For the spatial distribution, the values of ASDs on the ocean are usually smaller than those on land for a particular year. In the YRD, ASDs generally show a considerable regional heterogeneity, which is similar with the finding reported by Garrido-Perez et al. (2018) for Europe and Huang et al. (2017) for China.

EOF analysis of ASDs is carried out to investigate the above spatiotemporal variations. The first three EOFs explain 68.6, 11.3, and 7.1% of the total variations in ASDs, respectively. Figure 3 shows the results for the first EOF (EOF1). As shown in Figure 3A, the spatial pattern of EOF1 is characterized by the same phase over the entire YRD. Figure 3B further illustrates the atmospheric circulation condition causing this consistency. The time series of EOF1 is correlated with the South China Sea subtropical high area index and the South China Sea subtropical high intensity index, with the correlation coefficients of 0.40 and 0.42, respectively (statistically significant at the 90% confidence level), implying that the consistency over the YRD may depend on the area and the intensity of the South China Sea subtropical high. When the South China Sea subtropical high is large and strong, this spatial pattern is usually typical, which means it is likely to appear air stagnation in the YRD. As the main component of the East Asian summer monsoon, the subtropical high acts on the monsoon anomaly through its location, area and intensity, which can directly affect the summertime precipitation, rain belt distribution, and drought/flood anomalies in eastern China (Ding and Chan, 2005). Figures 3C–E illustrate the regression coefficients between the time series of EOF1 and the three meteorological variables used to define ASD. The role of each meteorological variable is different in EOF1. Wind speed at 500 hPa (Figure 3C) and precipitation (Figure 3E) are the main affecting factors. Generally, with the north jump of the subtropical high, its intensity increases and it controls wide areas of eastern China. In this case, the YRD region is dominated by the subtropical high. The weather is usually stable, and the upper and lower wind fields are generally weak (Figures 3C,D). In addition, once the subtropical high is northward, the rain belt generated by the cold and the warm air also moves to north, which can result in less precipitation in midsummer of the YRD as well (Figure 3E). Both weak wind and little precipitation in EOF1 are favorable to cause more air stagnation days.

Figure 4 illustrates the results for the second EOF (EOF2). As shown in Figure 4A, the spatial pattern of EOF2 shows significant maritime-continental contrasts, suggesting opposite changes between ocean and land, which can be explained by the fact that the surface wind is usually stronger on sea than on land (Figure 4D). Figure 4B presents that the time series of EOF2 cannot be linked to existing atmospheric circulation indexes. However, as shown in Figure 4C, it can be found that the north of the YRD is affected by the subtropical jet while the south is affected by the easterlies, corresponding to high wind speeds in the north and low in the south of the YRD at the 500 hPa layer. This may be the factor leading to the spatial pattern of EOF2 on land.

Figure 5 presents the results for the third EOF (EOF3). With regard to the spatial pattern, the ASDs in the north and the south of the YRD has the opposite phase, with a dividing line of approximately 31°N (Figure 5A). Figure 5B shows that the time series of EOF3 has a good relation to the northern hemisphere polar vortex (the tropospheric polar vortex) area index and the northern hemisphere polar vortex intensity index, with the correlation coefficients of −0.42 and −0.42, respectively (statistically significant at the 90% confidence level). When the northern hemisphere polar vortex is small and weak, the north part of the YRD has high ASDs values while the south has low ones. The specific physical process is probably as follows: with the small and weak tropospheric polar vortex, an abnormal high-pressure center is prone to occur above the mid-high latitudes of Asia, which is conductive to the weakening of the westerly wind at 500 hPa (Figure 5C). The atmospheric circulation shows a meridional distribution and the position of the upper-level jet stream is southward. Thus, the cold air can move southward easily, which can form precipitation in the south of the YRD with the warm and wet flow conveyed by the subtropical high (Figure 5E). For another, the subtropical high is usually strong while the tropospheric polar vortex is weak. The near-surface wind speed is low in the areas dominated by the subtropical high (Figure 5B). The tropospheric polar vortex controls the semi-permanent atmospheric center of action, which has important effects on atmospheric circulation at high latitudes. These effects are usually continuous more than seasonal scales. Furthermore, these effects are commonly expressed by significant atmospheric oscillations in climatology, such as the Northern Atlantic Oscillation, the Arctic Oscillation and the East Atlantic-West Russia Pattern (Waugh et al., 2017), which provides a way to link the variances of air stagnation to climate anomalies.

Previous studies, as well as the result of Spatiotemporal Distribution of O ₃ Section, show a positive trend of O₃ concentrations in the YRD during 2013–2017 (Lu et al., 2020; Wang et al., 2020; Zhan et al., 2021). However, the O₃ pollution days account for 17.0, 16.1, 14.7, 15.3, and 20.2% (from 2013 to 2017 respectively) of the summer days of the YRD. It seems that O₃ pollution episodes are more likely to occur in the summer of 2013 than 2015. Additionally, Figure 2 shows that there are more air stagnant days in 2013 than 2015. Figures 3–5 also illustrate that the time series of all three EOFs show opposite phase in 2013 and 2015. Therefore, comparing the meteorological fields between 2013 and 2015 may reveal how air stagnation affects O₃ pollution in the YRD from the perspective of meteorological conditions.

Figure 8 shows the summertime averages and the differences between 2013 and 2015 for the three meteorological variables identifying air stagnation day. In the 500 hPa weather map of 2013 (Figure 8A), the tropospheric polar vortex is large and strong. The flat isopleths indicate that there is a small geopotential height gradient, which can result in weak wind in the YRD. However, in 2015, the tropospheric polar vortex is small and weak. There are obvious troughs and ridges in the upper reaches of the YRD (Figure 8B), and thereby the cold air can move southward. These differences are visually shown in Figure 8C. Compared to 2015, there is a significant positive potential height perturbation in 2013, which can lead to easterly wind component over the YRD and result in a smaller wind speed at 500 hPa atmospheric layer. Furthermore, the westerly component appearing at 45°N confirms the fact that the westerly jet lies more to the north. For 10 m wind speed, it is higher on the ocean than on the land in both 2013 (Figure 8D) and 2015 (Figure 8E), which partly explains the reason why the ocean usually has a smaller ASDs. Figure 8F shows that the 10 m wind speed has a southwest component in 2013 (compared with 2015), which can block the prevailing southeasterly wind.

As for precipitation, there is a significant difference in the spatial distribution between 2013 and 2015. In the summer of 2013, there are more rainfalls in the south of China, and there is less precipitation in the YRD and the north part of China (Figure 8G). In the summer of 2015, however, more rains fall in the YRD, especially in the south of Anhui province (Figure 8H). The difference of precipitation between 2013 and 2015 (precipitation in the summer of 2013 minus that in 2015) are shown in Figure 8I. Obviously, the precipitation over the YRD is much higher in 2015, implying that the processes of the dilution and the removal of air pollutants are more intensive in 2015. Less precipitation in 2013 also means there are more intensive solar radiation reaching the ground and higher air temperature near the surface over the YRD, as shown in Figure 9.

In summary, larger and stronger of the tropospheric polar vortex and the subtropical high result in more air stagnation days over the YRD in the summer of 2013. Under this circumstance, the upper- and lower-level winds are weak, and more air pollutants including O₃ and its precursors are trapped near surface over the YRD. Moreover, less precipitation, higher air temperature and stronger solar radiation that are associated with air stagnation facilitate the photochemical reactions of O₃ formation. Then, more O₃ pollution days occur in the summer of 2013. The weather systems and the changes in relevant meteorological factors related to air stagnation can affect the physical and chemical processes of O₃ formation, and thereby severe O₃ pollution tends to form under air stagnant conditions.