Ozone Trends during 1979–2019 over Tibetan Plateau Derived from Satellite Observations

Ozone (O3) is a very important atmospheric component related to many atmospheric processes. The Tibetan Plateau (TP) has played an important role to weather and climate in South Asia and was found to modulate the variation of ozone amount, leading to a low-ozone event over the TP in both summer and winter. To better understand the changing trend of ozone, a novel statistic method, that is, the ensemble empirical mode decomposition method, is used to analyze the trends of total column O3 (TCO) and the ozone deficit after removing the corresponding zonal mean of TCO over the TP. Data of TCO over the TP from satellite observations by Ozone Monitoring Instrument/Aura from 2004 to 2019 and Ozone Mapping and Profiler Suite/S-NPP from 2013 to 2019 are used, and, for comparison, the long-term merged ozone data from NASA are used to analyze the trend of TCO in three latitude zones of 25°N–30°N, 30°N–35°N, and 35°N–40°N near the TP. In consistent with some other studies, a slight recovery of ozone around 1996–1997 is evident, with the peak occurring in 2000, but the rate is near zero in 2005. The mean annual trends over 30°N–35°N are −0.836 ± 0.233 DU/yr during 1979–1996 and 0.021 ± 0.124 DU/yr during 1997–2018. Based on the Ozone Monitoring Instrument data, TCO over the TP [27.5°N–37.5°N, 75.5°E–105.5°E] shows a consistent increase trend from 2005 to 2019, and the ozone deficit starts to increase since 2009, indicating that the ozone increase rate over TP is larger than the zonal average. The ozone deficit over the TP shows an even more rapid increase trend after 2015, with a larger increase during winter than summer. Further analysis of the low-ozone event in winter and early summer (MJJ) from 2004 to 2019 shows that the low-ozone event in May is much deeper than that in winter and also lasts longer than that in February. This study confirms the mechanism of the winter low-ozone event over the TP is more complex than that in summer. This finding of more and faster increase in ozone over the TP relative to the zonal mean suggests the possibility of the decrease in the ultraviolet radiation over the TP, which will benefit human skin health and biosphere over the TP; therefore, it is important to continue to monitor the ozone trend and study its mechanism.


