A DFA-based bivariate regression model for estimating the dependence of PM2.5 among neighbouring cities

TL;DR

This paper introduces a DFA-based bivariate regression model that estimates multi-scale dependence between variables, demonstrated through PM2.5 pollution data among neighboring cities, revealing seasonal and scale-dependent influences.

Contribution

The paper presents a novel DFA-based regression model that captures multi-scale dependence, providing new estimators and indices for analyzing complex variable relationships.

Findings

The model accurately depicts dependence between PM2.5 variables.

Dependence varies across seasons and time scales.

New indices effectively evaluate scale-dependent relationships.

Abstract

On the basis of detrended fluctuation analysis (DFA), we propose a new bivariate linear regression model. This new model provides estimators of multi-scale regression coefficients to measure the dependence between variables and corresponding variables of interest with multi-scales. Numerical tests are performed to illustrate that the proposed DFA-based regression estimators are capable of accurately depicting the dependence between the variables of interest and can be used to identify different dependence at different time scales. We apply this model to analyze the PM2.5 series of three adjacent cities (Beijing, Tianjin, and Baoding) in Northern China. The estimated regression coefficients confirmed the dependence of PM2.5 among the three cities and illustrated that each city has different influence on the others at different seasons and at different time scales. Two statistics based on the scale-dependent $t$-statistic and the partial detrended cross-correlation coefficient are used to demonstrate the significance of the dependence. Three new scale-dependent evaluation indices show that the new DFA-based bivariate regression model can provide rich information on studied variables.

Tables

Table 1. Table 1: Partial correlation coefficients and t 𝑡 t -statistics between per two cities of Beijing, Tianjin, and Baoding in four seasons.

	winter	spring	summer	fall
Beijing vs. Tianjin	0.3048	0.3072	0.2625	0.1660
$t$-statistics	${25.8011}^{\ast }$	${26.2714}^{\ast }$	${22.1366}^{\ast }$	${13.2672}^{\ast }$
Beijing vs. Baoding	0.2745	0.4545	0.4815	0.4468
$t$-statistics	${23.0124}^{\ast }$	${41.5216}^{\ast }$	${44.7079}^{\ast }$	${39.3711}^{\ast }$
Tianjin vs. Baoding	0.5570	0.4517	0.3461	0.5992
$t$-statistics	${54.0752}^{\ast }$	${41.1961}^{\ast }$	${30.0154}^{\ast }$	${58.9941}^{\ast }$

Equations

$$Y={\beta }_0+{\beta }_1{X}_1+{\beta }_2{X}_2+ϵ,$$

