Land use and anthropogenic heat modulate ozone by meteorology: a perspective from the Yangtze River Delta region

With the rapid advance in urbanization, land use and anthropogenic heat (AH) dictated by human activities significantly modify the urban climate and in turn the air quality. Focusing on the Yangtze River Delta (YRD) region, a highly urbanized coastal area with severe ozone (O3) pollution, we estimate the impacts of land use and AH on meteorology and O3 using the Weather Research and Forecasting model coupled to Chemistry (WRF-Chem). These results enhance our understanding of the formation of O3 pollution in rapidly developing city clusters with place-specific topography, as most of our results can be supported by previous studies conducted in other regions around the world. Regional O3 pollution episodes occurred frequently (∼ 26 times per year) in the YRD from 2015 to 2019. These O3 pollution episodes are usually in calm conditions characterized by high temperature (over 20 C), low relative humidity (less than 80 %), light wind (less than 3 m s−1) and shallow cloud cover (less than 5 okta). In this case, O3 pollution belts tend to appear in the converging airflows associated with the sea and the lake breezes. On the other hand, rapid urbanization has significantly changed land use and AH in this region, which subsequently affects meteorology and O3 concentration. The largest change in land use comes from urban expansion, which causes an increase in 2 m temperature (T2) by a maximum of 3 C, an increase in planetary boundary layer height (PBLH) by a maximum of 500 m, a decrease in 10 m wind speed (WS10) by a maximum of 1.5 m s−1 and an increase in surface O3 by a maximum of 20 μg m−3. With regard to the sea and lake breezes, the expansion of coastal cities, like Shanghai, can enhance the sea breeze circulation by ∼ 1 m s−1. During the advance of the sea breeze front inland, the updraft induced by the front causes strong vertical mixing of O3. However, once the sea breeze is fully developed in the afternoon (∼ 17:00 LT), further progression inland will stall. Then O3 removal by the low sea breeze will be weakened, and surface O3 can be 10 μg m−3 higher in the case with cities than in the case with no cities. The expansion of lakeside cities, such as Wuxi and Suzhou, can extend the lifetime of lake breezes from noon to afternoon. Since the offshore flow of the lake breeze transports high O3 from the land to the lake, the onshore flow brings high O3 back to the land. Surface O3 in lakeside cities can increase by as much as 30 μg m−3. Compared to land use, the effects of AH are relatively small. The changes mainly appear in and around cities where AH fluxes are large. There are increases in T2, PBLH, WS10 and surface O3 when AH fluxes are taken into account, with increments of approximately 0.2 C, 75 m, 0.3 m s−1 and 4 μg m−3, respectively. AH contributes largely to the urban environment, altering meteorological factors, O3 concentration and urban breeze circulation, but its effect on the sea and the lake breezes seems to be limited. Published by Copernicus Publications on behalf of the European Geosciences Union. 1352 C. Zhan and M. Xie: Land use and anthropogenic heat modulate ozone by meteorology