INTRODUCTION
As one of the most important trace gases in the atmosphere, ozone plays a controlling role in the oxidation capacity of the atmosphere. Ozone is also a greenhouse gas owing to its strong absorption, centered at 9.6 micron, in the atmospheric window region, which is particularly important in the upper troposphere where the temperature is low.
However, one of the most important parts that have long been recognized is the absorption to ultraviolet (UV) radiation by the total amount of ozone in the atmosphere, which protects the biosphere from harmful UV radiation reaching the ground. Most of the atmospheric ozone exists in the lower stratosphere, and since 1970s, the ozone loss was found in the stratosphere. Many studies found the main cause for the severe ozone depletion is related to the man-made chemicals, such as halocarbon refrigerants and chlorofluorocarbons (CFCs), and these sufficiently long-lived gases reach the stratosphere and release halogens, which destroys ozone (Solomon, 1999). One major milestone of the intergovernmental action to protect the atmospheric ozone layer is the Montreal Protocol and its amendments initiated in 1986, an international binding agreement on phasing out the use of ODS. The global stratospheric ODS peak occurred in approximately 1996.
The use of the ground-based Dobson spectrophotometers to observe stratospheric ozone started in the mid-1920s (Dobson, 1968;Staehelin et al., 1998). The number of ground stations with Dobson spectrophotometer observations increased significantly after the International Geophysical Year (IGY) 1957/1958 (Dobson, 1968). The first measurement of ozone from space occurred in 1970 with the launch of the backscatter UV (BUV) spectrometer, and the follow-on missions include the solar backscatter ultraviolet (SBUV) and total ozone mapping spectrometer (TOMS) instruments . The European Global Ozone Monitoring Experiment (GOME)-type instruments were launched since 1995 (Burrows et al., 1999;Wagner et al., 2001;Richter et al., 2005). The Ozone Monitoring Instrument (OMI) onboard on the NASA EOS Aura satellite was launched on July 15, 2004, which was followed by the Ozone Mapping and Profiler Suite (OMPS) on the Suomi NPP since November, 2011. Up to present, global continuous ozone observations from space span a time period of more than 40 years. From satellite and ground-based observations, a dramatic total ozone column decline of about −3 to −6% per decade (dependent on latitude) throughout the 1980s until the middle 1990s (Pawson et al., 2014) was found. A study by Harris et al. (2008) found a rapid increase in the annual mean total ozone occurred in the late 1990s in the northern hemisphere. A recent analysis to the ozone trends by Weber et al. (2018), using a few different merged data sets from satellite and ground-based observations for the period from 1979 to 2016, suggested that we were about to emerge into the expected recovery phase. However, Weber et al. (2018) found the trends in most data sets and regions were not significantly different from zero since the stratospheric halogen reached its maximum (∼1996 globally and ∼2000 in polar regions). In the southern hemisphere, extratropics, and northern hemisphere subtropics, most data sets showed small positive trends of slightly below +1% per decade, but in the tropics, some data sets showed significant trends of +0.5 to +0.8% per decade and the others showed near-zero trends.
The Tibetan Plateau (TP) is located in the southwest of China and is the highest and biggest plateau in the world. Many studies suggest that the TP acts as a very strong heat source in summer and has a significant impact on the Asian summer monsoon (ASM). The intense convective activities generated at the TP and the large-scale vertical motion associated with the ASM result in the transport of large amounts of sensible heat, moisture, and chemical pollutants from the near-surface layers to upper layers (Ye and Wu, 1998). Using observations by AIRS on EOS/Aqua CH4 products, Xiong et al. (2009) showed an enhanced CH4 plume over the TP during the ASM with the largest plume occurring in the end of August to early September and then dissipated with the withdrawal of ASM. This enhancement of CH4 is mainly associated with strong convection in summer over the TP and the emission from rice paddies, which is a secondary factor. The study of Guo et al. (2017) showed strong ozone valley index over the TP from 1979 to 2009 in summer. Based on the study on ozone variation over the TP, Zhou et al. (2006) pointed out that the TP was a pathway of mass exchange between the troposphere and stratosphere, which was the major factor for the formation of the lower ozone center in summer over the TP. Zhou et al. (2006) also pointed out that the total column ozone (TCO) decrease trend over the TP was one of the strong centers of TCO decrease trend in the same latitude. Tian et al. (2008) analyzed the cause of the significant low TCO over the TP from spring to summer and attributed it to combined effects of high surface altitude and thermal-dynamical forcing. Zhang et al. (2014) pointed out that the increasing surface temperature was one of the important causes for low TCO during winter and summer over TP. Case study from Liu et al., 2009 revealed close relation of low ozone in winter with the variation in upper troposphere-lower stratosphere (UTLS) dynamics over TP.
Using Nimbus-7 TOMS data from 1978to 1991, Zou (1996 found a decrease trend of the ozone, −0.79 ± 0.82 DU/yr (−2.7 7 ± 2.8 percent/decade), with the monthly trends ranging from −0.17 DU/yr (−0.6 percent/decade) to −1.79 DU/yr (6.0 percent/decade). Using TOMS and OMI data, Chen et al. (2016) verified a decrease trend of the TCO over TP from 1978 to 1993, and found the decrease rate slowed down since 1996, but since 2003 the TCO trend increased gradually. Liu et al. (2001) predicted the recovery of ozone trend over TP after 1995 according to simulation from a 2-D chemistry model. Using ground-based observation from 1994 to 2013 at the Global Atmospheric Watch (GAW) station in the north-eastern Tibetan Plateau region (Mt Waliguan, 36°17′ N, 100°54′ E, 3,816 m a.s.l.), Xu et al. (2016) and Xu et al. (2018) showed significant positive trends of surface ozone over the TP, with the largest increase occurring during spring and autumn. The enhanced stratosphere-to-troposphere exchange associated with the strengthening of the mid-latitude jet stream contributed to the observed high ozone anomalies at WLG during the springs of 1999 and 2012. The period with the largest increase occurred around May 2000 and October 2010.
Frontiers in Earth Science | www.frontiersin.org October 2020 | Volume 8 | Article 579624 A mini-hole in winter of 2003(December 14-17, 2003 was first found over the Qinghai-Tibet Plateau by Bian et al. (2006). Chen et al. (2016) found that the low-ozone event during the winter months (November, December and January) had a decreased trend from 1978 to 2015. Using the ozone data sets from the Copernicus Climate Change Service (C3S) to analyze the TCO trends over the TP, Li et al. (2020) found a consistent decreasing EESC (equivalent effective stratospheric chlorine loading) -based trends of −0.56 ± 0.21 DU/yr in winter and −0.30 ± 0.11 DU/yr in summer over TP from 1979to 1996. However, from 1997 to 2017 the EESC-based trends showed increase trends of 0.21 ± 0.08 DU/yr in winter and 0.11 ± 0.04 DU/yr in summer. The decreasing trend of TCO over the TP from 1979 to 1996 was less than the trend of the zonal mean, and from 1997 to 2017, the increase trend of TCO over the TP was also slower than that of the zonal mean. Compared to the low-ozone event from spring to summer, the mechanism of the winter lowozone event is more complex and can be impacted by many factors, like El Niño-Southern Oscillation (ENSO Modoki), the quasi-biennial Oscillation (QBO), 11-year solar cycle, and volcanic aerosol. From the correlation between TCO and different factors, Li et al. (2020) pointed out that fluctuations of geopotential height at 150 hPa (GH150) played a key role in controlling the DJF mean TCO variability over the TP and highlighted the influence of wintertime GH150 variations on summertime TCO trends.
We can see the trends from different studies do not agree with each other, mainly due to the use of different data sets and different analysis methods. This motivates us to use an adaptive time-frequency data analysis method, that is, the ensemble empirical mode decomposition (EEMD), to analyze the trends. In this study, we first analyze the trend of the zonal mean using a long-term merged ozone data from NASA. More focus of this study is to use the OMI-observed TCO data from 2004 to 2019 to study the TCO trend over TP and the ozonedeficit trend over the TP. We also examine the low-ozone event in winter and compare the trends of ozone deficit in winter and summer seasons. A brief introduction of the three measurements' data sets and the EEMD method is presented in Materials and Methods. Results presents the results with the discussion of the ozone and ozone-deficit trends in different seasons. Finally, Conclusion is the conclusion.

NASA Solar Backscatter Ultraviolet Merged Ozone Data Sets v8.6
Through the NASA program MEaSURES (an acronym for Making Earth System Data Records for Use in Research Environments), a monthly mean zonal and gridded average ozone products were constructed by merging data from the nine SBUV instruments: the Nimbus 4 BUV, the Nimbus 7 SBUV, and SBUV/2 on NOAA 11,14,16,17,18, and 19 (McPeters et al., 2013;Frith et al., 2014). This product is referred to as the SBUV merged ozone data sets (MOD). Data from NOAA-9 SBUV-2 and data taken as the equator crossing time as the satellite approaches the terminator are of lesser quality and are excluded from the MOD composite (DeLand et al., 2012;Kramarova et al., 2013). Those instruments are of similar design, and data from each instrument have been reprocessed using the version 8.6 retrieval algorithm . The version 8.6 designation indicates the version 8 algorithm has been used with new calibration. The consistency of calibrations from instrument to instrument is analyzed over the 41-year time series, and then the ozone is retrieved using the Rodgers (1976) optimal estimation approach. The V8.6 data contain ozone profiles in mixing ratio on pressure levels and in Dobson units on layers, and the total ozone is computed as the sum of the layer data. To maintain consistency over the entire time series, the individual instrument records are analyzed with respect to each other, and absolute calibration adjustments are applied as needed based on comparison of radiance measurements during periods of instrument overlap (DeLand et al., 2012). A detail description of the MOD data can be referred to Frith et al. (2014). For total ozone, no additional adjustments to the individual instrument measurements are applied since the differences between SBUV measurements computed during the overlap periods are typically less than the differences between any given instrument and external data sources McPeters et al., 2013;Frith et al., 2014). As there is no physical rationale to identify one instrument as better than the others, MOD comprises all available data, and data from multiple instruments are averaged during the periods of overlap. There are two kinds of data sets, 5 degree zonal mean of TCO and daily merged overpass data for downloading covering period from January 1970 to December 2018. The daily merged overpass data supply ozone profiles at different in situ stations around the world. Comparison of MOD column ozone with measurements from the Dobson/Brewer network shows good agreement within ±1% over the past 40 years . More detailed information and MOD data can be found at https://acd-ext.gsfc.nasa.gov/Data_services/ merged/.

Ozone Monitoring Instrument Total Column Ozone
Ozone monitoring instrument (OMI) is a nadir-viewing pushbroom UV/visible instrument which measures backscattered radiances in three channels covering the 270-to 500-nm wavelength range (UV-1: 270-310 nm; UN-2: 310-365 nm; visible: 350-500 nm) with spectral resolution of 0.42-0.63 nm (Schoeberl et al., 2006;Liu et al., 2010). It was carried on the NASA EOS Aura satellite and launched on July 15, 2004 into a sun-synchronous polar orbit with a local passing time of ∼13: 45. The OMI has a wide field of view (114°) with a cross-track swath of 2,600 km, which makes a daily global coverage with a spatial resolution of 13 km × 24 km (along × cross-track) under nadir-observing model for UV-2 and visible channels, and 13 km × 48 km for the UV-1 channel. Since the OMI measures the entire Huggins ozone absorption bands, it makes it possible to use differential optical absorption spectroscopy Frontiers in Earth Science | www.frontiersin.org October 2020 | Volume 8 | Article 579624 (DOAS) to retrieve the total ozone column. The OMI total ozone DOAS algorithm consists of three steps (Veefkind et al., 2006). First, the DOAS method is used to fit the differential absorption cross section of ozone to the measured sunnormalized radiance, to obtain the slant column density. Then, the slant column density is translated into vertical column density using the calculated air mass factor. This factor depends on the viewing and solar geometry, the fit window, surface albedo and pressure, the actual ozone profile, clouds, and aerosols. The third step consists of a correction for cloud effects, to account for ozone obscured by clouds. Validation of the OMI TCO product indicates good agreement with airborne measurements within less than 1.5% with a standard deviation of nine DU (Kroon et al., 2008). The level 2 ozone total column daily product from October 2004 to present can be achieved at https://giovanni.gsfc.nasa.gov/ giovanni/. . Both OMPS-NP and OMPS-LP are ozone profile sensors but lack sensitivity to the troposphere due to the spectral coverage and observing geometry, respectively. The retrieval algorithm used for the OMPS-NM ozone total column is similar to the NASA version 8 TOMS algorithm used by the OMI. Compared to ground TCO measurements, there is a 2% bias over most latitudes and viewing conditions for the OMPS-NPP level 2 NM TCO product (Bai et al., 2015). The OMPS TCO product can be achieved at https://disc.gsfc.nasa. gov/datasets?keywords omps&page 1&measurement Ozone.

Ensemble Empirical Mode Decomposition Analysis
Ensemble empirical mode decomposition (EEMD) is used to analyze the trend of TCO over the TP. As to trend analysis, there are many curve fits with a priori determined functional forms. Generally, these fits are subjective, and there is no foundation to support that the selected functional forms can express the underlying mechanisms of various data sets, except for the cases where physical processes are completely known (Wu et al., 2007). EEMD is based on EMD (Huang et al., 1998;Huang and Wu, 2008), an adaptive time-frequency data analysis method. Comparing to other decomposing methods, EMD is empirical, intuitive, direct, and adaptive, without requiring any predetermined basis functions. EMDbased decomposition is designed to seek intrinsic modes of oscillation in any data according to the principle of local scale separation. Each acquired intrinsic mode is defined as an IMF (intrinsic mode function). The IMFs are extracted step by step. First, the most high-frequency local oscillations riding on the lower frequency part of the data are extracted, leaving the residual of the original data, and then the next level of highfrequency local oscillations of the residual are extracted, leaving another residual. Extractions go until no complete oscillation can be identified. EMD-based decomposition of data x(t) can be expressed as below: In Eq. 1, C j is IMF; R n represents the residual component, and it could be a constant, a monotonic function, or a function that contains only a single extremum, from which no more oscillatory IMFs can be extracted. The total number of IMFs of a data set is close to log2N, with N being the number of total data points. In practice, each extraction for an IMF uses a sifting process as follows: i) identify all of the local extrema (maxima and minima) and connect these extrema with a cubic spline as the upper/lower envelope; ii) obtain the first component h by taking the difference between the data and the local mean of the two envelopes; iii) treat h as the data and repeat former two steps until the envelopes are symmetric with respect to zero mean under certain criteria, and the final h is the IMF for one