$$\{\begin{array}{rl}{\hat{\beta }}_1& =\frac{({\sum }_{t=1}^N{x}_{1t}{y}_t)\cdot ({\sum }_{t=1}^N{x}_{2t}^2)-({\sum }_{t=1}^N{x}_{2t}{y}_t)\cdot ({\sum }_{t=1}^N{x}_{1t}{x}_{2t})}{({\sum }_{t=1}^N{x}_{1t}^2)\cdot ({\sum }_{t=1}^N{x}_{2t}^2)-({\sum }_{t=1}^N{x}_{1t}{x}_{2t}{)}^2}\sim \frac{\widehat{{\sigma }_{X_{1}Y}}\cdot \widehat{{\sigma }_{X_2}^2}-\widehat{{\sigma }_{X_{2}Y}}\cdot \widehat{{\sigma }_{X_1{X}_2}}}{\widehat{{\sigma }_{X_1}^2}\cdot \widehat{{\sigma }_{X_2}^2}-{\widehat{{\sigma }_{X_1{X}_2}^2}}^2},\\ {\hat{\beta }}_2& =\frac{({\sum }_{t=1}^N{x}_{2t}{y}_t)\cdot ({\sum }_{t=1}^N{x}_{1t}^2)-({\sum }_{t=1}^N{x}_{1t}{y}_t)\cdot ({\sum }_{t=1}^N{x}_{1t}{x}_{2t})}{({\sum }_{t=1}^N{x}_{1t}^2)\cdot ({\sum }_{t=1}^N{x}_{2t}^2)-({\sum }_{t=1}^N{x}_{1t}{x}_{2t}{)}^2}\sim \frac{\widehat{{\sigma }_{X_{2}Y}}\cdot \widehat{{\sigma }_{X_1}^2}-\widehat{{\sigma }_{X_{1}Y}}\cdot \widehat{{\sigma }_{X_1{X}_2}}}{\widehat{{\sigma }_{X_1}^2}\cdot \widehat{{\sigma }_{X_2}^2}-{\widehat{{\sigma }_{X_1{X}_2}^2}}^2},\end{array}$$

$$\{\begin{array}{r}var({\hat{\beta }}_1)=\frac{({\sum }_{t=1}^N{x}_{2t}^2)\cdot \frac{{\sum }_{t=1}^N{\hat{e}}_t^2}{N-3}}{({\sum }_{t=1}^N{x}_{1t}^2)\cdot ({\sum }_{t=1}^N{x}_{2t}^2)-({\sum }_{t=1}^N{x}_{1t}{x}_{2t}{)}^2}\sim \frac{1}{N-3}\cdot \frac{\widehat{{\sigma }_{X_2}^2}\cdot \widehat{{\sigma }_{ϵ}^2}}{\widehat{{\sigma }_{X_1}^2}\cdot \widehat{{\sigma }_{X_2}^2}-{\widehat{{\sigma }_{X_1{X}_2}^2}}^2},\\ var({\hat{\beta }}_2)=\frac{({\sum }_{t=1}^N{x}_{1t}^2)\cdot \frac{{\sum }_{t=1}^N{\hat{e}}_t^2}{N-3}}{({\sum }_{t=1}^N{x}_{1t}^2)\cdot ({\sum }_{t=1}^N{x}_{2t}^2)-({\sum }_{t=1}^N{x}_{1t}{x}_{2t}{)}^2}\sim \frac{1}{N-3}\cdot \frac{\widehat{{\sigma }_{X_1}^2}\cdot \widehat{{\sigma }_{ϵ}^2}}{\widehat{{\sigma }_{X_1}^2}\cdot \widehat{{\sigma }_{X_2}^2}-{\widehat{{\sigma }_{X_1{X}_2}^2}}^2}.\end{array}$$

$$R^2=1-\frac{{\sum }_{t=1}^N{\widehat{e_t}}^2}{{\sum }_{t=1}^N{y}_t^2}=1-\frac{\widehat{{\sigma }_{ϵ}^2}}{\widehat{{\sigma }_Y^2}},$$

$$\beta_{j}^{*}=\hat{\beta_{j}}\sqrt{\frac{\sum^{N}_{t=1}x_{jt}^{2}}{\sum^{N}_{t=1}y_{t}^{2}}}\qquad\mbox{for $j=1$ and $2$}$$

$$\eta_{j}=\hat{\beta_{j}}\frac{\langle X_{j}\rangle}{\langle Y\rangle}\qquad\mbox{for $j=1$ and $2$,}$$

$$f_Z^2(n,j)=\frac{1}{n}\sum _{k=1}^n(Z_{j,k}-\widehat{Z_{j,k}}{)}^2.$$

$$F_Z^2(n)=\frac{1}{2N_n}\sum _{j=1}^{2N_n}{f}_Z^2(n,j),$$

$$f_{Z_1{Z}_2}^2(n,j)=\frac{1}{n}\sum _{k=1}^n(Z_{1j,k}-\widehat{Z_{1j,k}})(Z_{2j,k}-\widehat{Z_{2j,k}}),$$

$$F_{Z_1{Z}_2}^2(n)=\frac{1}{2N_n}\sum _{j=1}^{2N_n}{f}_{Z_1{Z}_2}^2(n,j).$$

$${\rho }_{DCCA}(Z_1,Z_2,n)=\frac{F_{Z_1{Z}_2}^2(n)}{\sqrt{F_{Z_1}^2(n)F_{Z_2}^2(n)}}.$$

$${\rho }_{PDCCA}(Z_1,Z_2,n)=-\frac{C_{j_1,j_2}(n)}{\sqrt{C_{j_1,j_1}(n)C_{j_2,j_2}(n)}},$$

$$\left\{\begin{array}[]{rcl}\hat{\beta_{1}}^{DFA}(n)=\frac{F_{X_{1}Y}^{2}(n)\cdot F_{X_{2}}^{2}(n)-F_{X_{2}Y}^{2}(n)\cdot F_{X_{1}X_{2}}^{2}(n)}{F_{X_{1}}^{2}(n)\cdot F_{X_{2}}^{2}(n)-[F_{X_{1}X_{2}}^{2}(n)]^{2}},\\ \hat{\beta_{2}}^{DFA}(n)=\frac{F_{X_{2}Y}^{2}(n)\cdot F_{X_{1}}^{2}(n)-F_{X_{1}Y}^{2}(n)\cdot F_{X_{1}X_{2}}^{2}(n)}{F_{X_{1}}^{2}(n)\cdot F_{X_{2}}^{2}(n)-[F_{X_{1}X_{2}}^{2}(n)]^{2}}.\end{array}\right.$$

$$\widehat{e_t}(n)=Y_t-{\hat{\beta }}_1^{DFA}(n)X_{1t}-{\hat{\beta }}_2^{DFA}(n)X_{2t}-⟨Y_t-{\hat{\beta }}_1^{DFA}(n)X_{1t}-{\hat{\beta }}_2^{DFA}(n)X_{2t}⟩$$

$$\left\{\begin{array}[]{rcl}var(\hat{\beta_{1}}^{DFA}(n))=\frac{1}{N-3}\cdot\frac{F_{X_{2}}^{2}(n)\cdot F_{\epsilon}^{2}(n)}{F_{X_{1}}^{2}(n)\cdot F_{X_{2}}^{2}(n)-[F_{X_{1}X_{2}}^{2}(n)]^{2}},\\ var(\hat{\beta_{2}}^{DFA}(n))=\frac{1}{N-3}\cdot\frac{F_{X_{1}}^{2}(n)\cdot F_{\epsilon}^{2}(n)}{F_{X_{1}}^{2}(n)\cdot F_{X_{2}}^{2}(n)-[F_{X_{1}X_{2}}^{2}(n)]^{2}}.\end{array}\right.$$

$$R_{DFA}^2(n)=1-\frac{F_{ϵ}^2(n)}{F_Y^2(n)},$$

$$\beta_{j}^{*DFA}(n)=\hat{\beta_{j}}^{DFA}(n)\sqrt{\frac{F_{X_{j}}^{2}(n)}{F_{Y}^{2}(n)}}\qquad\mbox{for $j=1$ and $2$,}$$

$$\eta_{j}^{DFA}=\hat{\beta_{j}}^{DFA}(n)\frac{\langle X_{j}\rangle}{\langle Y\rangle}\qquad\mbox{for $j=1$ and $2$.}$$

Full text

Keywords: bivariate DFA regression; time-scale; regression coefficient; partial DFA coefficient.

1 Introduction

In recent years, air pollution has become a more and more serious problem around the world. The new air quality model presented by the World Health Organization in 2016 confirmed that 92% of the world’s population lives in areas where air quality levels exceed their limits [1]. Fortunately, more and more governments have realized the importance of managing air pollution and some actions have been placed. Nowadays, a common topic around the world is the governance of the air pollution source such as smog (the main ingredient is fine particulate matter). Many researchers have been involved in the study on the cause and propagation of smog [2, 3, 4, 5, 6, 7, 8, 9]. Modern statistic methods provide some new perspectives to assess smog trends and propagation characteristics. Among them, most studies have focused on studying the correlations among various air pollution indicators including air pollution index (API), air quality index (AQI), fine particulate matter of PM2.5 (diameter $\leq 2.5\mu m$) concentrations, and PM10 (diameter $\leq 10\mu m$) concentrations, and very limited studies considered the correlations among neighboring areas. A common sense is that smog produced at one source place can spread to surrounding areas [6, 7, 8, 9, 10]. Therefore, it is more practical to explore the dependence of air pollution indicators among adjacent cities as it helps assess the causes of local smog and its spread behavior. It has been found by a newly proposed time-lagged cross-correlation coefficient in Ref. [10] that there are different degrees of correlation for PM2.5 series between four neighboring cities in Northern China. However, what has not been investigated is how the PM2.5 series of one city depends on those of the neighbouring cities. In this work, we will develop a detrended fluctuation analysis (DFA)-based bivariate regression model to investigate this dependence.

The simplest and maturest method to describe the dependence of variables is the linear regression. However, the information gained from the traditional linear regression cannot fully meet our need of investigation on the dependence among different variables at different time periods. On the other hand, note that the DFA proposed in 1990s [11, 12] performs excellently in analyzing the long-range correlations [13] of a nonstationary series with fractality and multifractality [14, 15] at different time-scales. To obtain the cross-correlation between two nonstationary series, DFA was extended to the detrended cross-correlation analysis (DCCA) [16]. By defining scale-dependent detrended fluctuation functions, the methods of DFA and DCCA together with their extensions have been applied in a wide range of disciplines [17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30]. Since the ordinary least squares (OLS) method expresses the estimated parameters of standard regression framework as a form of variances and covariances, it builds a bridge between the regression framework and the family of DFA/DCCA as the latter can also produce variances and covariances. Then, the idea of estimating multiple time scale regression coefficients can be achieved by the DFA/DCCA. Recently, Kristoufek [31] constructed a simple DFA-based regression framework exactly by this bridge. The selected examples show the relationship between the pair of variables varies strongly across scales.

In this work, we focus on the interaction of PM2.5 series of three adjacent cities in Northern China, namely, Beijing, Tianjin, and Baoding. The three cities form a triangle shape in the map. The distances between Beijing and Tianjin, Beijing and Baoding, and Tianjin and Baoding are about 115km, 140km, and 150km, respectively. All three cities have a population of more than 10 million and have been greatly affected by heavy smog in recent years. The real-time data of PM 2.5 series of these three cities from December 1, 2013 to November 30, 2016 are chosen for our consideration, which are taken from the Ministry of Environmental Protection of the People’s Republic of China (http://datacenter.mep.gov.cn). The original data show an obvious periodic characteristic and roughly similar trends among the three cities, which imply that there is a possible relevance between per two cities of them. To verify that, the partial correlation technique is employed to get the intrinsic relations between two cities by deleting the interference from the third variable. Four seasons, classified as winter (December, January, and February), spring (March, April, and May), summer (June, July, and August), and fall (September, October, and November), are considered. The results are listed in Table 1.

In Table 1, we also list the $t$-statistics ($t=r_{12,3}\sqrt{\frac{N-3}{1-r_{12,3}^{2}}}$, where $r_{12,3}$ denotes the partial correlation coefficient between the first and second variables eliminating the effects of the third one, $N-3$ is the degree of freedom) of the partial correlation coefficients to access the statistical significance at the given significance level. Unsurprisingly, Table 1 shows that the correlations of PM2.5 between both per two cities are of statistical significance. It explains that the air quality in one city of Northern China cannot be irrelevant to that of its neighbouring cities, which implies potential dependence among the three cities. However, we also note in Table 1 that the degree of relevance is different among different cities and in different seasons though all of them are significant.

To fully detect and quantify the dependence among the PM2.5 series of the above-mentioned three cities, in this work, we construct a new bivariate regression framework which prevails the DFA method and allows us to investigate the dependence of three nonstationary series with multiple time scales. With the DFA-based variance instead of the standard variance, this new DFA bivariate regression model provides more information on the dependence among variables at different time scales. We organize the rest of this paper as follows. The performance of the proposed DFA regression model and the results on the application to PM2.5 series analysis are reported and discussed in the following section, which is followed by our conclusions. The methodologies including the standard regression method, the DFA method, and the DFA-based regression method are introduced at the end of this paper.

2 Results and Discussions

Performance of DFA estimators

The bivariate DFA-based regression model produces two time scale-based regression coefficients. This allows us to detect the dependence of a response variable and two independent variables at different time scales. In order to examine the validity of the model and show its advantages, in this section, we perform two numerical tests on the non-stationary bivariate regression frameworks $Y=\beta_{0}+\beta_{1}X_{1}+\beta_{2}X_{2}+\epsilon.$

In the first test, we investigate the performance of the DFA estimators under different levels of long-term dependence in $X_{1}$, $X_{2}$, and $Y$. According to  [31], the setting I is given as below: two artificial series $X_{1}$ and $X_{2}$ with length $10000$ are generated by ARFIMA($0,d,0$) process with identical fractional integration parameter ($d$) and independent Gaussian noises ($\xi_{i}(t),i=1$ and 2 ) as $X_{i}(t)=\sum_{n=0}^{\infty}a_{n}(d)\xi_{i}(t-n)$. The quantity $a_{n}(d)$ is defined by $a_{n}(d)=\Gamma(n-d)/[\Gamma(-d)\Gamma(n+1)]$, where $\Gamma(\cdot)$ is the Gamma function. The error-term $\epsilon$ is set as a standard Gaussian noise so that the response variable $Y$ has the same parameter $d$ as the two independent variables. The regression coefficients are set as $\beta_{0}=\beta_{1}=1$ and $\beta_{2}=2$. Fig. 1 shows mean values and standard deviation of the two DFA estimators $\beta_{i}^{DFA}$ ($i=$1 and 2) for the generated series with $d$ ranging from $-0.5$ to $0.5$ (at the step size of $0.1$). The estimators are averaged over scales between $10$ and $1000$ with a logarithmic isometric step. Each case is run $1000$ times to eliminate the noise interference. It is clear that the two estimators locate the two given regression coefficients of $1$ (Fig. 1a) and 2 (Fig. 1b) unbiasedly, and are independent of the value of $d$. In addition, the standard deviations of both estimators decrease with the increasing memory. The good performance shows that the method is feasible. On the other hand, to investigate the performance of the DFA estimators faced with a long-range dependent error-term $\epsilon$, we use setting II given as: the memory parameter $d$ is fixed at 0.4 for both $X_{1}$ and $X_{2}$, and the $\epsilon$ is produced by an ARFIMA process with $d_{\epsilon}$ varying from $-0.5$ to $0.5$. Other settings are as those in setting I. Fig. 2 records similar information as that in Fig. 1. Although the fluctuation of DFA estimators increases with $d_{\epsilon}$, which is expected due to an increasing weight of the error-term in the dynamics of $Y$ with the increasing memory of the error-term, we are satisfied to find that the two estimators are still unbiased pointing to the given values with a narrow range for each level of memory of the error-terms.

Our second numerical test aims to show that the DFA estimators are able to identify the dependence of studied variables at different time scales whereas the classical method cannot. To this end, a binomial multifractal series (BMFs) is employed to be regarded as the independent variable $X_{1}$, which is constructed as $X_{1}=p^{n-n[k-1]}(1-p)^{n[k-1]}$, $k=1$, $2$, $\ldots$, $2^{n}$, where the parameter $p\in(0,0.5)$ (We take $p=0.3$ in our test), $n[k]$ denotes the number of digit $1$ in the binary representation of the index $k$. The variable $X_{2}$ is a Gauss variable with [math] mean and $0.0001$ standard deviation. Both $X_{1}$ and $X_{2}$ are of length $2^{15}$. The bivariate regression framework $Y=\beta_{0}+\beta_{1}X_{1}+\beta_{2}X_{2}+\epsilon$ is set with the same coefficients as the first test ($\beta_{0}=\beta_{1}=1$ and $\beta_{2}=2$). The error-term $\epsilon$ is the Gauss noise of the same strength as $X_{2}$. For the BMFs $X_{1}$, we remove all values smaller than $0.00001$ so that only a few of the largest elements are left. In their places, we substitute Gaussian distributed random numbers with [math] mean and $0.0001$ standard deviation. Then we obtain a binomial cascade series embedded in random noise. We analyze the dependence between the response variable $Y$ and two independent variables and find that the estimated $\beta_{2}^{DFA}$ is unbiased at $2$ with a desirable error bar for every time scale, as shown in Fig. 3. However, the performance of $\beta_{1}^{DFA}$ has changed a lot. The dependence between $Y$ and $X_{1}$ is obviously less than the given value at the smaller scales contrary to the larger ones. This is because in the smaller scales, the dependency has been destroyed by the random noise. Our DFA estimators have the capability to recognize this effect while the classical estimators fail to do so (see the errorbar with circle symbol in Fig. 3).

Performance of the three models’ regression coefficients

As mentioned above, air pollution in Northern China is very serious in recent years. Fine particulate matter from industrial exhaust and smoke dust forms smog to fill in the air. We now apply our DFA regression model to investigate the dependence of PM2.5 series in these three cities. We build three bivariate models for Beijing, Tianjin, and Baoding, respectively. In Model I, the dependent variable ($Y$) is the PM2.5 of Beijing while the two independent variables are the PM2.5 of Tianjin ($X_{1}$) and Baoding ($X_{2}$); in Model II, $Y$ is the PM2.5 of Tianjin, $X_{1}$ is the PM2.5 series of Beijing and $X_{2}$ is the PM2.5 series of Baoding; in Model III, $Y$ is the PM2.5 of Baoding, $X_{1}$ and $X_{2}$ stand for the PM2.5 series in Beijing and Tianjin, respectively. In this section, we first show the performance of the regression coefficients at different scales in the three models and then make two statistical tests for the two regression coefficients in each model. Some evaluations for the DFA-based regression and the standard regression are conducted at the end of this section.

The two regression coefficient estimators together with their standard deviations of the three models are sketched in Figs. 4–6, respectively. As expected, the effect is obviously positive. However, a strong variation across scales is found in different seasons. More specifically,

• (a)

In the Beijing’s model, Tianjin ($X_{1}$) has strongly positive effect in every season, especially for the larger time scales. On the contrary, Baoding ($X_{2}$) has different effects on Beijing. Compared to spring and summer, the effect is quite weak in the other two seasons, especially in winter, $\beta_{2}^{DFA}(n)$ is nearly [math] when the scale is more than $800$ hours.

• (b)

In the Tianjin’s model, Baoding ($X_{2}$) presents more unstable effect at different scales. Particularly in summer, $\beta_{2}^{DFA}(n)$ is close to [math] from the smaller scale to the larger scale at about $50$ days ($1200$ hours), which implies that the positive correlation between Tianjin and Baoding can last less than $50$ days. In addition, the two coefficients are less than $0.5$ in most days, which indicates that Beijing and Baoding have little impact on the PM2.5 in Tianjin.

• (c)

For the model of Baoding, the effect of Tianjin ($X_{2}$) to Baoding is similar to that of Baoding to Tianjin in model II. However, the fact that after approximately $17$ days ($408$ hours) the effect reaches the value greater than $1$ indicates that an increase of $1$ unit PM2.5 concentration of Tianjin will lead to the increase of more than $1$ unit PM2.5 concentration in Baoding. In this regard, Tianjin has more impact on Baoding. In addition, the narrow confidence intervals and low standard deviations (less than $0.02$) shown in all sub-plots suggest satisfied reliability of the estimates.

Statistic significance tests of regression coefficients

As mentioned above, the estimated $\hat{\beta}^{DFA}(n)$ is able to theoretically describe the dependence between the impulse variables and the response variables at different time scales. In theory, as long as $\hat{\beta}^{DFA}_{j}(n)$ is not equal to zero, the independent variable $X_{j}$ will affect $Y$. However, for finite time series, $\hat{\beta}^{DFA}_{j}(n)$ is not always equal to [math] even in the absence of relationship between $X_{j}$ and $Y$ due to the size limitation. Therefore, we perform a hypothesis test for the estimated $\hat{\beta}^{DFA}(n)$ to ensure the significance. The standard regression analysis provides a so-called $t$ statistic defined as $t_{j}=\frac{\hat{\beta}_{j}-\beta_{j}}{\sqrt{var(\hat{\beta}_{j})}}$ ($j=1$, 2) for this purpose. We have $t_{j}\sim t(N-3)$ for the bivariate regression model as $\hat{\beta}_{j}\sim N(\beta_{j},var(\hat{\beta}_{j}))$. In general, if $\mid t_{j}\mid>t_{1-\alpha/2}(N-3)$ with a given $\alpha$, we should reject the null hypothesis of $\beta_{j}=0$ and the dependence between $X_{j}$ and $Y$ is considered to be statistically significant. However, since lots of time scales are taken accounted in the DFA regression model, using a single critical value of $t_{1-\alpha/2}(N-3)$ is inappropriate. A correct way is to generate a critical value $t^{c}(n)$ for each time scale. To this end, inspired by the idea proposed by Podobnik et al. [32], we shuffle the considered PM2.5 series and repeat the DFA regression calculations for $10,000$ times. Then let the integral of probability distribution function (PDF) from $-t^{c}(n)$ to $t^{c}(n)$ be equal to $1-\alpha$ (here, we take $\alpha=0.01$). As an example, we show the PDF of $t^{c}(n)$ with five given $n$’s produced by the shuffled PM2.5 series of fall in Fig. 7.

As expected, the symmetrical PDF of $t^{c}(n)$ converges to a Gaussian distribution according to the central limit theorem. In addition, the critical value increases as $n$ increases. This implies that large time scale may strengthen dependence between two variables. By using $t^{c}(n)$, we can determine whether the dependence between the impulse variable and the response variable is significant or not. In practice, the dependence between $X_{j}$ and $Y$ is present when $t_{j}(n)(=\frac{\hat{\beta}_{j}^{DFA}(n)-\beta_{j}}{\sqrt{var(\hat{\beta}_{j}^{DFA}(n))}})$ is larger than $t^{c}(n)$. For the four seasons, the scale-dependent $t$-statistics of regression coefficient together with the scale-dependent critical value $t^{c}(n)$ are presented in Fig. 8.

Note that in Model I (for Beijing), the $t(n)$-statistics of $\beta_{1}^{DFA}(n)$ (Tianjin’s coefficient) is equal to that of $\beta_{1}^{DFA}(n)$ (Beijing’s coefficient) in Model II (for Tianjin), the $t(n)$-statistics of $\beta_{2}^{DFA}(n)$ (Baoding’s coefficient) is equal to that of $\beta_{1}^{DFA}(n)$ (Beijing’s coefficient) in Model III (for Baoding), and in Model II, the $t(n)$-statistics of Baoding’s coefficient $\beta_{2}^{DFA}(n)$ is equal to that of Tianjin’s coefficient $\beta_{2}^{DFA}(n)$ in Model III (for Baoding). Here the three colored lines with different symbols represent the $t(n)$-statistics between each per two cities while the black dashed line stands for $t^{c}(n)$. The partial DCCA coefficient $\rho_{PDCCA}(n)$ is recently developed to uncover the intrinsic relation for two nonstationary series at different time scales. We also calculate the partial DCCA coefficients $\rho_{PDCCA}(n)$ of Beijing and Tianjin, Beijing and Baoding, and Tianjin and Baoding, respectively, and present the results in Fig. 9. For the same purpose of testing the statistical significance, we also produce a critical value for the four seasons. Similarly, the PM2.5 data are shuffled $10,000$ times in the PDCCA calculations repeatedly, and thus $\rho_{PDCCA}^{c}(n)$ for 99% confidence level is obtained, which is also shown in Fig. 9.

Comparing results in Fig. 8 and Fig. 9 gives amazing similarities, which are also in agreement with the results shown in Figs. 4–6. Based on the results, we can draw the following three main points.

• (a)

The dependence between Beijing and Tianjin (the blue square line) gradually increases with the increasing time scales in all seasons. However, the dependence between the two cities is lower than other cities. This finding uncovers that the reason for the serious air pollution in these two cities are mainly due to their own heavy smog or are impacted by other cities.

• (b)

The dependence between Beijing and Baoding (the green triangle line) is significant in spring, summer, and fall. In winter, however, the dependence disappears at long time scale, which implies that the two cities can only affect each other at relatively short term. Moreover, compared to winter and fall, the dependence is much stronger in spring and summer, especially at long time scales, which indicates that they affect much longer in warm weather.

• (c)

In spring and summer, the $t(n)$-statistics and $\rho_{PDCCA}(n)$ of Tianjin vs. Baoding (the red circle line) go down through the critical lines of $t^{c}(n)$ and $\rho^{c}_{PDCCA}(n)$, respectively at about $800$ hours. This suggests that the dependence between Tianjin and Baoding will disappear when it’s more than one month. However, the exact opposite occurs in winter and fall. In these two seasons, both $t(n)$-statistic and $\rho_{PDCCA}(n)$ increase with the increasing time scales, which demonstrates that the interaction of bad air quality between the two cities will last longer in cold days.

Evaluations of DFA-based regression model

To evaluate our estimated DFA-based bivariate regression model, we plot the scale-dependent determination coefficient $R^{2}_{DFA}(n)$, and the beta coefficient $\beta^{*DFA}(n)$ and the average elasticity coefficient $\eta^{DFA}(n)$ in Fig. 10, and Figs. 11–13, respectively.

To show the new model provides more information than the standard regression model does, we also include the three corresponding coefficients of standard bivariate regression model in these figures. As seen from Fig. 10 that $R^{2}_{DFA}(n)$ is superior to the standard $R^{2}$ at most time scales. The good performance illustrates that one will gain richer information in explaining the response variable when using our DFA-based regression model. On the other hand, we can conclude from Figs. 11–13 that (1) Baoding has more influence than Tianjin on Beijing in all seasons except for winter. (2) Tianjin is more sensitive to Baoding’s changes in air quality than Beijing’s in winter and fall. (3) Tianjin affects Baoding more than Beijing does in winter and fall, but less in the other two seasons. In addition, Figs. 11–13 illustrate that the standard $\beta^{*}_{j}$ and $\eta_{j}$ can be seen as the mean values of the DFA-based $\beta^{*DFA}_{j}(n)$ and $\eta^{DFA}_{j}(n)$, respectively. This means that $\beta^{*DFA}_{j}(n)$ and $\eta^{DFA}_{j}(n)$ are able to measure the dependence degree of the studied independent variable on the dependent variable in all directions. Thus one can access the measurement according to his/her needs. For example, in winter of Model I, we find that the $\beta^{*DFA}_{2}(n)$ and $\eta^{DFA}_{2}(n)$ are larger than $\beta^{*DFA}_{1}(n)$ and $\eta^{DFA}_{1}(n)$, respectively, at smaller scales but much smaller at larger scales, which shows that the sensitivity of $Y$ to $X_{2}$ (Baoding) is greater than that of $Y$ to $X_{1}$ (Tianjin) for short term ($\leq 300$ hours) but Tianjin is more sensitive to Beijing at the long term. This can help air quality inspectors make the correct analysis for Beijing’s PM2.5 at different periods.

3 Conclusions

The study of dependence between variables helps expose the causal relationship and correlation of the variables of interest in the real world. The linear regression model is undoubtedly one of the simplest methods among many approaches. However, single variety of regression coefficient and evaluation index cannot show all aspects of the dependence between independent variables and dependent variable. As a meaningful extension, we design a new framework for bivariate regression model using the prevailing DFA method. The proposed bivariate DFA regression model allows us to estimate multi-scale regression coefficients and other corresponding scale-dependent evaluation indicators. It has been shown via two artificial tests that these DFA-based regression coefficients are able to describe the dependence between the response variable and two independent variables exactly; and can capture different dependence at different time scales.

An application of the new model to the study of dependence of PM2.5 series among three heavily air polluted cities in Northern China unveils that huge difference of the dependence exists in per two cities in different seasons and at different periods. Three new indicators of the scale-dependent determination coefficient, the scale-dependent beta coefficient, and the scale-dependent elasticity coefficient are proposed, which turned out to be more practical than those in standard regression models. Three main points can be concluded as (1) Beijing and Baoding have little impact on the PM2.5 in Tianjin while Tianjin takes more impact on Baoding and the air quality of Beijing is more sensitive to the changes in Baoding. (2) In contrast, the air quality in Beijing and Tianjin is not significantly relevant, while the air quality in Tianjin and Baoding has a very significant impact on each other especially in the cold weather. (3) In comparison, the fluctuation of PM2.5 in Baoding has the greatest impact on the other two cities in most days. While Baoding’s air quality is more sensitive to Beijing’s changes in spring and summer, and is more sensitive to Tianjin’s changes in winter and fall. These findings may provide some useful insights on understanding air pollution sources among cities in Northern China.

4 Methods

The standard bivariate regression model

To study the dependence of air quality among three neighboring cities, we consider a bivariate linear regression model as

[TABLE]

where $Y$ is a dependent variable, $X_{1}$ and $X_{2}$ are two independent variables, $\epsilon$ is a Gaussian error term with zero mean value, and $\beta_{j}$ ($j=1$, $2$) is the partial regression coefficient characterizing the dependence on $X_{j}$. The most critical work in empirical studies is to estimate $\beta_{1}$ and $\beta_{2}$. The OLS method gives

[TABLE]

where $\langle\cdot\rangle$ denotes the mean value of the whole time period, $x_{1t}=X_{1t}-\langle X_{1}\rangle$, $x_{2t}=X_{2t}-\langle X_{2}\rangle$, and $y_{t}=Y_{t}-\langle Y\rangle$. Then the estimator of residuals can be determined by $\hat{e}_{t}=Y_{t}-\hat{\beta}_{1}X_{1t}-\hat{\beta}_{2}X_{2t}-\langle Y_{t}-\hat{\beta}_{1}X_{1t}-\hat{\beta}_{2}X_{2t}\rangle$. With it one can obtain the estimators of variance of the two regression coefficients as

[TABLE]

The variance illustrates the accuracy of the estimated parameters. The estimated regression coefficients together with their variances can be further employed for hypothesis test and model evaluation. As an important indicator to evaluate the regression model, the determination coefficient $R^{2}$ is defined by

[TABLE]

with the range of $[0,1]$. $R^{2}$ measures a proportion of variance of $Y$ explained by $X_{1}$ and $X_{2}$ and higher value of $R^{2}$ implies better model interpretation ability. Moreover, to quantify sensitivity of explained variable to each explaining variable, two quantities, namely, the beta coefficient (denoted as $\beta_{j}^{*}$) and the average elasticity coefficient (denoted as $\eta_{j}$), are defined

[TABLE]

and

[TABLE]

which can explain the relative importance of variables $X_{1}$ and $X_{2}$ to $Y$. According to [31], the advantage of translating the standard notation into variance and covariance shown on the right-hand side of Eqs. (4.2)-(4.5) is available to use the DFA/DCCA methods based on the same idea.

The DFA-based variance and DCCA-based covariance functions

DFA and DCCA methods are described as follows. For a time series $\{z_{t}\}$, $t=1$, 2, $\cdots$, $N$, we split its profile $Z_{t}=\sum_{i=1}^{t}(z_{i}-\langle z\rangle)$ into $N_{n}=[N/n]$ nonoverlapping segments of equal length $n$, denoted as $Z_{j,k}$, $k=1$, 2, $\cdots$, $n$. The same procedure is repeated starting from the opposite end to avoid disregarding a short part of the series in the end and thus $2N_{n}$ segments are obtained altogether. In the $j^{th}$ segment, we have $Z_{j,k}=Z_{(j-1)n+k}$ for $j=1$, 2, $\cdots$, $N_{n}$ and $Z_{j,k}=Z_{N-(j-N_{n})n+k}$ for $j=N_{n}+1$, $N_{n}+2$, $\cdots$, $2N_{n}$, where $k=1$, 2, $\cdots$, $n$. In each segment, the local linear (or other) trend  [33, 34] can be fitted as $\widehat{X_{j,k}}$ (in our work, we use $2^{nd}$ order polynomial to fit the trend in each segment). Fluctuation function $f_{Z}^{2}(n,j)$ is then defined for each segment as

[TABLE]

Averaging the fluctuation $f_{Z}^{2}(n,j)$ over all segments yields

[TABLE]

which is the so-called DFA-based scale-dependent variance function. To obtain the scale-dependent covariance of two equal length series $\{z_{1t}\}$ and $\{z_{2t}\}$, $t=1$, 2, $\cdots$, $N$, we only need to translate the univariate fluctuation function in each segment and average fluctuation into the bivariate case directly,

[TABLE]

[TABLE]

The scale-characteristic fluctuation $F_{Z_{1}Z_{2}}^{2}(n)$ is the so-called DCCA-based covariance, which expresses the cross-correlation fluctuations between the series of $\{z_{1t}\}$ and $\{z_{2t}\}$. Thus we have obtained all objects to create the DFA-based regression model. But for purpose of testing, we need some accessories of the DFA process. The DCCA cross-correlation coefficient $\rho(n)$, proposed by Zebende [35], can measure the cross-correlation between two nonstationary series at multiple time scales, which is defined as

[TABLE]

To access intrinsic relations between the two time series on time scales of $n$, Yuan et al. [36] and Qian et al. [37] developed a so-called partial DCCA coefficient independently, which applies partial correlation technique to delete the impact of other variables on the two currently studied variables. This coefficient is defined as

[TABLE]

where $C$ is the inverse matrix of the cross-correlation matrix produced by $\rho_{DCCA}(n)$ of $Z_{1}$, $Z_{2}$, $\cdots$, and subscripts $j_{1}$ and $j_{2}$ stand respectively for the row and column of the location of $\rho_{DCCA}(Z_{1},Z_{2},n)$.

The DFA-based bivariate regression model

We now translate the standard bivariate regression process described above into the DFA-based bivariate regression model. The two estimators in Eq. (4.2) can be extended to the scale-dependent estimators in the following way using the scale-dependent variance and covariance defined in Eqs. (4.7) and (4.9),

[TABLE]

Similarly, the scale-dependent residuals are

[TABLE]

with zero mean value. Inserting the calculated $\hat{e_{t}}(n)$ into the DFA process, we obtain the fluctuation $F_{\epsilon}^{2}(n)$ to estimate the variances of $\hat{\beta_{1}}^{DFA}(n)$ and $\hat{\beta_{2}}^{DFA}(n)$ via Eq. (4.12) as

[TABLE]

Then Eqs. (4.3)-(4.5) can be translated into the DFA regression form as

[TABLE]

[TABLE]

and

[TABLE]

Comparing to the standard $R^{2}$, $\beta^{*}$, and $\eta$, the scale-dependent $R^{2}_{DFA}(n)$, $\beta^{*DFA}(n)$, and $\eta^{DFA}(n)$ express more abundant information on model interpretation from multiple time scales.

Acknowledgements

This work was partially supported by National Natural Science Foundation of China (No.31501227, No. 11401577), and NSERC of Canada. The authors would like to thank two anonymous reviewers and the handling editor for their comments and suggestions, which led to a great improvement to the presentation of this work.

Author Contributions

F.W. designed the framework and performed the statistical analysis. L.W. wrote the manuscript and checked the analysis. Y.C. reviewed the analysis.

Additional Information

Competing Interests: The authors declare no competing interests.