Introduction
Tropospheric O 3 is a secondary pollutant formed by a series of complex chemical reactions (Chameides and Walker, 1973;Xie et al., 2014) of precursor gases such as nitrogen oxides (NO x , which is NO + NO 2 ) and volatile organic compounds (VOCs) in combination with sunlight. It has received continuous attention over the last few decades due to its negative effects on the human respiratory system (Jerrett et al., 2009) and the growth of vegetation (Mills et al., 2011). The global average lifetime of tropospheric O 3 is 20 to 25 d, and it will decrease to 5 d in the boundary layer (Young et al., 2013). The relatively long lifetime of O 3 favors regional transport and brings huge challenges to its control (Bergin et al., 2005). O 3 concentrations considerably depend on weather conditions because they play an important role in determining the chemistry, dispersion and removal of O 3 (Jacob and Winner, 2009). Elevated O 3 usually occurs in warm dry weather with strong sunlight, high temperature, low relative humidity and light wind speed (Zhang et al., 2015). Furthermore, weather conditions have many similarities in terms of certain weather patterns (Buchholz et al., 2010;Zhan et al., 2019), and the main weather patterns associated with O 3 pollution episodes in China are tropical cyclones and continental anticyclones .
O 3 concentration and meteorology in urban areas are of great concern because urban areas have huge populations. A report from the United Nations has pointed out that 69.6 % of the world's population will live in cities by 2050 . Because of historical and cultural factors, many cities have similar topography, usually along the coasts, close to mountains or in basins. For these cities, the local circulations induced by the thermal contrast of the topography, such as sea-land breezes, mountain-valley breezes and lake-land breezes, have an important impact on urban air quality under calm conditions . Examples can be found around the world. Ding et al. (2004) simulated the main features of sea-land breezes during a multiday episode in the Pearl River Delta (PRD) region and found that the sea-land breezes can transport air pollutants between inland and coastal cities. Miao et al. (2015) studied the effects of mountain-valley breezes on the boundary layer structure in the Beijing-Tianjin-Hebei (BTH) region, suggesting that mountain-valley breezes are vital to the vertical transport and distribution of air pollutants in Beijing. Wentworth et al. (2015) identified a causal link between lake breezes and O 3 in the Greater Toronto Area, finding that the daytime O 3 maxima were 13.6-14.8 ppb higher on lake breeze days than on no-lake-breeze days.
Ongoing urbanization, including changes in land use and anthropogenic heat (AH), can affect meteorology and atmospheric compositions at local, regional and even global scales (Fu and Liao, 2014;Park et al., 2014;Oke et al., 2017). Land use changes via urban expansion (typically from vegetation to impervious surfaces) directly alter the surface physical properties (e.g., albedo, surface moisture and roughness), subsequently affecting the exchange of energy, moisture and momentum and hence impacting the urban climate and air quality (Jiang et al., 2008;Wang et al., 2009). Li et al. (2019) found that increases in thermal inertia, surface roughness and evapotranspiration due to urban expansion can increase O 3 by 5.6 ppb in southern California. AH is an important waste by-product of urban metabolism. Nearly all energy consumed by human activities will be dissipated as heat within Earth's land-atmosphere system (Flanner, 2009;Sailor, 2011), which is then "injected" into the energy balance processes. Ryu et al. (2013) reported that AH affects the characteristics/structures of the boundary layer and local circulations, resulting in an increase in O 3 by 3.8 ppb in the Seoul metropolitan area.
Previous studies have usually investigated the impacts of topography, land use and AH on meteorology and air quality separately and have mainly focused on a specific megacity. However, these factors can work together in near-calm conditions. Furthermore, complex interactions exist widely among thermally driven circulations, and the effects can even spread from one city to nearby areas. For example, Zhu et al. (2015) demonstrated that the meteorological conditions and air quality over Kunshan are significantly affected by Shanghai urban land-surface forcing (Kunshan is located downstream of Shanghai with a straight-line distance of approximately 50 km). Given the increasing prevalence of cities, cities gradually appear in the form of clusters. Therefore, assessing the effects of land use and AH (the topography rarely changes) in city clusters is meaningful, which also helps us to understand interactions between the urban environment and human activities.
The Yangtze River Delta (YRD) region, located on the western coast of the Pacific Ocean (Fig. 1a), has undergone accelerated urbanization processes and rapid economic development over the past few decades. It is now one of the largest economic zones in the world. The YRD region consists of the southern part of Jiangsu Province, the northern part of Zhejiang Province and the eastern part of Anhui Province, including 26 megacities/large cities such as Shanghai, Hangzhou and Nanjing (Fig. 1b). With dense population and huge energy consumption, this area is now suffering from air quality deterioration (Xie et al., 2017;Zhan et al., 2020), especially increasingly severe O 3 pollution in recent years . Furthermore, cities with hot spots of O 3 are usually concentrated in the central YRD region surrounding Lake Tai (Zhan et al., 2021). Numerous cities, unique topography and severe O 3 pollution make the YRD an ideal study place.
In this study, the impacts of land use and AH on meteorology in the central YRD region and how these effects further modulate O 3 are investigated using the Weather Research and Forecasting model coupled to Chemistry (WRF-Chem). These results fill a knowledge gap about the formation of O 3 pollution in this region and provide valuable insight for other rapidly developing regions with complex topography in the world. The remainder of this paper is organized as follows. Section 2 gives a detailed description of the observation data, the model setup and the experimental design. The main results, including the characteristics of O 3 pollution episodes, the model evaluation, and the changes in meteorology and O 3 caused by land use and AH, are presented in Sect. 3. A summary and conclusions are given in Sect. 4.