Seasonal Variation of Total Ozone over the Tibetan Plateau and the Ozone Deficits After Removing the Zonal Mean
One important characteristic of the ozone over the TP is the lowozone event in the summer, which is usually referred as ozone deficit and defined by the ozone over the TP minus the zonal mean in the same latitudes (Zou et al., 1996;Zhou et al., 2006). Bian et al. (2006) reported an low-ozone event occurring in winter 2003 but is not so deep as compared to the summer. In our study period using OMI data, we found an ozone deficit over the TP in the June of 2016, which is about −30 DU, or ∼10% ( Figure 1A). And a similar extremely low-ozone event ( Figure 1B) has also been observed in the winter of 2016 as shown in Figure 1A, using 16-day OMI data from December 4 to December 19. By selecting a region at 27.5-37.5°N, 75.5-105.5°E over the TP (see the red box in Figure 1), we computed the monthly average of TCO from 2004 to 2019 ( Figure 2A) and the corresponding seasonal variation of the ozone deficit over the TP ( Figure 2B). The peak of monthly mean ozone appears in March, and the minimum occurs in October (Figure 2A). However, after removing the latitudinal zonal mean, the minimum of ozone deficit appears in the end of spring to summer ( Figure 2B). The reason for the summer low-ozone event over the TP is linked to the strong summer convection over the TP that brings the tropospheric air with low ozone to lower stratosphere .

