**Xin Gao , Xuan Wang, Hongbin Shi**

School of Urban-rural Planning and Landscape Architecture, Xuchang University, Henan 46100, China

Received:
October 31, 2019

Revised:
May 1, 2019

Accepted:
May 10, 2019

Download Citation:
||https://doi.org/10.4209/aaqr.2018.10.0364

Cite this article:

Gao, X., Wang, X. and Shi, H. (2019). Multifractal Cascade Analysis on the Nature of Air Pollutants Concentration Time Series over China. *Aerosol Air Qual. Res.* 19: 2100-2114. https://doi.org/10.4209/aaqr.2018.10.0364

**Highlights**

- Statistical and spectral analysis of the air pollution data are performed.
- The underlying physical meanings of multifractal parameters are explored.
- The random multiplicative cascade models are applied in the air pollution fields.

**ABSTRACT**

This study investigates the temporal multifractal cascade behaviors of air pollution in major cities of China by using 1-year data of the hourly time series of six pollutants (PM_{2.5}, PM_{10}, CO, NO_{2}, O_{3}, and SO_{2}) and the air quality index (AQI) as model inputs. Given that the air pollution time series generally exhibit positive skewness with a heavy tail, a stable distribution with four parameters, such as Levy *α*-stable distribution, reasonably fits the frequency distributions of the data over the entire range. The periodicity and nonstationarity shared by the series are also spectrally analyzed, and all of the pollutants display a daily and semi-daily cycle as well as two dynamical regimes with a cut-off scale around 11 days. Using universal multifractal analysis, the two parameters (*α* and *C*_{1}) in an air pollution system are measured, and their underlying significance is also addressed via analogues of the Monte Carlo simulation data. Moreover, self-organizing criticality analysis on the data of the air pollution index (API) indicates first-order multifractal phase transitions, with an estimated interval of the order of the moment (*q _{D}*) that varies between 1.95 and 2.98. Finally, downscaled performances of multifractal cascade models are numerically evaluated, demonstrating that log-normal and log-Poisson models, but not

*α*,

*β*, binomial, and uniform models, can effectively recover extreme values.

Keywords:
Cascade; Multifractal; Air pollutant; Time series; Downscaling.

**INTRODUCTION**

Air pollution has become one of the most concerning environmental problems in the world because of its significant negative impacts on human health and restrictions on national development. However, environmental pollution is a complex system, which is essentially the result of a series of physical and chemical interactions between natural and anthropogenic environmental conditions (Mayer, 1999; Kinney, 2008). With the heavy dependence on the fluctuations of a variety of meteorological and emission variables, air pollutant concentrations randomly and irregularly varied with regard to time and space. Knowledge of the structure of the time series of the concentrations of major air pollutants is needed in order to make pollution prevention policies and better understand dynamical mechanisms governing their temporal variabilities (Xepapadeas, 1992; Liu *et al.*, 2003; Esposito *et al.*, 2016).

Time series of air pollution data measured at any location usually distributes as a strong right-skewed unimodal shape (Georgopoulos and Seinfeld, 1982). Identifying optimal parent probability functions for these distributions is fundamental for precisely predicting concentration related statistics (Sedek *et al.*, 2006; Yahaya *et al.*, 2011; El-Shanshoury *et al.*, 2017). Several probability density functions (log-normal, gamma, Weibull, beta, log-logistic, Pearson type V, etc.) have recently been frequently and successfully used to fit the measured data of air pollution (Karaca *et al.*, 2005; Gavriil *et al.*, 2006; Leiva *et al.*, 2010). Some of these studies estimated the exceedance frequencies of air pollutant concentrations over the air quality standard (AQS) and also predicted the needed emission reduction in order to meet the AQS (Marchant *et al.*, 2013). Among all the applied parent probability models, the log-normal distribution is used the most; however, results from this distribution, while fitting well with the observed data with medium level of pollution, do not fit well to heavy pollution cases with high concentrations. In the latter case, the type I two-parameter exponential distribution is a better choice (Maciejewska *et al.*, 2015). Earlier studies have also suggested that the optimal model depends on pollutant species of interest, e.g., the lognormal distribution for particulate matter and the majority of the nitric oxide, oxides of nitrogen, and sulphur dioxide; the gamma distribution for carbon monoxide, nitrogen dioxide, and ozone; and the Weibull distribution for carbon monoxide and ozone (Taylor *et al.*, 1986). In addition, power law distributions have also been frequently employed to fit concentration values of air pollutants (Kütchenhoff and Thamerus, 1996).