Surface observations
Hourly O 3 concentrations monitored by the National Environmental Monitoring Center (NEMC) of China are used in this study. These data strictly follow the national monitoring standards HJ 654-2013and HJ 193-2013(http://www.cnemc. cn/jcgf/dqhj/, last access: 18 January 2022 and are available at https://quotsoft.net/air/ (last access: 18 January 2022). The nationwide observation network initially operated in 74 major cities from 2013, and it had grown to more than 1500 stations covering 454 cities by 2017 (Lu et al., 2018). The urban hourly O 3 concentrations are the average results of measurements at all monitoring sites in each city. The daily maximum 8 h average (MDA8) O 3 concentrations are calculated using the hourly O 3 concentrations for days with more than 18 h of measurements (Liao et al., 2017).
Meteorological data are taken from the National Climatic Data Center (NCDC), including 2 m air temperature (T 2 ), relative humidity (RH), 10 m wind speed (WS 10 ) and direction (WD 10 ), and cloud cover (CC). These data as well as the technical documents recording the quality control, data collection and archive are available at ftp://ftp.ncdc.noaa.gov/ pub/data/noaa/isd-lite/ (last access: 18 January 2022). Locations of the weather stations are shown in Fig. 1b. Specifically, the weather stations in the innermost domain are Pudong (Pd), Shanghai (Sh), Hongqiao (Hq) and Xiaoshan (Xs).

MODIS-based and USGS land use classifications
To explore the effects of land use, two land use categories set as default in WRF are used to set up the first two scenario simulations ( Table 2). The MODIS-based land cover product was created from 500 m MODIS Terra and Aqua satellite imagery (Friedl et al., 2010) and has replaced USGS as the default settings in WRF since version 3.8. The USGS data were primarily derived from the Advanced Very High Resolution Radiometer (AVHRR) from 1992 to 1993 at a 1 km spatial resolution (Loveland et al., 2000), which reflects the distribution of cities in the late 1980s. Figure 2 presents the land cover maps in the innermost domain. The most obvious difference between MODIS and USGS is the urban fraction, which is related to the urban expansion caused by rapid urbanization. In addition, the Finer Resolution Observation and Monitoring of Global Land Cover in 2015 (GLC), one of the latest (2015) and finest (30 m) land cover datasets (Gong et al., 2019), is quite consistent with the performance of MODIS in this region. This confirms that the urban fraction in MODIS is close to reality. Thus, the MODIS data can generally refer to today's distribution of cities.

Anthropogenic heat flux modeling
Another scenario simulation incorporates the urban canopy model with gridded AH fluxes to diagnose the impact of AH. The AH fluxes were calculated based on the statistical data of energy consumption of China in 2016 and then were gridded as 144 rows and 144 columns with a resolution of 2.5 • using population density downloaded from Columbia University's Socioeconomic Data and Applications Center. AH fluxes with their diurnal variations are considered by adding them to the sensible heat flux from the urban canopy layer within the single-layer urban canopy model (SLUCM). The AH fluxes for each grid are determined by the fixed AH value for the urban land use category, the urban fraction value on each grid and the fixed temporal diurnal pattern. Details on the calculation of AH fluxes and how to add AH fluxes into the model can be found in Xie et al. (2016a, b). Figure 3 gives the spatial distribution of AH fluxes in the innermost domain. In urban areas, the AH fluxes usually exceed 20 W m −2 . Some megacities, like Shanghai, can have an AH flux value as high as 200 W m −2 . Except for the urban areas, the AH fluxes are generally less than 5 W m −2 in most parts of the YRD region. In particular, in places where there is no human activity, the AH flux is 0.

Model setup and experimental designs
In this study, WRF-Chem version 3.9.1 was applied. The WRF-Chem model is a fully coupled online numerical weather prediction model with chemistry components , in which chemical and meteorological variables use the same coordinates, transport schemes and physics schemes in space and time. There are 32 vertical levels extending from the surface to 100 hPa with 12 levels located below 2 km to resolve the boundary layer processes. The domain and options for physical and chemical parameterization schemes are summarized in Table 1. The initial and boundary conditions of meteorological fields are from the National Centers for Environmental Prediction (NCEP) global final analysis fields every 6 h at 1 • × 1 • resolution. The anthropogenic emissions are provided by the Multi-resolution Emission Inventory for China (MEIC) in 2017 with a resolution of 0.25 • (http://meicmodel.org/), which includes 10 air pollutants and CO 2 from power, industry, residential, transportation and agriculture sectors. The biogenic emissions are calculated online using the Model of Emissions of Gases and Aerosols from Nature (MEGAN) available in WRF-Chem (Guenther et al., 2006). As our main objective is to explore  the response of O 3 to the meteorological changes induced by land use and AH, we use the same surface biogenic emission rates for different land use scenarios (Li et al., , 2017a. Further studies will be carried out to quantify the contribution of biogenic volatile organic compounds changed by meteorological conditions to O 3 . As shown in Table 2, three numerical experiments are performed. The MODIS_noAH experiment is a control simulation with commonly used settings. Compared with MODIS_noAH, USGS_noAH selects the USGS land use classification at runtime through the geogrid program. Thus, the difference between the modeling results of MODIS_noAH and USGS_noAH can illustrate the changes caused by land use. The impact of AH can be identified by comparing the modeling results of MODIS_AH and MODIS_noAH because AH is only added in MODIS_AH. To exclude the uncertainty conceivably caused by different configurations, all three simulations use the same emission inventory and physical and chemical parameterization  (Grell and Dévényi, 2002) RADM2 (Stockwell et al., 1990) Fast-J photolysis (Fast et al., 2006) MADE/SORGAM (Schell et al., 2001)

Model evaluation
To verify the model performance, the simulation results in the innermost domain, including O 3 concentration, T 2 , RH, WS 10 and WD 10 , are examined against the surface observations described in Sect. 2.1. The statistical metrics, including the mean bias (MB), root mean square error (RMSE) and correlation coefficient (COR), are also calculated. They are defined as follows: where S i and O i are the simulations and observations, respectively. N is the total number of valid data, and S and O represent the average of simulations and observations, respectively. Generally, the model performance is acceptable if the values of MB and RMSE are close to 0 and that of COR is close to 1 (Xie et al., 2016a, b;Zhan et al., 2020).

Results and discussions
3.1 Regional O 3 pollution episodes in the YRD On cloudless sunny days, regional O 3 pollution episodes occur frequently in the YRD (Gao et al., 2020;Zhan et al., 2021), affecting an area of 3.5 × 10 6 km 2 and harming more than 200 million people. Regional O 3 pollution is generally defined as when more than half of the 26 typical cities in the YRD fail to meet the national O 3 standard (in China, the national ambient air quality standard for MDA8 O 3 is 160 µg m −3 ). Based on the surface O 3 observations, we sort out all regional O 3 pollution episodes and the correspond-ing weather patterns from 2015 to 2019 (Table S1). There were 20, 19, 34, 28 and 30 regional O 3 pollution cases in the YRD from 2015 to 2019, respectively. These cases mainly occurred in April to October of each year and were usually related to high-pressure, uniform-pressure fields and typhoon activity. Figure 4 further displays the monthly distribution of meteorological factors during the daytime (from 08:00 to 20:00 local time) when regional O 3 pollution occurs in the YRD. All the variables show significant monthly variations. The highest (lowest) temperature is found in July (April), and the relative humidity is highest in June. As for the cloud cover, the sky is covered with fewer clouds in October than in other months. In addition, southeast winds prevail in the YRD from April to October under the influence of the monsoon climate. The correlation coefficients between temperature, relative humidity, cloud cover and wind speed and MDA8 O 3 are 0.12, −0.34, −0.15 and 0.04, respectively. Therefore, O 3 pollution episodes tend to occur on days characterized by high temperature, low relative humidity, cloudless sky and light wind (the weak correlation between wind speed and MDA8 O 3 is due to the small change in light wind). More specifically, on days when the temperature exceeds 20 • C (Fig. 4b), the relative humidity is less than 80 % (Fig. 4c), the cloud cover is less than 5 okta (Fig. 4d) and the wind speed is less than 3 m s −1 (Fig. 4e) in the YRD. On the other hand, local circulations are clearest when there is an absence of clouds, radiative heating is strongest and wind is weakest. In this case, local circulation will inevitably have an impact on the distribution of O 3 .

Case for O 3 pollution episode
For simplicity but without loss of generality, the longestlasting regional O 3 pollution episode is selected to investigate the impacts of land use and AH on meteorology and O 3 in the YRD. This 10 d regional O 3 pollution episode occurred from 25 May to 3 June 2017 (Table S1). Dominated by high-pressure/uniform-pressure fields (Fig. S1), high O 3 concentrations are accompanied by high air temperature, low relative humidity, light wind and shallow cloud cover during this smog episode. An average of 18 out of the 26 cities experienced O 3 pollution every day with MDA8 O 3 concentrations ranging from 80.0 to 269.0 µg m −3 in the YRD. With regard to the meteorological factors, T 2 ranged from 12.9 to 33.5 • C, with an average of 26.4 • C; RH ranged from 26.6 % to 99.4 %, with an average of 58.6 %; WS 10 ranged from 0.5 to 10.8 m s −1 , with an average of 2.8 m s −1 ; CC ranged from 0 to 8.4 okta, with an average of 4.2 okta ( Table 3). The values of these meteorological factors meet the standards in the previous section and can cover both the whole YRD and the central YRD region (the innermost domain). Therefore, this O 3 pollution episode meets the requirements of not only high O 3 concentrations but also calm weather conditions. The long duration also provides relatively universal results.

Evaluation of model performance
To evaluate the model performance, the simulation results in the innermost domain are validated by comparing them with the observational data. Table 4 presents the statistical metrics in meteorological factors. Figure 5 further illustrates the time series of these meteorological factors and their modeling results. T 2 is reasonably well simulated as the CORs (the mean of all the sites) are 0.87, 0.86 and 0.86 in MODIS_noAH, USGS_noAH and MODIS_AH, respectively. The small negative MBs at all sites suggest that our simulations underestimate T 2 to some extent, although this light underestimation is acceptable because of the small RMSEs (2.3, 3.1 and 2.3 • C). The MBs for T 2 in USGS_noAH, MODIS_noAH and MODIS_AH are −2.4, −1.0 and −0.8 • C. The smaller RMSEs and MBs indicate that the T 2 simulation is improved when new land use and AH are taken into account, which may be related to the improved latent and sensible heat fluxes in models (De Meij and Vinuesa, 2014). The improvement can also be confirmed by Fig. 5a. With respect to RH, the CORs are 0.82, 0.75 and 0.83 in MODIS_noAH, USGS_noAH and MODIS_AH, respectively. Thus, all three simulations can capture the diurnal variation in RH well but have different performances at different sites. In USGS_noAH, RH is overestimated at all sites, especially at the Pudong site with an MB of 11.2 %, while RH is overestimated at the two coastal sites (Pudong and Shanghai) but underestimated at two other sites (Hongqiao and Xiaoshan) in MODIS_noAH and MODIS_AH. Moreover, USGS_noAH has the highest RMSEs of RH (16.3 %), followed by MODIS_AH (12.4 %) and MODIS_noAH (12.1 %). Among the three experiments, the simulation result of RH in USGS_noAH is the worst (Fig. 5b). As for WS 10 , the modeling values are slightly overestimated at all sites in all three simulations. The overestimation of WS 10 may partly be attributed to the unresolved terrain features by the default surface drag parameterization causing an overestimation of wind speed especially at low values (Jimenez and Dudhia, 2012). In particular, WS 10 in USGS_noAH is the most overestimated, followed by MODIS_AH and MODIS_noAH, with MBs of 1.2, 1.0 and 0.8 m s −1 , respectively. In addition, high MBs of WS 10 correspond to high RMSEs (1.9, 1.8 and 1.7 m s −1 ) in our simulations. The MODIS-based land cover shows a larger urban fraction and thereby leads to more friction and higher roughness than the USGS dataset. In this case, the overestimation of WS 10 is somewhat neutralized to fit the observations. In terms of WD 10 , the model captures well the shift in wind direction during the study period (Fig. 5d). In summary, both the statistical metrics in Table 4 and the time series in Fig. 5 indicate that all three numerical experiments can capture the major changes in meteorological factors dur-  ing this O 3 pollution episode. Nevertheless, updating the land use and adding AH can reduce the underestimation of T 2 and the overestimation of RH and WS 10 in models. Benefitting from the development of computers, numerical models grow more complex, and the resolution requirements are higher. This puts forward demands on many basic parameters. An updated land cover dataset, like the GLC, contains useful information and can improve the model results (Chen et al., 2016;Wang et al., 2018). Some factors that are generally ignored in models, like the AH fluxes, should also be estimated and evaluated. Table 5 lists the statistical metrics in O 3 , and Fig. 6 gives the hourly variations in O 3 for observations and simulations in different cities. With high CORs (the CORs are 0.80, 0.81 and 0.80 in MODIS_noAH, USGS_noAH and MODIS_AH, respectively), all three simulations reproduce well the diurnal variation in O 3 , which is that O 3 concentration reaches its maximum in the afternoon and decreases to its minimum in the morning. The magnitudes of O 3 modeling results are reasonable (Fig. 6), but the peak and valley values of O 3 simulations sometimes differ from the observations, especially the peak values, like at Huzhou. Considering the relatively low MBs (6.9, −1.6 and 9.0 µg m −3 ) and RMSEs (49.3, 46.2 and 49.0 µg m −3 ), the modeling results of O 3 are generally reasonable and acceptable.
Above all, the WRF-Chem model using our configuration has a good capability in simulating the meteorological factors and O 3 over the studied region. In addition, it is noteworthy that the object of inter-comparisons between the three numerical experiments is not to determine which setting is the most skillful in reproducing the observations. Rather, it is to diagnose and understand the changes induced by land use and AH and the response of O 3 to these changes.

Overall behaviors of O 3 and local circulations
Based on the results of the control simulation (MODIS_noAH), we first give an overall behavior of O 3 and local circulations during the study period. Then the changes induced by land use and AH are discussed via inter-comparisons between different numerical experiments. Therefore, only difference plots between USGS_noAH/MODIS_AH and MODIS_noAH are shown in this paper, and the corresponding original plots for USGS_noAH/MODIS_AH can be found in the Supplement (Figs. S2-S7).

Spatiotemporal variations in O 3
As shown in Fig. 7, O 3 concentration began to rise at around 08:00 local time (LT is UTC+8 h) and became noticeable after only 3 h ( Fig. 7a and h). During this stage, the nocturnal residual layer vanished due to the development of the convective boundary layer (Fig. 8). The O 3 -rich air mass in the residual layer was mixed with the O 3 -poor air mass on the ground, which enhanced the surface O 3 in the morning (Hu et al., 2018). At around 11:00 LT, a convective boundary layer was established (Fig. 8), and high O 3 produced by photochemical reactions appeared over the central YRD and persisted until 18:00 LT (Fig. 7b and c). After sunset, surface O 3 concentration generally decreased due to nitrogen oxide (NO) titration and reached its minimum in the early morning ( Fig. 7f and g). In general, O 3 has a typical diurnal variation with high concentrations in the daytime and low concentrations at night. This is consistent with the observations in Fig. 6, and this rule of O 3 can be applied to most parts of the world. Therefore, the situation during the daytime (we select 11:00, 14:00, 17:00 and 20:00 LT in this study) should be considered in regard to O 3 pollution.

The sea and lake breezes
As shown in Fig. 7a and b, high O 3 concentrations in the central YRD tended to appear in the converging airflows associated with the sea breeze, the lake breeze and the background southeast wind (Fig. 4e). The sea breeze and the background southeast wind were usually in the same direction, and thereby the sea breeze affected a wide area and lasted a long time. The sea breeze was obvious at around 14:00 LT, matured at around 17:00 LT and continuously transported high O 3 from coastal to inland areas during this period (Fig. 7bd). Compared with the sea breeze, the lake breeze had a much smaller influencing area and a shorter duration. At around 11:00 LT, the lake breeze was established. It reached its maximum intensity at around 14:00 LT and then disappeared sharply due to the predominant southeast wind (Fig. 7b and    c). Both the sea and the lake breezes are the typical local circulations in the YRD and play important roles in the horizontal distributions of O 3 over this region. Since the coastline is generally north-south (Fig. 1b), the cross sections along the blue line AB (east-west) depicted in Fig. 2a are illustrated to show representative examples of the vertical structure of the sea breeze (Fig. 9a-c). By 11:00 LT, the sea breeze below 500 m had already developed. The sea breeze front was found in front of Shanghai (∼ 121.6 • E), with a height of 1.5 km (Fig. 9a). At around 14:00 LT, the speed of the sea breeze increased, exceeding 5 m s −1 . The sea breeze front moved inland for a distance of 20-30 km (∼ 121.4 • E) and was elevated to ∼ 2 km (Fig. 9b). At around 17:00 LT, the sea breeze matured, reducing the O 3 concentration near the surface in Shanghai. However, the O 3 in the mixed layer still maintained a high level, which may result in an O 3 -rich reservoir forming in the nocturnal residual layer (Fig. 9c). The penetration of the sea breeze front and its effect on surface O 3 are also observed in other coastal regions, such as Taiwan (Lin et al., 2007), the Athens basin (Mavrakou et al., 2012) and São Paulo (Freitas et al., 2007).
As for the lake breeze, since the lake is inland, the lake breeze has different directions. Thus, the cross sections along the green line CD (south-north) and black line EF (eastwest) in Fig. 2a are given. The lake breeze was established by 11:00 LT ( Fig. 9d and g) although it was shallow at that time. At around 14:00 LT, the lake breeze strengthened. The extension of the lake breeze circulation zone can reach up to 2 km in the vertical dimension (Fig. 9e). The offshore flow of the lake breeze circulation (∼ 2 m s −1 ) transported high O 3 concentrations from the urban areas to the lake, while the onshore flow blew the O 3 back to urban areas ( Fig. 9e and h). Thus, the net effect of the lake breeze is to accelerate the vertical mixing of O 3 in the boundary layer, resulting in high surface O 3 in the lakeside cities. This has also been reported in other lakeside cities, such as Lake Michigan (Lennartson and Schwartz, 2002), the Great Lakes (Sills et al., 2011) and the Great Salt Lake (Blaylock et al., 2017). By 17:00 LT, the lake breeze had disappeared. Figure 10 presents the spatial differences of the main factors (T 2 , PBLH, WS 10 and O 3 ) between MODIS_noAH and USGS_noAH. Land use changes via urban expansion can enhance surface heating through upward sensible heat fluxes so that T 2 will increase. As shown in Fig. 10a-d, the spatial pattern of the remarkable warming effect for T 2 was consistent with the urban fraction change associated with urbanization ( Fig. 2a and b), which is the positive temperature anomaly that mainly appeared in cities and their surrounding areas. In megacities like Shanghai, T 2 increased by as much as 3 • C. The change in PBLH was similar to that in T 2 because the warming up of T 2 was conducive to vertical movement in the boundary layer, which increased the PBLH (Fig. 10eh). The maximum positive change in PBLH can reach up to 500 m at noon but down to 100 m after sunset. With regard to WS 10 , it decreased in MODIS_noAH (Fig. 10i-l), with a maximum decrease up to 1.5 m s −1 in Hangzhou at around 17:00 LT. This is because the roughness of cities and forest is larger than that of cropland ( Fig. 2a and b). Apart from the abovementioned meteorological factors, urban expansion also increased the surface O 3 concentration (Fig. 10m-p). The largest increment of O 3 appeared in the afternoon, with a value of 20 µg m −3 at around 17:00 LT in Changzhou. In addition to these results, it is noteworthy that there were confusing "false" changes at the junction of land and sea/lake, especially for meteorological factors, such as T 2 and WS 10 . This was caused by the different treatments of the MODISbased and USGS land use classifications at the boundary conditions of land versus water instead of urban expansion.

The changes in vertical direction
Urban expansion alters not only the meteorological factors but also the local circulations. As shown in Fig. 11a-c, the sea breeze below 500 m increased by ∼ 1 m s −1 due to the enhanced temperature contrast between the land and the sea induced by the expansion of Shanghai. During the advance of the sea breeze front inland, the updraft induced by the sea breeze front promoted the vertical mixing of O 3 , elevating the surface O 3 concentration in Shanghai ( Fig. 11a  and b). When the sea breeze matured at around 17:00 LT, its transport effect reduced the surface O 3 concentration of the coastal cities (Fig. 9c). However, this "transport effect" was weakened because the sea breeze near the surface was slowed, affected by the rough urban surface. Finally, surface O 3 of ∼ 10 µg m −3 was left compared to the scenario without cities. In contrast to the onshore flow, the offshore flow transported O 3 to the sea, which may be an important source of O 3 in the nocturnal residual layer. Influenced by the strong background southeast wind, the offshore flow was imperceptible during the daytime (Fig. 9), but it was enhanced by urban expansion (Fig. 11c).
As for the lake breeze, it was also enhanced by ∼ 1 m s −1 because of the larger temperature contrast resulting from the expansion of lakeside cities (Fig. 11e and h). Moreover, the life of the lake breeze was extended to 17:00 LT (Fig. 11f and i). Since the lake breeze circulation was conducive to the vertical mixing of O 3 in the boundary layer (Sect. 3.3.2), the O 3 concentration will increase in the lakeside cities, with a maximum of 30 µg m −3 in Wuxi at 14:00 LT.

The mechanism of land use modulating O 3
Changing land use from USGS to MODIS leads to an increase in T 2 by a maximum of 3 • C, an increase in PBLH by a maximum of 500 m and a decrease in WS 10 by a maximum of 1.5 m s −1 in the YRD, which is comparable to those in the BTH region (Yu et al., 2012), the PRD region  and the National Capital Region of India (Sati and Mohan, 2017). These changes are particularly evident in and around cities as the largest change in land use is in urban surface fractions. Elevated air temperature is conducive to the photochemical production of O 3 , and the well-developed convective boundary layer favors the vertical mixing of O 3 . These changes in meteorological factors eventually increase the surface O 3 concentration by a maximum of 20 µg m −3 in the YRD. Furthermore, local circulations, including the sea and the lake breezes, are also influenced by urban expansion, which further alters O 3 in the vertical direction. For coastal cities, like Shanghai, the larger temperature contrast caused by cities enhances the sea breeze. As the sea breeze front moves inland, it enhances upward movement, which is conducive to the mixing of O 3 in the boundary layer. However, the movement of the sea breeze is slowed due to the rough urban surface after the sea breeze matures. The removal of O 3 via the sea breeze near the surface is then weakened. A similar response of the sea breeze to cities as well as its impact on O 3 has also been reported in the PRD region (You et al., 2019) and São Paulo (Freitas et al., 2007). For the lakeside cities, like Wuxi and Suzhou, the lifetime of lake breezes is extended to the afternoon due to the expansion of cities. The offshore flow of the lake breeze transports high O 3 concentrations in the middle of the boundary layer from the land to the lake, while the onshore flow brings the O 3 back to the land, which accelerates the vertical mixing of O 3 and increases the surface O 3 . Thus, high surface O 3 usually appears when the lake breeze is established. This has also been observed in the Greater Toronto Area (Wentworth et al., 2015) and Lake Michigan (Abdi-Oskouei et al., 2020).

Impacts of anthropogenic heat on meteorology
and O 3

Horizontal changes
Compared with land use, the changes caused by AH are much smaller (Fig. 12). Furthermore, these changes are effective in and around cities as they usually have large AH flux densities (Fig. 3). By adding more surface sensible heat into the atmosphere, the AH fluxes led to an increase in T 2 of 0.2 • C in urban areas, with a typical value of 0.42 • C in Shanghai (Fig. 12a-d). Vertical air movement in the boundary layer was then enhanced by the warming of T 2 , and the PBLH increased as well. According to the simulations, the PBLH increased by ∼ 75 m in the urban areas ( Fig. 12e-h). Contrary to the decrease in WS 10 caused by urban expansion, WS 10 increased by ∼ 0.3 m s −1 in the urban areas when AH fluxes were taken into account (Fig. 12i-l). This is ascribed to the strengthened urban breeze circulations, which are conducive to the transmission of momentum from the upper layer to the surface. With regard to surface O 3 concentration, it increased by ∼ 4 µg m −3 in the simulation with adding AH. What is more, the increases in T 2 , PBLH, WS 10 and O 3 were clearer after sunset as the solar shortwave radiation disappeared.

Vertical changes
The phenomenon that cities are almost always warmer than their surroundings is widely known as the urban heat island (UHI), and this urban-rural difference in temperature can further induce urban heat island circulation (UHIC). There was enhanced UHIC driven by AH in the megacity Shanghai at around 14:00 LT (Fig. 13b). This circulation extended horizontally 20-30 km from the city center to the urban edge and vertically to ∼ 2 km from the ground to the top of the urban boundary layer, resulting in a small increase (4-6 µg m −3 ) in surface O 3 . However, for the lakeside cities, the enhanced UHIC was not perceptible. The O 3 concentration in urban areas was reduced on average, with a maximum of 16 µg m −3 in Wuxi at around 14:00 LT (Fig. 13e). The decrease in O 3 may be related to the increased wind ( Fig. 12i-k), which was beneficial to the diffusion and dilution of O 3 . Furthermore, AH has a limited effect on the sea and the lake breezes as it cannot affect any branch of the two as significantly as urban expansion.
3.5.3 The mechanism of anthropogenic heat modulating O 3 AH and land use play different roles in meteorology and O 3 . AH allows the atmosphere to reserve more energy via the ad- Figure 11. Same as Fig. 9 but for the differences between MODIS_noAH and USGS_noAH (MODIS_noAH − USGS_noAH).
ditional sensible heat fluxes, which increases T 2 by ∼ 0.2 • C. Higher temperature induces stronger upward air movement to the development of the convective boundary layer, raising the PBLH by ∼ 75 m. In the convective boundary layer, the atmosphere is associated with turbulent motions and is unstable. Together with the enhanced urban breeze caused by AH, WS 10 increases by ∼ 0.3 m s −1 . These findings are comparable to the values estimated in other cities around the world, such as Philadelphia in the United States (Fan and Sailor, 2005), Winnipeg in Canada (Ferguson and Woodbury, 2007) and Berlin in Germany (Menberg et al., 2013). These changes in meteorological factors eventually lead to an increase in surface O 3 by ∼ 4 µg m −3 . It is noteworthy that the effects of AH are usually clearer at night in urban areas.
In addition, although AH plays an important role in urban breeze circulations, it may not be powerful enough to affect local circulations such as the sea and the lake breezes.

Summary and conclusions
Land use change via urban expansion and the increase in AH release are important manifestations of urbanization. They can alter the regional meteorology and thereby affect O 3 concentrations in and around cities. In this study, the YRD region, a highly urbanized coastal area with severe O 3 pollution, is selected to investigate this issue. Firstly, the basic characteristics of O 3 pollution in the YRD are provided based on the surface observations. Secondly, a representative case is selected for further study using the WRF-Chem model, and the model performance is evaluated by comparing it with the observations. Finally, the response of O 3 to changes in meteorology caused by land use and AH are discussed via the model inter-comparisons. The main findings are listed below: 1. Regional O 3 pollution occurred frequently in the YRD (∼ 26 times per year) from 2015 to 2019. These O 3 pollution episodes mainly occur in calm conditions char-acterized by high temperature (over 20 • C), low relative humidity (less than 80 %), light wind (less than 3 m s −1 ) and shallow cloud cover (less than 5 okta). In this case, the sea and the lake breezes tend to develop and have an important impact on the distribution of O 3 in the YRD.
2. By updating the land use dataset from USGS to MODIS, we find an increase in T 2 by a maximum of 3 • C, an increase in PBLH by a maximum of 500 m and a decrease in WS 10 by a maximum of 1.5 m s −1 in the YRD. The higher temperature and PBLH elevate the surface O 3 concentration by a maximum of 20 µg m −3 via stronger photochemical reactions and vertical mixing processes. Apart from meteorological factors, the sea and the lake breezes are also modified by urban expansion associated with urbanization. The sea breeze is enhanced due to the expansion of coastal cities. Nevertheless, further progression inland of the sea breeze is stalled in the afternoon on account of the rough urban surface, reducing the transmission of O 3 from the coast to the land. The expansion of lakeside cities extends the lifetime of the lake breeze from noon to afternoon. Since the lake breeze accelerates the vertical mixing of O 3 in the boundary layer, the surface O 3 in lakeside cities can increase as much as 30 µg m −3 .
3. When the AH fluxes are taken into account, T 2 , PBLH, WS 10 and O 3 will increase by about 0.2 • C, 75 m, 0.3 m s −1 and 4 µg m −3 , respectively, in and around cities. These changes are relatively small compared to urban expansion and mainly appear around cities where the AH fluxes are large. In addition, unlike urban expansion, AH may have a limited impact on the sea and the lake breezes. But the urban breeze circulations are found to be sensitive to AH inputs.
Studying the impacts of land use and AH forced by human activities on the urban environment is fundamental in improving urban air quality. Although this study only focuses on the YRD region, most of the results can be supported by previous studies conducted in other regions around the world. As an increasing number of city clusters composed of largeand medium-sized cities are built, this work can provide valuable insight into the formation of O 3 pollution in rapidly developing regions with unique geographical features. Data availability. Air quality monitoring data were acquired from a mirror of data from the official NEMC real-time publishing platform (https://quotsoft.net/air/, National Environmental Monitoring Centre of China, 2022). Meteorological data were issued by the NCDC (ftp://ftp.ncdc.noaa.gov/pub/data/noaa/isd-lite/, National Climatic Data Center, 2016). The FNL meteorological data were taken from the NCEP (https://doi.org/10.5065/D6M043C6, National Centers for Environmental Prediction/National Weather Service/NOAA/U.S. Department of Commerce, 2000). The MEIC data were accessible at http://meicmodel.org/ (Li et al., 2017b;Zheng et al., 2018). These data can be downloaded for free as long as you agree to the official instructions. Model output data used for analysis and plotting, as well as the code used for simulations, can be made available upon request.
Author contributions. CZ and MX had the original ideas, designed the research, collected the data and prepared the original draft. CZ did the numerical simulations and carried out the data analysis. MX acquired financial support for the project leading to this publication.

Competing interests.
The contact author has declared that neither they nor their co-author has any competing interests.
Disclaimer. Publisher's note: Copernicus Publications remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Acknowledgements.
We are grateful to the MEPC for the air quality monitoring data, to the NCDC for the meteorological data, to the NCEP for global final analysis fields and to Tsinghua University for the MEIC inventories. The numerical calculations were performed on the blade cluster system in the High Performance Computing and Massive Data Center (HPC&MDC) of the School of Atmospheric Sciences, Nanjing University. We also thank the anonymous reviewers for the constructive comments and suggestions. Review statement. This paper was edited by Andrea Pozzer and reviewed by two anonymous referees.