Trend Analysis of Total Column Ozone in Different Latitudinal Zones near Tibetan Plateau Using Ensemble Empirical Mode Decomposition Method
MOD is a monthly zonal mean and gridded product constructed by merging individual SBUV/SBUV/2 (total and profile ozone)  Frontiers in Earth Science | www.frontiersin.org October 2020 | Volume 8 | Article 579624 data sets, and it covers a period of ∼48°years since 1972. Considering MOD is a gridded data in 5 degree of zonal mean, we use these data over 30°N-35°N latitudinal zone in the TP and its trend is analyzed using the EEMD method. For comparison, the trends of TCO in the south of TP over 25°N-30°N and north of the TP over 35°N-40°N are also analyzed. Figure 3 shows the EEMD decomposition analysis of monthly TCO over the TP, with the right panels showing the components and the left panel showing the residual components of decomposition. To get a better view, R i is plotted by adding (5 − i) × 205 (i 1, 2, 3, 4, 5); the original monthly mean MOD data (green line) and final residual R 5 (red) are plotted by adding (5 × 205 + 80). Similarly, in the right panel, different components C i are plotted by adding (5 − i) × 29 (i 1, 2, 3, 4, 5). IMFs C i represent the TCO variations with frequency from high to low. The amplitudes of IMFs represent the strength of TCO variations. The EEMD method can give IMF components, but the individual component cannot guarantee a well-defined physical meaning which is true for all decomposing methods, especially for the methods with an a priori basis (Huang et al., 1998). The rule for interpreting the physical significance of IMFs is based on the time scales. C 1 shows a fast oscillation while low amplitude with a time scale less than a month. C 2 and C 3 show strong oscillations with time scales of several months, and these two components can reflect the seasonal behavior of TCO. C 4 shows an oscillation with time scale of years, and C 5 shows little oscillation with even longer time scale. After IMFs are extracted, the final residual R 5 represents the TCO trend from EEMD analysis, and the monthly increase rates are calculated from the derivative of residual R 5 . As shown in Figure 4, before 1996, the TCO all have negative trends, but the decrease rate gets smaller and then starts to increase after 1996. The peak of the positive increase trends occurs in ∼2000 over the latitudinal zone 30°N-35°N, which is slightly behind the faster Frontiers in Earth Science | www.frontiersin.org October 2020 | Volume 8 | Article 579624 7