Compared to the complicated mathematical models with huge computational cost, the interpretation and regulation about the distributions of concentration values for air pollution time series through some easy-to-operate physical models are obviously inspiring and encouraging. In order to give a physical explanation of the log-normality of pollutant concentrations, Ott (1990) proposed a successive random dilution (SRD) theory which refers to processes of dilution and mixture at a series of stages by constantly adding fresh air into newly added containers, and the resultant concentrations are approximately lognormally distributed. As a counterpart, this phenomenological idea can be traced back to Richardson’s celebrated poem on self-similar cascades describing atmospheric dynamics as a cascade process (Richardson, 1965; Schertzer and Lovejoy, 2011). The idea was strictly conceptualized with help of fractal theory until 1980s, and then it evolved into multifractals. Since then, several probabilistic cascade models have been built and widely applied, such as in the scaling analysis of fully developed turbulence, downscaling of soil moisture, and evaluation of mineral resources (Schertzer *et al.*, 1997; Agterberg, 2007; Mascaro *et al.*, 2010), with the simplest monofractal and multifractal models being *β*-model and *α*-model (Schertzer and Lovejoy, 2011). Moreover, inspired by simple cascade models, the general forms of multifractal, also called multifractal formalisms, were deduced when scales tend to be infinitesimal. Of the two most famous formalisms, the multifractal spectrum *f*(*α*) are adopted for *α*-*f*(*α*) model, while *γ*-*c*(*γ*) model lies in the codimension function *c*(*γ*) respectively (Evertsz and Mandelbrot, 1992; Schertzer and Lovejoy, 1987). In the past 30 years, techniques of specific calculation and implementation for multifractal formalisms have been proposed, which include box-counting method, histogram method, probability distribution, wavelet analysis and detrended fluctuation analysis (Lovejoy *et al.*, 1987; Meneveau and Sreenivasan, 1989; Lazarev *et al.*, 1994; Muzy *et al.*, 1994; Kantelhardt *et al.*, 2002). Indeed, multifractal analysis has gradually become a standard technique for quantitatively delineating the nonlinear evolution of a complex system and the multiscaling characteristics of physical quantity and is used to aid the understanding of the intrinsic regularity and the mechanism of physical changes (Shen *et al.*, 2015). For instance, Shi *et al.* (2015) studied the dynamic characteristics of PM_{2.5} series during a heavy haze episode of Chengdu by constructing a self-organized critical model and showed that the multifractal features of the series were highly correlated with long-range correlation and non-Gaussian distribution. Gao *et al.* (2016) studied the nonlinear characteristics of the spatial distribution of particulate matter in China, resulting in a power law distribution with a spatial scale index of *β* = 0.45 and a mesoscale range of 90–550 km.

To date, most studies focused on characterizing multifractal properties that reflect tempo-spatial structure complexities inherent in the air pollution systems and investigating relationships between air pollutants and meteorological factors for seeking the direct evidences of mechanical variations in fractality. Thus, more research is needed on the multifractal formalisms based on the probabilistic cascade models representing physical processes in air pollution systems. For example, scale invariance mechanism provides a creative idea for straightly establishing a relationship to connect the concentration values at different scales, making it possible to easily obtain the concentration values at infinite finer scales. Such a tool is useful for calculating frequencies at finer scales exceeding air quality standards so as to estimate more accurate reduction in favor of providing better and more effective decision-making for local government. Besides, in-depth discussions about the physical meanings of multifractal parameters in air pollution systems are still lacking, thus limiting their applications in pollution assessment studies. In addition, the continuity of a multifractal spectrum also needs to be addressed, which may provide some useful information for extreme analysis. This is especially the case in the determination of critical values of the moment order *q* for representing typical concentrations when Monte Carlo simulation is used.

The present study tries to fill some knowledge gaps discussed above. Firstly, the basic statistical analysis of the confirmation of right-skewed, fat-tailed, long-range dependent and nonstationary features for air pollution time series is made with spectral analysis and frequency distribution fitting methods. Secondly, the feasibility of using multifractal parameters as an evaluation tool of air pollution is discussed, i.e., using the calculated multifractal parameters to evaluate the pollution levels and identify the main pollutants in seven cities; in addition, the continuities of codimension functions of air pollution time series are discussed to determine the estimation intervals for the moment order *q* and the types of discontinuity points. Finally, the downscaling capability in air pollution is discussed based on several probabilistic cascade modes and multifractal characteristics.

METHODOLOGY

METHODOLOGY

Data Collection and Quality Assurance and Control

Data Collection and Quality Assurance and Control

Data Collection and Quality Assurance and Control

Six air pollutants, including CO (mg m^{–3}), NO_{2} (µg m^{–3}), O_{3} (µg m^{–3}), PM_{10} (µg m^{–3}), PM_{2.5} (µg m^{–3}), and SO_{2} (µg m^{–3}), were measured using automated monitoring systems widely distributed in cities at county level and above. SO_{2}, NO_{2}, and O_{3} were measured by the ultraviolet ﬂuorescence method, the chemiluminescence method, and the UV spectrophotometry method, respectively; CO was measured by the nondispersive infrared absorption method and the gas ﬁlter correlation infrared absorption method; and PM_{2.5} and PM_{10} were measured by the micro-oscillating balance method and the *β* absorption method.

AQI is a dimensionless index to quantitatively describe air quality. An air quality index is defined as the maximum of the individual air pollutant quality indices standardized to the corresponding concentrations of the six pollutants separately based on AQI Equation and air quality standards. Seven cities in China, including Beijing, Shenyang, Shanghai, Wuhan, Chengdu, Guangzhou, and Xi’an, representing different regions as well as climatic and topographic conditions are selected for investigation (Fig. 1). Time series of the hourly concentration of six pollutants and AQI from January 2016 to December 2016 were used in the analysis (available from the China National Environmental Monitoring Center, http://106.37.208.233:20035/). In addition, 13 years (2000–2012) of daily API data are also used for phase transitions and self-organized criticality analysis. Linear piecewise interpolation method is used for filling the data gaps caused by power outages, equipment failures, and other situations.

**Fig. 1.** The distribution of the cities in China that are used for the study.

Methods

Methods

Methods

For a scale invariant measure, *ε*, which can be a value of rainfall intensity or pollutant concentration, a fractional Brownian motion (fBm) can be constructed by incremental variances of the measures at different scales. Moreover, the parameter *H _{u}*, also called the Hurst index, is commonly used to express the intermittency and variation for a fractional Brownian motion. Obviously, a fractional Brownian motion is an addition process, in which each point is generated by the addition of the variances of two adjacent points and a random number, and

*H*can be calculated using the probability density function method and R/S analysis. Corresponding to the addition process is the multiplication process, also called cascade process. For example, if we consider the 1D cascade process, it is assumed that there is an amount,

_{u}*ε*

_{1}, with a certain support of unit length. Firstly, the unit is divided into two subunits of equal length, which are each multiplied by a random variable,

*W*

_{1}, that satisfies a certain probability distribution; these values are denoted as

*ε*

_{2,1}and

*ε*

_{2,2}with <

*W*

_{1}> = 1, where

*W*is the weight and <> means expectation. Secondly, according to the same operation, the resulting subunits of the segmentation in the previous step are further divided, multiplied by the same distributed random variable,

*W*

_{2}, and this process is repeated

*n*times to obtain a series,

*ε*. For

_{n}*ε*, the geometrical supports become

_{n}*r*= 1/2

*n*, meaning the finest scale is 1/2

*n*, and 1/

*r*is

*λ*, where

*λ*is the scale ratio. From the generation process shown in Fig. 2, it can be seen that as the scale ratio increases,

*ε*gradually becomes more intermittent and has a greater oscillation.

_{n}**Fig. 2.** The cascade generation process.

A parameter is needed to describe the fBm scale invariance and intermittence, while the latter requires multiple fractal dimensions to be described. Assuming there is a real value *γ*, as *γ* increases, there is a relationship between the probability distributions of the values larger than *λ ^{γ}* and the scale rations: Pr (

*ε*>

_{n}*λ*) =

^{γ}*λ*

^{–c(γ)}, where

*c*(

*γ*) is a function of

*γ*and is the characterization of the support sets corresponding to each value

*γ*for

*ε*. The function

_{n}*c*(

*γ*) is also called the codimension function and corresponds to the dimension function

*f*(

*α*) proposed by Evertsz and Mandelbrot (1992), where

*γ*and

*α*are singularity indexes. The

*γ*-

*c*(

*γ*) model, a multifractal formalism widely used in the fields of turbulence, rainfall, weather, stocks, and pollution, was first proposed by Schertzer and Lovejoy (1987). To obtain

*c*(

*γ*), we need to calculate the statistical properties of the higher moments <

*ε*> of

_{n}^{q}*ε*, and there is a relationship between <

_{n}*ε*> and the scales: <

_{n}^{q}*ε*> =

_{n}^{q}*λ*

^{K}^{(q)}, where

*K*(

*q*) is the scaling function and there are the following relationships:

*K*(

*q*) = min[

*qγ*– c(

*γ*)] and

*c*(

*γ*) = min[

*qγ*–

*K*(

*q*)], with

*c*(

*γ*) determined by the Lengdre transformation. The above

*γ*

*-c*(

*γ*) model is based on probability theory, while the

*α*-

*f*(

*α*) model is based on measurement theory, and there is also the following relation between the two models:

*f*(

*α*) = 1 –

*c*(

*γ*) and

*α*= 1 –

*γ*. More comparisons between the two models can be found in the research by Cheng and Agterberg (1996). Compared to one parameter of monofractals,

*K*(

*q*) and

*c*(

*γ*) contain an infinite number of parameters. For simplicity, Schertzer derived the universal multifractal model by nonlinear mixing and densification of different cascade processes with the following formulas:

where *α* denotes the multifractal intensity and *C*_{1} is the codimension of the information dimension (Schertzer *et al.*, 1987).

RESULTS AND DISCUSSION

RESULTS AND DISCUSSION

The Statistical Features of the Air Pollution Time Series

The Statistical Features of the Air Pollution Time Series

The Statistical Features of the Air Pollution Time Series

The possible divergence of high moments for statistical quantities about heavy-tailed random variables indicates the inefficiency and instability on the estimation or analysis by using traditional probability theory, and thus a new probabilistic model will be adopted to characterize the statistical distributions of the air pollutants. For the extremes or outliers in heavy tail, it is theoretically considered that their variances are infinite, and the sums of them will no longer be normally distributed; instead, they will be random variables of a stable distribution. The probability density function of this model commonly is defined based on four parameters, including stability parameter, *α _{s}*, allowing the tail shape of PDF to be adjusted; skewness,

*β*; scale parameter,

*c*; and position parameter,

*µ*. For a stable distribution, there are three known special cases: When

*α*= 2, it reduces to a normal distribution, while for

_{s}*α*= 1 and

_{s}*β*= 0, and

*α*= 1/2 and

_{s}*β*= 1, it corresponds to a Cauchy distribution and Levy distribution, respectively. Taking the PM

_{2.5}time series of Beijing as an example, the fitting results by using log-normal distribution, stable distribution, normal distribution, and kernel function models are shown in Fig. 3, which exhibit good fittings for log-normal and stable models over the whole range and also reveal that the time series is right-skewed with heavy tail. The features can be extended to the rest of the pollution series of seven cities, and Tables 1 and 2 list the estimates of

*α*and

_{s}*β*for all the pollution time series.

**Fig. 3.** The fitting results of the probability density function of the PM_{2.5} series observed in Beijing using log-normal distribution, normal distribution, kernel function, and Levy distribution.

For a relation of a spectrum versus wave numbers in a periodogram of a time series, which can be built under the Fourier transformation, the curve can reflect a distribution of spectral energy at different scales and also allows to detect periodic characteristics for the time series. Generally speaking, in the curve, an energy peak represents one cycle, clearly having higher value than surrounding areas. Like Fig. 4, the daily and semi-daily cycles are experimentally confirmed to be commonly existing in all of the time series, and like Fig. 5, quite a few time series are also characterized by the presence of lower energy peaks at 8 h, 6 h, and even earlier. It is evident that the daily and semi-daily cycles are closely associated with the forces driven by the sun, which can deeply affect human activities, such as various anthropogenic sources, and meteorological conditions, such as diffusion and deposition factors. However, it is difficult to explain the other energy peaks, which may be involved with the other influencing factors, such as pollutants’ properties (physical or chemical) or complex interactions between natural and anthropogenic environmental conditions, and it is noted that the cycles should not be ignored from their comparison of amplitudes with the daily or semi-daily cycles.

**Fig. 4.** The spectral density function of the PM_{2.5} series observed in Xi’an versus wave numbers.

**Fig. 5.** The spectral density function of the O_{3} series observed in Beijing versus wave numbers.

For a time or space series, the nonstationarity means the values differ with respect to time and location and requires that when making a forecasting or an analysis, the underlying trend need to be removed. The multifractal cascade framework provides a new effective idea to distinguish and identify the feature. As argued by Davis *et al.* (1993), when *β* < 1, it indicates that the process is stationary, and *β* > 1 means the process is non-stationary, while for 1 < *β* < 3, the increment of the process is stationary. Here, using the PM_{2.5} series of Beijing as an example, it can be seen from Fig. 6 that the spectrum can be divided into two parts at the point of 11 days, which is consistent with the conversion scale between low-frequency weather and the weather with a slope, *β*, that is close to 0 at large scales and 5/3 at small scales in the turbulence field (Lovejoy and Schertzer, 2012). Meanwhile, the valley value appears at approximately 60 days, and from the performance of *β*, the series is found to be non-stationary. Through the further investigation, non-stationarity and the crossover between the two regimes corresponding to different scales also generally exist in the rest of pollutants series. Another parameter related to *β* is the Hurst index, *H _{u}*, and there is a relation of

*β*= 2

*H*+ 1, where

_{u}*H*can be used to determine the positive or negative long-range correlation; specifically, when 0 <

_{u}*H*< 1/2, the series is positive, and when 1/2 <

_{u}*H*< 1, it is negative.

_{u}**Fig. 6.** The spectral density function of the PM_{2.5} series in Beijing, averaged over logarithmically spaced bins, 10 per order of magnitude.

Multifractal Analysis of the Air Pollution Time Series

Multifractal Analysis of the Air Pollution Time Series

Multifractal Analysis of the Air Pollution Time Series

Multifractal analysis has gradually become a standard high moment statistical method, which can relate the values of different scales through a scaling function and describe the scale invariance for a time series through multiple dimensions. For the estimation of the multifractal formalism, here we adopted a method called double-trace moment involving two moment operations to a scaling quantity, which has been proven to be more effective and stable than moment method (Lavallée, 1991). As proposed by Schertzer (1987), firstly, a non-stationary series must be converted to a stationary series using the relation Δ*v* ≈ *ε*^{1/3}*κ*^{1/3}, where Δ*v* is a difference in velocities at two arbitrary points with a certain interval in a turbulence field, and *ε* and *κ* are the conservative energy and scaling variables, respectively. Here, the concentration values are observed directly and act as substitutes for the velocity *v*. The specific estimation steps are as follows: Firstly, first-order differences are performed on the minimum scale for a time series. Then, the absolute values are used instead and divided by the mean value. Thus, the processed time series is normalized to obtain a new and conservative time series, *ε _{λ}* = |Δ

*µ*|/<Δ

*µ*>, where

*µ*denotes the concentration values of an air pollution time series. Finally, a series of statistical quantities at various scales defined based on the average for the above result,

*ε*, are subjected to power operations two times with the orders of

_{λ}*η*and

*q*, respectively, to obtain Tr

*λ*(

*ε*)

_{λ}^{η}*, which are related to*

^{q}*l*by Tr

*λ*(

*ε*)

_{λ}^{η}*≈*

^{q}*λ*

^{K}^{(q,η)}, where the functions

*K*(

*q*,

*η*) are the scaling functions of double moments and are related to the moment scaling function

*K*(

*q*) as

*K*(

*q*,

*η*) =

*η*(

^{α}K*q*). The double trace moment estimators are functions with three variables, and thus, the quantities are all three-dimensional arrays. Fig. 7 shows the relations of Tr

*λ*(

*ε*)

_{λ}^{η}*versus*

^{q}*λ*of the PM

_{10}series of Beijing, and we can see that there are good linear relationships between them; in other words, the scale invariance remains within certain intervals, and the values of

*K*(

*q*,

*η*) are accurately estimated by the least square method. Once the values of

*K*(

*q*,

*η*) are known, then

*α*is estimated by the relation

*K*(

*q*,

*η*) =

*η*(

^{α}K*q*), and

*C*

_{1}can also be estimated by Eqs. (1)–(2).

**Fig. 7.** The relation of the logarithm moments of the orders of *q* with the values from 0.1 to 3 with the step of 0.1 versus the logarithm scales for the Beijing PM_{2.5} series.

In the universal multifractal formalism, a continous fractal or multifractal cascade process can be simulated by an exponential operation to a weighted noise with the prior konwn value for levy variable*. *For a process, four types are categorized with respect to *α.* When *α* = 0, it is a mono-fractal; when 0 < *α* < 2 except *α* = 1, it is a log-Levy multifractal; and in the other cases when *α* = 2 and *α* = 1, it corresponds to a log-normal and log-Cauchy multifractal, respectively. In the formalism, *α* is commonly used as a characterization for multifractal intensity and reflect hierarchical structure for a process, while *C*_{1} is a first-order statistic, also called information co-dimension, which reflects sparseness of mean process. For *C*_{1}, having a relation *C*_{1} = 1 – *D*_{1} where D_{1 }is information dimension, which is associated with the uniformity for a distribution of a process over the whole range A more intuitive understanding of the physical meanings for the two parameters can be obtained from Fig. 8, which presents four simulated universal multifractal process with (a)* α* = 0.6 and *C*_{1} = 0.12, (b) *α* = 1.6 and *C*_{1} = 0.12, (c) *α* = 1.2 and *C*_{1} = 0.12, and (d) *α* = 1.2 and *C*_{1} = 0.52 (the specific simulation method: http://www.physics.mcgill.ca/~gang/so

ftware/index.html). According to Figs. 8(a) and 8(b), compared to the lower *α* having more areas with high values, the higher *α* owns fewer but stronger upward peaks, and it indicates that extreme pollution risk will increase. According to Figs. 8(c) and 8(d), compared with lower *C*_{1}, for higher *C*_{1}, the degree of sparseness increases, and some extreme outliers will occur; it seems to make a compression and stretch for the data as a whole.

**Fig. 8.** The comparison of cascade simulations with different values of *α* and *C*_{1}: (a)* α* = 0.6, *C*_{1} = 0.12; (b) *α* = 1.6, *C*_{1} = 0.12; (c) *α* = 1.2, *C*_{1} = 0.12; (d) *α* = 1.2, *C*_{1} = 0.52.

With the underlying significance in a process for the two parameters, next let us discuss the overall characteristics of all the air pollution time series in the seven cities. As seen from Figs. 9 and 10, in addition to Chengdu, the values of *C*_{1} and *α* of the two pollutants PM_{2.5}, PM_{10}, and AQI observed in the northern cities are generally greater than those in the southern cities due to the heating and easily calmed weather in the northern China during winter and also indirectly indicate that the north is the main polluted area of China in winter. It is noted that the maximum values of *α* for PM_{2.5}, PM_{10}, and AQI occur in a western city, Chengdu, which imply that there are stronger pollution peak areas than the other cities during 2016. This result may be caused by Chengdu’s unique natural and geographical environment. Chengdu is located in the Sichuan Basin in western China, a very closed geographic unit, and is perennially characterized by low wind speed, high cloudiness and humidity, more neutral and stable weather conditions, and so forth; therefore, the inversion occurs frequently with high intensity, and the pollutants diffuse poorly. With the higher values of *C*_{1,PM2.5}, *C*_{1,PM10}, and *C*_{1,AQI} for Beijing, it can be concluded from the time series that the extreme pollutions yielded from the sparseness occurred on February 8 and May 6, 2016. The former is caused by the calm weather during the heating period, which results in a worsening of pollutant diffusion conditions and is afterwards coupled with the migration from other areas, yielding the extreme pollution; moreover, the latter is caused by the dust storms and reaches extreme values in a short time. The pollution for the three pollutants in Shanghai, Guangzhou, andWuhan are relatively less serious; however, Guangzhou and Wuhan need to face certain pollution risks from PM_{2.5} and PM_{10}, respectively. Shenyang has the largest *α*_{SO2} and the lowest *C*_{1,SO2}, i.e., it has the strongest pollution peaks and the smallest sparse degree for SO_{2}, indicating that the pollutant consistently maintain at a high level of pollution, and an irrational energy structure, especially coal-fired heating, is the main reason. Wuhan has the highest level of O_{3} pollution, which is also identified from the measured series with high concentration values in July and August, and it is mainly caused by automobiles and marine exhausts, and volatile organics used in spray paintings, dyes, and petrochemical productions, which are processes that also lead to the larger value for *C*1,_{CO}, while the pollution level of O_{3} in Beijing is relatively low. The maximum values of *α*_{CO} and *α*_{NO2} occur in Xi’an, and the maximum value of *C*_{1,NO2} occurs in Shenyang, and the pollutants mainly come from vehicle exhaust, coal-fired boiler heating, and thermal power generation. Contrary to the other pollutants, *α*_{CO} and *C*_{1,CO} roughly show a negative correlation, indicating that there are fewer pollution peaks in southern cities, such as Shanghai, Wuhan, and Guangzhou, but they largely deviate from normal levels, i.e., meaning a rapid increase of CO in a short time against a stable background, and the case could also cause greater harm to human health and ecosystem. This is also the reason why some southern cities have listed automobile exhaust as their primary pollutant.

**Fig. 9.** The comparison of the *α* values of six pollutants and AQI series in seven cities.

**Fig. 10.** The comparison of the *C*_{1} values of six pollutants and AQI series in seven cities.

Although the above scaling behaviors of the pollutants were observed through the two parameters *α* and *C*_{1}, among which *C*_{1} is a first-order moment parameter, and although *α* reflects a high-order moment characteristic of the time series, it does not intuitively provide the scaling features of all the moments. It should be noted that the codimension function *c*(*γ*) refers to structures that can be described as a collection of interwoven fractal subsets that exhibit power-law behavior with a range of scaling exponents. A comparison of the estimates of *c*(*γ*) for the AQI series among different cities is shown in Fig. 11, where the axis of abscissa denotes singularity index *γ*, representing a sample space with the concentration values constantly increasing from left to right; on the other hand, the ordinate axis shows the estimates of codimension *c*(*γ*), i.e., the higher it is, the larger the sparseness degree of the time series is. From the figure, it can be concluded that the highest peak of AQI occurs in Beijing by the maximum estimates of *γ*, indicating that it is exposed to the most serious extreme pollution episode, and in contrast, the risks of extreme pollution in Guangzhou, Shanghai, and Wuhan maintain at a relatively low level. Similarly, through the comparisons of the other pollutants among different cites, more major and prominent pollutants can be identified and deduced for each city: Beijing, Chengdu, and Xi’an each have three pollutants, i.e., PM_{2.5}-PM_{10}-SO_{2}, NO_{2}-PM_{10}-SO_{2}, and O_{3}-PM_{10}-SO_{2,} respectively; Wuhan has high extreme pollution risks of CO and O_{3}, while the risks in Guangzhou and Shenyang are mainly associated with CO and PM_{10} pollution, respectively.

**Fig. 11.** The estimated *c*(*γ*) functions of the AQI series in the seven cities.

Due to limitations regarding resolutions and sample sizes, *c*(*γ*) is not always continuous in the domain of definition, and there are two kinds of discontinuities. According to the analogies between multifractal formalisms and statistical thermodynamics, the discontinuities are connected with the phase changes in the physical processes, and a first order and a second order phase transition were proposed by Schertzer and Lovejoy (1993) through the study of the cascade processes. For a cascade process, when the sample size is large enough, the discontinuity point, *q _{d}*, can be revealed, which is also called the first order phase transition, and the corresponding theoretical limit values,

*γ*

_{s}_{,d}and

*c*(

*γ*

_{s}_{,d}), can also be obtained. When the sample size is not sufficient, the discontinuity,

*q*, cannot be reached, but the sample discontinuity,

_{d}*q*, can be obtained, which is called the second order phase transition. Therefore, the sample size is very important in the determination of the theoretical limits of

_{s}*γ*,

*c*(

*γ*), and other parameters, especially the moment order,

*q*. Next, the quantitative analysis of the API data of seven cities within the period from 2000 to 2012 is conducted to observe the multifractal phase transitions in air pollution systems.

Before calculation, each series needs to be divided into two time segments, which are set to be 4 years and 13 years in the scenario. The specific calculation process is as follows: Taking a series as an example, firstly, the long-tailed behavior of the survival function for a conservative measures differenced from the series at the minimum scale is fitted to obtain *q _{d}*, and then the scaling moment function,

*K*(

*q*), is fitted using the formula

*K*(

*q*) =

*γ*

_{max}

*q*– Δ

*s*, where

*γ*denotes the maximum of singular values and Δ

_{max}*s*is the effective dimension, which correspond to the slope and intercept of the equation, respectively. Secondly, the parameters

*γ*,

_{D}*q*, and

_{s}*γ*are calculated using the formulas

_{s}*γ*=

_{D}*K*ʹ(

*q*),

_{d}*q*= (Δ

_{s}*s*/C

_{1})

^{1/α}, and

*γ*=

_{s}*K*ʹ(

*q*). The results analyzed as a representative of five cities selected for the two different lengths of data are listed in Table 3. For each city,

_{s}*q*>

_{s}*q*and

_{D}*γ*>

_{s}*γ*indicate that the temporal variations of the API series observed in air pollution systems belong to the first multifractal phase transitions. Theoretically, the difference of both

_{D}*q*and

_{D}*γ*under the different sample sizes are quite small, and

_{D}*γ*increases as sample sizes increase. However, for

_{s}*q*and

_{D}*γ*, a remarkable change appears in Xi’an and Shanghai, which means that the concentrations vary drastically with respect to the changes in time; in contrast,

_{D}*q*and

_{s}*γ*only have a few minor changes in all cities. Moreover, the theoretical extreme

_{s}*γ*

_{s}_{,d}values of these data are equal to the fitted

*γ*

_{max}values of the empirical data. The behaviors on the intervals [

*γ*,

_{D}*γ*] between the first discontinuity points and theoretical maximums are fully characterized, and it also has been proven that

_{s,d}*c*(

*γ*

_{s}_{,d}) ≈ Δ

*s*and oscillates around the theoretical value of 1 + [log(13)/log(365)] = 1.435.

The Downscaling Behaviors of Multifractal Cascade Analysis Applied in the Air Pollutants’ Concentration Data

The Downscaling Behaviors of Multifractal Cascade Analysis Applied in the Air Pollutants’ Concentration Data

The Downscaling Behaviors of Multifractal Cascade Analysis Applied in the Air Pollutants’ Concentration Data

Generally, the pollutants’ concentration time series are always measured at fixed resolutions, and normally, the concentration over a sampling interval varies with resolutions. However, people are often more interested in observing and obtaining the infinite fine-scale concentrations. In general, the adjacent concentrations at finer resolutions are averagely combined to generate the corresponding values at coarser resolutions, where the smoothness undoubtedly results in the loss of some extreme information, but the extremes are not easily recovered and regained by using the traditional imputation and interpolation methods. As mentioned earlier, fractals or multifractals have natural advantages in describing the extreme values based on a hierarchically generating mechanism realized by a power law relation of scaling moment function, and it indicates that they can provide a scale-independent random downscaling tool. From the generation process of a one-dimensional cascade, the parent node is multiplied by a random variable to obtain twice its number of child nodes. With the different theoretically distributed forms of random variables, multiple cascade models were proposed with the specific formulas of the function *K*(*q*) shown in Table 4 (Gupta *et al.*, 1993). Obviously, the key to downscaling is to calibrate the parameters in the theoretical models through the empirical data and then generate the downscaled series using the most appropriate one from the well calibrated models.

In order to examine the dependence on the distributions of the simulations using the six models, Fig. 12 shows a comparison of the fittings to the histogram of the Shanghai PM_{10} series for hourly simulated time series of one year using Monte Carlo method. It can be seen that the log-normal and log-Poisson models, especially the latter, relatively agree well with the original series in the whole, but there are also some deficiencies. For example, the heavy-tailed portrayal of the log-normal distribution model is clearly insufficient, while the slight oscillation appears over the tail for the log-Poisson model. The *α* and binomial models show greater oscillation than the log-Poisson model over almost the whole range, but their envelope curves are approximately consistent with the frequency distribution of the empirical data. For the uniform model, most of the simulation data are close to 0, but a few extreme outliers will occur occasionally. However, due to the limitation of the assignment of the values for the parameters, the *β* model, which is characterized by two peaks, deviates largely from the original series.

**Fig. 12.** The comparisons of fittings of PDF to the Shanghai PM_{10} series by six kinds of simulations generated from different cascade models.

Downscaling is a method as same as an interpolation method with the common purpose of constructing new data, and to be different, the former mainly focuses on inferring high-resolution information from known low-resolution data while the latter is more about the points in the range of a discrete set of known data points. However, generally the technique of the former is mainly based on dynamical and statistical approaches, while the latter mostly depends on the linear or nonlinear function relations built upon the known neighboring data. Next, let us examine the evolution and changes of distributions of air pollutant concentrations across scales based on multifractal downscaling tools and a comparison in relation to the interpolation method, and the Shanghai PM_{10} series is also employed. Firstly, the data with four coarser resolutions of 16 h, 8 h, 4 h, and 2 h are generated by merging adjacent sample points and averaging; then, the above six models and piecewise linear method are used to generate hourly data, and the estimated probability density functions of them using kernel function are shown in Fig. 13. The goodness of fit of the fitted distributions can quantitatively show the similarities between the original and downscaled series, namely, it can determine which type of the cascade models is the most appropriate to downscale the air pollutants’ time series, and the common test methods include mean bias error (MBE), root mean square error (RMSE), and Kolmogorov-Smirnov (*K*-*S*) and χ^{2} methods (Gavriil *et al.*, 2006). Here, the root mean square error is used to test the accuracy of the downscaled series, and the specific implementation are as follows: Firstly, for each time series at the studied resolutions, a quantile of each empirical concentration is obtained from its estimated probability density function; secondly, the concentration values of the downscaled series corresponding to the quantile values in the first step are calculated based on their estimated probability density functions; then, the RMSE results are calculated by using the two concentration pairs. As shown in Table 5, the piecewise linear interpolation has the best overall performance compared with the six cascade models; however, from Figs. 13(a) and 13(b) it has a poor fitting on the range of larger values, while the six models except the uniform model can behave quite well. Therefore, it can be inferred that the fittings of the range of the high concentrations determine the sizes of the errors. Of the six cascade models, the log-normal model has the highest accuracy, although the *α*, *β*, log-Poisson, and binomial models have the similar accuracy. Except the above comparison, as seen in Figs. 14 and 15, from *γ* – *c*(*γ*) functions and the overlay show, we also get the linear interpolation cannot reach the information of extreme values or outliers at the small scales, but these values can be obtained from the downscaling methods. Therefore, we conclude that the multifractal cascade model are ideal and effective downscaling models, and when people need the singular information eliminated because of combining and averaging, or unobservable due to limitations of laboratory equipments experimental conditions, it can be inferred or recovered using the scaling invariance tools.

**Fig. 13.** The comparison of PDF estimated by kernel function for the downscaled series based on the Shanghai PM_{10} series by using the six cascade models and piecewise linear interpolation method.

**Fig. 14.** The comparison of *c*(*γ*) functions of the original series (the Shanghai PM_{10} series), averaged series, and downscaled series at four scales based on the log-normal random cascade model.

**Fig. 15.** The overlay display of the downscaled series using the log-normal random cascade model, original series (the Shanghai PM_{10} series), and the interpolated series.

**CONCLUSIONS**

The universal multifractal model is a framework that is defined by the infinite densifications and multiplications of many independent cascade processes on various scales, which are assumed to follow a stable distribution in an identical manner; hence, we selected a stable distribution as a candidate for fitting air pollution time series. The presence of the lower energy peaks at 8 h, 6 h, and even earlier is worth studying on account of the unignorable amplitudes; however, it is difficult to explain and may be related to factors that influence air pollution, such as physical or chemical pollutant properties or complex interactions between natural and anthropogenic environmental conditions. Two regimes corresponding to low-frequency weather, and weather scales with a slope (*β*) close to 0 and 5/3 for larger and smaller scales, respectively, were observed, confirming self-similarity and nonstationarity. As the mean and the variance, the parameters *α* and *C*_{1} serve as standard statistical quantities that comprehensively describe the properties of the air pollution series, accurately capturing extreme values and peak clustering, and can be used to evaluate air quality. Exhibiting excellent downscaled performance, multiplicative cascade processes can realistically represent air pollution time series, enabling the further exploration of temporal variations and scaling evolution on finer scales.

D**ISCLAIMER**

D

The authors declare that they have no competing financial interests.

ACKNOWLEDGEMENT

ACKNOWLEDGEMENT

The study is supported by Key Science and Technology Program of Henan Province, China (152102310365, 182102310950).