Zou
TCO Trend Analysis over TP recoveries over latitudinal zones of 25°N-30°N and 35°N-40°N. Such a recovery of ozone continues to 2011 in the latitudinal zone 30°N-35°N over the TP, but since 2012, the TCO shows a slightly decrease again. Since most studies found the pinpoint of the trend change occurred in 1996 (e.g., Li et al., 2020), the MOD data from 1979 to 1996 are kept as one period, and the period from 1997 to 2018 is treated as another one. Meanwhile, we broke the data from 1997 to 2018 into two other subperiods: 1997-2003 and 2004-2018, as the OMI TCO data date from 2004, so they can be compared with each other. In the following analysis, the annual trends will be calculated based on EEMD decomposition of MOD data during four periods (1979-1996, 1997-2018, 1997-2003, and 2004-2018). Using the monthly trend shown in Figure 4 to compute the annual mean ozone trend in the latitudinal zone 30°N-35°N, we got the mean annual trends with confidence interval at 95% confidence level which are −0.836 ± 0.233 DU/   Figures 5A,B). Therefore, we can do similar FIGURE 6 | EEMD decompositions of ozone deficit over the TP (A) and the monthly increase rate derived from R5 (B). For comparison, the increase rate from TCO over the TP is also shown in the lower panel (B).

Zou
TCO Trend Analysis over TP analysis to the trend of TCO (and ozone deficit as well) as we already did using the MOD data, but with a shorter time series. Below we just focus on the trend of ozone deficit over the TP. Figure 6A shows the decompositions of ozone deficit over the TP. Since the IMFs C 2 and C 3 from EEMD-based decomposition reflect the seasonal behavior, we found a combination of C 2 and C 3 decomposed from OMI measurements can clearly present the ozone deficit seasonal variation. Compared to TCO, the combination of C 2 and C 3 from ozone deficit shows larger seasonal variation. The increase rates of TCO derived using residual R 5 are shown in Figure 6B. For comparison, the increase rate from TCO over the TP is also shown in Figure 6B. The trend of the TCO over the TP is similar to the trend of the zonal mean TCO in Figure 4 but has slight difference. TCO trend shows a continual recovery over the TP since 2005, with the peak of the increase appearing in 2010. The ozone deficit has a positive trend since 2009 and increases more rapidly after 2015. This indicates that the ozone increase rate over the TP is larger than the zonal mean since 2009. Using the monthly trend of TCO and ozone deficit ( Figure 6B), we got the mean annual increase rate of TCO over the TP with confidence interval at 95% confidence level which is 0.377 ± 0.148 DU/yr during 2004-2019, while the mean annual increase rate of ozone deficit is 0.263 ± 0.437 DU/yr.
To better understand the difference of ozone trend over the TP, Figure 7 shows the ozone deficit in the early summer and winter, and the corresponding increase rates are shown in Figure 8. The low-ozone event usually occurs in May, but in the winter, it occurs mostly in December and February. It is evident that there is relatively faster ozone recovery over the TP than its zonal mean since 2009, with a faster increase in speed after 2015. The recovery of ozone in the winter season is faster than that in the summer season, and December takes the lead in this recovery.    Bian et al. (2006) reported the finding of low-ozone events in the December of 2003. From Figure 7, we found a few large ozone-deficit events in February, May, and December 2004 and February 2019. Considering the ozone deficits in December 2016 and February 2017 are extremely lower than others, as well as OMI and OMPS both have observations in these two periods, we plot the daily TCO change over the TP in these two periods (two left panels in Figure 9) using measurements from the OMI and OMPS. The daily ozone products from the OMI and OMPS agree very well, and both shows a lasting low-ozone event of about two weeks in the December of 2016. By removing the zonal mean ozone, the daily TCO deficit (two right panels in Figure 9) more clearly show the low-ozone days in these two months. Further analysis of low-ozone events in February ), May (2006), and December (2005 has been done based on the daily TCO deficit. According to the statistics of monthly ozone deficit (shown in Figure 2B), we choose the 10 percentile as reference, which means −34 DU for May, −16 DU for February, and −17 DU for December, to count for the extreme low-ozone days below these reference values. Selection of such thresholds is little arbitrary, and the number of the days representing how long the low-ozone event lasts is sensitive to these thresholds. The detailed statistics of these selected months are listed in Table 1. The low-ozone event in May is deeper than that in February and December. The number of extreme low-ozone days is around 17 in May which is equivalent to that in December but longer than that in February; 2015 is the only year in which the trends of ozone deficit in all selected six months are positive. Compared to the low-ozone days in the summer that are mainly controlled by the transport from the troposphere to stratosphere over the TP (e.g., Zhou et al., 2006;Tian et al., 2008), the mechanism of low-ozone days is the winter is more complex (e.g., Li et al., 2020).

CONCLUSION
To better understand the trend of ozone over the Tibetan Plateau, we used ozone measurement from OMI/Aura from 2004 to 2019 and the EEMD method to analyze the trends of total column O 3 and the ozone deficit after removing the corresponding zonal mean TCO over the TP. Analysis of the zonal mean TCO trends using the merged ozone data from NASA from 1978 to 2019 shows a slight recovery of ozone during 1996-1997; however, after 2013, only the north of the TP shows an increase trend, but the zonal mean between 30°N and 40°N shows a decrease trend. The mean annual increase rates over the latitudinal zone 30-35°N are −0.836 ± 0.233 DU/yr during 1979-1996, 0.291 ± 0.056 DU/ yr during 1997-2003, and −0.104 ± 0.14 DU/yr during 2004-2018. The trend of 0.021 ± 0.124 DU/yr during 1997-2018 confirms a positive but relatively small rate of increase after 1997.
Comparison of the zonal mean ozone from 2004 to 2019 from the OMI product with the merged ozone data from NASA, as well as comparison of TCO over the TP from the OMI and the most recent measurements from OMPS on S-NPP, shows that they agree with each other quite well. Using OMI data, we examined the trend of ozone over the TP and the trend of ozone deficit in different seasons. We found that the TCO over the TP shows a consistent increase trend from 2005 to 2019, and the ozone deficit over the TP has a positive trend since 2009 and an even more rapid increase trend after 2015, with a larger increase in the winter season summer. It indicates that the ozone increase rate over the TP is larger than the zonal mean since 2009. This finding of more and faster increase of ozone over the TP relative to the zonal mean suggests the possibility of the decrease of the UV radiation to the ground surface over the TP, which will benefit to human skin health and biosphere there. Further analysis of the low-ozone events in winter and early summer from 2004 to 2019 shows that the low-ozone event in May is much deeper and lasts longer than that in February. Overall, trend analysis from OMI measurements shows positive increase rates during the past 2 decades over the TP. The mean annual increase rates over TP are 0.377 DU/yr for TCO and 0.263 DU/yr for ozone-deficit events during 2004-2019. There are many causes for this trend, such as a local photochemical reaction and transformation. Simulation from the atmospheric chemical model and in situ observations will be used to learn the mechanism of ozone variation in the future study.

DATA AVAILABILITY STATEMENT
The data sets generated for this study are available on request to the corresponding author.

AUTHOR CONTRIBUTIONS
MZ and XX wrote the paper. MZ performed the data curation and visualization. XX and MZ did formal analysis, the result investigation, and paper editing. ZW supplied code of EEMD. CY was responsible for data preprocessing.