# A Novel Approach for Nonstationary Time Series Analysis with Time-Invariant Correlation Coefficient

Authors: Chengrui Liu; Zhihua Wang; Huimin Fu; Yongbo Zhang
Journal: Mathematical Problems in Engineering
 (2014)

---

## Abstract

We will concentrate on the modeling and analysis of a class of nonstationary time series, called correlation coefficient stationary series, which commonly exists in practical engineering. First, the concept and scope of correlation coefficient stationary series are discussed to get a better understanding. Second, a theorem is proposed to determine standard deviation function for correlation coefficient stationary series. Third, we propose a moving multiple-point average method to determine the function forms for mean and standard deviation, which can help to improve the analysis precision, especially in the context of limited sample size. Fourth, the conditional likelihood approach is utilized to estimate the model parameters. In addition, we discuss the correlation coefficient stationarity test method, which can contribute to the verification of modeling validity. Monte Carlo simulation study illustrates the authentication of the theorem and the validity of the established method. Empirical study shows that the approach can satisfactorily explain the nonstationary behavior of many practical data sets, including stock returns, maximum power load, China money supply, and foreign currency exchange rate. The effectiveness of these processes is addressed by forecasting performance.


---

## Body

## 1. Introduction

Time series methods have been generally accepted as one of the most important means in an increasing number of real-world applications including finance. In the past several decades, considerable efforts have been made for time series analysis and prediction [1–3]. Time series approaches [4], regression models [5], artificial intelligence method [6], and Grey theory [7] are the commonly used techniques [8]. Many analyses are based on the assumption that the probabilistic properties of the underlying process are time invariant; that is, the series to be analyzed is covariance stationary. Modeling this stationary time series, one frequently chooses time series methods because of their high performance and robustness, which mainly include autoregressive (AR), moving average (MA), autoregressive moving average (ARMA), autoregressive integrated moving average (ARIMA), and Box-Jenkins models.Although the stationary assumption is very useful for the construction of simple models, it does not seem to be the best strategy in practice, and sometimes such stationarity assumptions are often questionable [9], because time series with time-varying means and variances are commonly seen in economic forecast [10], fault diagnosis [11], quality control [12], signal processing [13], performance test [14], automatic control [15], biopharmaceutical [16], and other fields. When the heteroscedasticity time series is processed by existing covariance stationary time series analysis method, the model parameters will lose the minimum variance property, and the variance estimator is no longer the unbiased estimation [17]. Referring to time series approaches and regression analysis, reasonable analysis and accurate prediction cannot be achieved for the nonstationary time series. Considering artificial intelligence, such as expert system and neural network, abundant prediction rule and practical experience from specific experts and large historical data banks are requisite for precise forecast. Although the Grey prediction model has been successfully applied in various fields and has demonstrated satisfactory results, its prediction performance still could be improved. The reason is that the Grey forecasting model is constructed of exponential function, and hence it may have worse prediction precise in the case of more random data sets.Meanwhile, great progress has been achieved related to process monitoring in industrial fields. To solve the multimode problem, illustrated by industrial process because of multiple production patterns in the same production line, various methods including partial least squares methods [18], model library-based methods [19], the Gaussian mixture model [20], the localized Fisher discriminant analysis approach [21], and the recent independent component analysis (ICA) based statistical processing methods [22, 23] have been constructed. The industrial process monitoring is of significant importance in the literature. However, the statistical method based on time series analysis is focused on in this study.In 1982, Engle [24] proposed the concept of conditional heteroscedasticity, with which they solved the conditional heteroscedasticity estimation problem for time series with constant unconditional variance. The proposed theory has been widely applied in financial risk evaluation. For his significant contribution, Engle gained 2003’s Nobel Economics Prize. However, the analysis problem for time series with unconditional time-varying variance still exists. It can be commonly seen in application [25, 26]. In addition, some hybrid models are also seen in the literature [27, 28], which combine dissimilar models or models that disagree with each other strongly to lower the generalization variance or error. Although hybrid models have shown advantages in some circumstances, there is no denying that they are much complicated for application.A simple nonstationary model contains a second-order stationary process modulated by a deterministic time-varying mean and a deterministic unconditional time-varying variance [29]. Let y
t, t
=
1,2
,
3
,
…, be a stationary process with zero mean and a simple nonstationary model can be given by x
t
=
μ
(
t
)
+
σ
(
t
)
y
t, where μ
(
t
) is the deterministic time-varying mean function and σ
(
t
) is the deterministic unconditional time-varying standard deviation function which is strictly positive. We can conclude that the non-stationarity of x
t is expressed by its evolving mean and unconditional variance. The efficient analysis of this nonstationary process is a substantial drawback in practice and has gradually attached great importance to the researchers.Based on systematic study of mass measured data, Fu and Liu [30] found that some common characteristics are shared by certain nonstationary time series. Like general nonstationary ones, these series exhibit time-varying mean μ
(
t
) and variance σ
2
(
t
), and their autocovariance function γ
(
t
,
t
+
τ
)
=
Cov
(
x
t
,
x
t
+
τ
)
=
E
{
[
x
t
-
μ
(
t
)
]
[
x
t
+
τ
-
μ
(
t
+
τ
)
]
} is no longer a univariate function of time interval τ, that is, γ
(
t
,
t
+
τ
)
≠
γ
(
τ
), while their correlation coefficient function ρ
(
t
,
t
+
τ
)
=
Cov
(
x
t
,
x
t
+
τ
)
/
[
σ
(
t
)
σ
(
t
+
τ
)
] is still a univariate function of time interval τ, that is, ρ
(
t
,
t
+
τ
)
=
ρ
τ. Accordingly, we can conclude that (i) they are not covariance stationary time series [31], whose autocovariance function γ
(
t
,
t
+
τ
) is a univariate function of time interval τ; that is, γ
(
t
,
t
+
τ
)
=
γ
(
τ
); (ii) they are a certain class of nonstationary time series and different from other nonstationary time series whose correlation coefficient function ρ
(
t
,
t
+
τ
) varies with time t. On this basis, Fu and Liu [30] proposed the concept of “correlation coefficient stationary process,” and discussed the establishment of the correlation coefficient autoregressive moving average (CCARMA) model.In this paper, we further study the nonstationary behavior of this correlation coefficient stationary series. First, characteristics of the variance function have been further studied, and a rigorous theorem was proposed, which can help not only determination of the standard deviation but also verification of the modeling process. Second, a rolling window determination scheme named moving multiple-point average method has been established to obtain the mean and standard deviation functions. This technology can enhance the accuracy under the same sample size, and the effect is more obvious in case of limited sample size. Third, we studied the scope of correlation coefficient stationary process, in which discussion can be helpful to better understand the concept of correlation coefficient stationary process. Finally, the correlation coefficient stationary test method has been investigated, which can assess the validity of the modeling process and make the modeling process a closed-loop system.In the next section, the concept of CCARMA process and its basic properties are introduced; we also discuss the CCARMA model. In Section3, we develop a method for determining the function forms for mean and standard deviation. Section 4 establishes the parameter determination method and a correlation coefficient stationary test method. Section 5 illustrates simulation studies to assess the validity of the approach. And Section 6 is devoted to the practical evaluation of the proposed method on several data sets, including daily returns to Shanghai composite index, Guangxi monthly maximum power load, China monthly money supply, and daily foreign exchange (FX) rate EUR/USD. A comparison between our forecasting results and ARIMA, variable differential, GARCH, GM
(
1,1
), and Modified GM
(
1,1
) models is also provided in this section. Finally, we conclude this paper with a discussion in Section 7.


## 2. Concept of CCARMA Process

### 2.1. Concept of CCARMA Process

Generally speaking, traditional stationarity means covariance stationarity [31]. Time series x
t, t
=
1,2
,
…
,
n, is a covariance stationary time series if the following two conditions are satisfied.(i)
The mean functionμ
(
t
)
=
E
(
x
t
) does not evolve through time; that is, μ
(
t
)
=
μ.
(ii)
The autocovariance functionγ
(
t
,
t
+
τ
)
=
Cov
(
x
t
,
x
t
+
τ
)
=
E
{
[
x
t
-
μ
(
t
)
]
[
x
t
+
τ
-
μ
(
t
+
τ
)
]
} is a univariate function of time interval τ; that is, γ
(
t
,
t
+
τ
)
=
γ
(
τ
).Letτ
=
0 and the autocovariance Cov
(
x
t
,
x
t
) equals the variance Var
⁡
(
x
t
) at time t. Consequently, variance of covariance stationarity time series dose not vary with time; that is, Var
⁡
(
x
t
)
=
σ
2. However, most time series, encountered in practice, cannot satisfy the above two requirements. To solve the analysis problem of a certain class of nonstationary time series, Fu and Liu [30] extended the above concept and proposed the following two concepts of correlation coefficient stationary time series.Concept 1. Let x
t, t
=
1,2
,
…
,
n, be a second-order moment time series, and let its correlation coefficient function ρ
(
t
,
t
+
τ
)
=
Cov
(
x
t
,
x
t
+
τ
)
/
[
σ
(
t
)
σ
(
t
+
τ
)
] be a univariate function of time interval τ; that is, ρ
(
t
,
t
+
τ
)
=
ρ
τ; then x
t, t
=
1,2
,
…
,
n, is called a correlation coefficient stationary time series.Concept 2. Correlation coefficient stationary time series x
t, t
=
1,2
,
…
,
n, is called CCARMA series, if y
t
=
[
x
t
-
μ
(
t
)
]
/
σ
(
t
), t
=
1,2
,
…
,
n, is an ARMA sequence, where E
(
x
t
)
=
μ
(
t
) and Var
⁡
(
x
t
)
=
σ
2
(
t
) are mean and variance functions of series x
t, t
=
1,2
,
…
,
n, respectively.The GaussianCCARMA
(
p
,
q
) model can be denoted by
(1)
Φ
(
L
)
x
t
-
μ
(
t
)
σ
(
t
)
=
Θ
(
L
)
ε
t
,
ε
t
~
NID
[
0
,
σ
ε
2
]
,
where Φ
(
L
)
=
1
-
φ
1
L
-
⋯
-
φ
p
L
p specifies the AR lag-polynomial, Θ
(
L
)
=
1
-
θ
1
L
-
⋯
-
θ
q
L
q specifies the MA polynomial, and ε
t
~
i
.
i
.
d
.
N
(
0
,
σ
ε
2
). As special cases of CCARMA
(
p
,
q
) model, the CCAR
(
p
) model and CCMA
(
q
) model can be given as
(2)
Φ
(
L
)
x
t
-
μ
(
t
)
σ
(
t
)
=
ε
t
,
ε
t
~
NID
[
0
,
σ
ε
2
]
,
x
t
-
μ
(
t
)
σ
(
t
)
=
Θ
(
L
)
ε
t
,
ε
t
~
NID
[
0
,
σ
ε
2
]
.Based on the above definitions, we can conclude the following.(i)
The statistical property difference between covariance stationarity series and correlation coefficient stationary series is whether its mean and variance vary with time. By plotting the sequence, this statistical property difference can be detected intuitively. And a quantitative method to determine the operation is the covariance stationary test method, introduced in Section4.2.
(ii)
The difference between the correlation coefficient stationary series and other nonstationary sequence lies in whether the correlation coefficientρ
(
t
,
t
+
τ
) of x
t and x
t
+
τ, τ
=
1,2
,
…, varies with time t, meaning that ρ
(
t
,
t
+
τ
) is only a univariate function of time interval τ. We can get a judgment by a simple way. First, divide the series into several subseries, which are considered as series with constant mean and variance. Second, calculate the correlation coefficients ρ
(
τ
), τ
=
1,2
,
…, of each subsequence. Generally speaking, it is enough to derive the first five-order correlation coefficients. Finally, examine whether the correlation coefficient ρ
(
τ
), τ
=
1,2
,
…, of each subsequence equals to each other. This can be completed by plotting or simple quantitative tests.


### 2.2. Scope of CCARMA Process

Basic properties show that the covariance stationary series is a special case of correlation coefficient stationary sequence. When the mean and variance do not vary with time; that is,E
(
x
t
)
=
μ, Var
⁡
(
x
t
)
=
σ
2, the correlation coefficient stationary series x
t, t
=
1,2
,
…
,
n, degenerates to covariance stationary sequence.Suppose thaty
t is a covariance stationary series and μ
(
t
) is a deterministic function; then we can conclude that series x
t
=
y
t
+
μ
(
t
), t
=
1,2
,
…
,
n, is a correlation coefficient stationary time series, with time-varying mean μ
(
t
) and constant variance Var
⁡
(
x
t
)
=
σ
2. For this series, y
t and μ
(
t
) represent the random part and deterministic part, respectively. Sequences of this kind are very common in practice, such as ground movement and deformation time series, meteorological data, and observation sequence in other fields. Wang [32] called it variance stationary sequence.Lety
t be a zero mean covariance stationary series, and σ
(
t
) is a deterministic positive function and then we can know that x
t
=
σ
(
t
)
×
y
t, t
=
1,2
,
…
,
n, is a correlation coefficient stationary time series with zero mean and time-varying variance. Amplitude modulation signal, commonly seen in radio communication, monitoring, and other fields, belongs to this case. In these signals, carrier signal is the zero mean covariance stationary series y
t, and modulated signal is the positive deterministic function σ
(
t
).Considering a more composite circumstance, we assume thaty
t is a zero mean covariance stationary series, μ
(
t
) is a deterministic function, and σ
(
t
) is a positive deterministic function; then, it can be inferred that x
t
=
σ
(
t
)
×
y
t
+
μ
(
t
), t
=
1,2
,
…
,
n, is a correlation coefficient stationary time series with time-varying mean and variance. Actually, this is a comprehensive result of the former two cases.Furthermore, correlation coefficient stationary series also includes the sequences which can satisfy the correlation coefficient stationary conditions. The correlation coefficient stability test method will be discussed in Section4.


## 2.1. Concept of CCARMA Process

Generally speaking, traditional stationarity means covariance stationarity [31]. Time series x
t, t
=
1,2
,
…
,
n, is a covariance stationary time series if the following two conditions are satisfied.(i)
The mean functionμ
(
t
)
=
E
(
x
t
) does not evolve through time; that is, μ
(
t
)
=
μ.
(ii)
The autocovariance functionγ
(
t
,
t
+
τ
)
=
Cov
(
x
t
,
x
t
+
τ
)
=
E
{
[
x
t
-
μ
(
t
)
]
[
x
t
+
τ
-
μ
(
t
+
τ
)
]
} is a univariate function of time interval τ; that is, γ
(
t
,
t
+
τ
)
=
γ
(
τ
).Letτ
=
0 and the autocovariance Cov
(
x
t
,
x
t
) equals the variance Var
⁡
(
x
t
) at time t. Consequently, variance of covariance stationarity time series dose not vary with time; that is, Var
⁡
(
x
t
)
=
σ
2. However, most time series, encountered in practice, cannot satisfy the above two requirements. To solve the analysis problem of a certain class of nonstationary time series, Fu and Liu [30] extended the above concept and proposed the following two concepts of correlation coefficient stationary time series.Concept 1. Let x
t, t
=
1,2
,
…
,
n, be a second-order moment time series, and let its correlation coefficient function ρ
(
t
,
t
+
τ
)
=
Cov
(
x
t
,
x
t
+
τ
)
/
[
σ
(
t
)
σ
(
t
+
τ
)
] be a univariate function of time interval τ; that is, ρ
(
t
,
t
+
τ
)
=
ρ
τ; then x
t, t
=
1,2
,
…
,
n, is called a correlation coefficient stationary time series.Concept 2. Correlation coefficient stationary time series x
t, t
=
1,2
,
…
,
n, is called CCARMA series, if y
t
=
[
x
t
-
μ
(
t
)
]
/
σ
(
t
), t
=
1,2
,
…
,
n, is an ARMA sequence, where E
(
x
t
)
=
μ
(
t
) and Var
⁡
(
x
t
)
=
σ
2
(
t
) are mean and variance functions of series x
t, t
=
1,2
,
…
,
n, respectively.The GaussianCCARMA
(
p
,
q
) model can be denoted by
(1)
Φ
(
L
)
x
t
-
μ
(
t
)
σ
(
t
)
=
Θ
(
L
)
ε
t
,
ε
t
~
NID
[
0
,
σ
ε
2
]
,
where Φ
(
L
)
=
1
-
φ
1
L
-
⋯
-
φ
p
L
p specifies the AR lag-polynomial, Θ
(
L
)
=
1
-
θ
1
L
-
⋯
-
θ
q
L
q specifies the MA polynomial, and ε
t
~
i
.
i
.
d
.
N
(
0
,
σ
ε
2
). As special cases of CCARMA
(
p
,
q
) model, the CCAR
(
p
) model and CCMA
(
q
) model can be given as
(2)
Φ
(
L
)
x
t
-
μ
(
t
)
σ
(
t
)
=
ε
t
,
ε
t
~
NID
[
0
,
σ
ε
2
]
,
x
t
-
μ
(
t
)
σ
(
t
)
=
Θ
(
L
)
ε
t
,
ε
t
~
NID
[
0
,
σ
ε
2
]
.Based on the above definitions, we can conclude the following.(i)
The statistical property difference between covariance stationarity series and correlation coefficient stationary series is whether its mean and variance vary with time. By plotting the sequence, this statistical property difference can be detected intuitively. And a quantitative method to determine the operation is the covariance stationary test method, introduced in Section4.2.
(ii)
The difference between the correlation coefficient stationary series and other nonstationary sequence lies in whether the correlation coefficientρ
(
t
,
t
+
τ
) of x
t and x
t
+
τ, τ
=
1,2
,
…, varies with time t, meaning that ρ
(
t
,
t
+
τ
) is only a univariate function of time interval τ. We can get a judgment by a simple way. First, divide the series into several subseries, which are considered as series with constant mean and variance. Second, calculate the correlation coefficients ρ
(
τ
), τ
=
1,2
,
…, of each subsequence. Generally speaking, it is enough to derive the first five-order correlation coefficients. Finally, examine whether the correlation coefficient ρ
(
τ
), τ
=
1,2
,
…, of each subsequence equals to each other. This can be completed by plotting or simple quantitative tests.


## 2.2. Scope of CCARMA Process

Basic properties show that the covariance stationary series is a special case of correlation coefficient stationary sequence. When the mean and variance do not vary with time; that is,E
(
x
t
)
=
μ, Var
⁡
(
x
t
)
=
σ
2, the correlation coefficient stationary series x
t, t
=
1,2
,
…
,
n, degenerates to covariance stationary sequence.Suppose thaty
t is a covariance stationary series and μ
(
t
) is a deterministic function; then we can conclude that series x
t
=
y
t
+
μ
(
t
), t
=
1,2
,
…
,
n, is a correlation coefficient stationary time series, with time-varying mean μ
(
t
) and constant variance Var
⁡
(
x
t
)
=
σ
2. For this series, y
t and μ
(
t
) represent the random part and deterministic part, respectively. Sequences of this kind are very common in practice, such as ground movement and deformation time series, meteorological data, and observation sequence in other fields. Wang [32] called it variance stationary sequence.Lety
t be a zero mean covariance stationary series, and σ
(
t
) is a deterministic positive function and then we can know that x
t
=
σ
(
t
)
×
y
t, t
=
1,2
,
…
,
n, is a correlation coefficient stationary time series with zero mean and time-varying variance. Amplitude modulation signal, commonly seen in radio communication, monitoring, and other fields, belongs to this case. In these signals, carrier signal is the zero mean covariance stationary series y
t, and modulated signal is the positive deterministic function σ
(
t
).Considering a more composite circumstance, we assume thaty
t is a zero mean covariance stationary series, μ
(
t
) is a deterministic function, and σ
(
t
) is a positive deterministic function; then, it can be inferred that x
t
=
σ
(
t
)
×
y
t
+
μ
(
t
), t
=
1,2
,
…
,
n, is a correlation coefficient stationary time series with time-varying mean and variance. Actually, this is a comprehensive result of the former two cases.Furthermore, correlation coefficient stationary series also includes the sequences which can satisfy the correlation coefficient stationary conditions. The correlation coefficient stability test method will be discussed in Section4.


## 3. Function Form Determination of Mean and Standard Deviation

We know that one can hardly efficiently obtain the mean and standard deviation functions when these two functions discontinuously vary with time. Consequently, in this study, we consider the general cases in which mean and standard deviation functions vary with time continuously and slowly. This is a common assumption in model constructing for time series analysis. In our theoretical study process, we found that some constraints have to be satisfied for rigorous derivation. Accordingly, we proposed the following theorem for determining the mean and standard deviation functions for nonstationary time series.Theorem 1.
Letx
t, t
=
1,2
,
…
,
n, be a Gaussian correlation coefficient stationary series with time-varying mean E
(
x
t
)
=
μ
(
t
) and variance Var
⁡
(
x
t
)
=
σ
2
(
t
) and its standard deviation function σ
(
t
) has the same form with the trend item of series |
∇
x
t
|
=
|
x
t
+
1
-
μ
(
t
+
1
)
-
x
t
+
μ
(
t
)
|, t
=
1,2
,
…
,
n
-
1, when Δ
σ
(
t
)
/
σ
(
t
) is a constant or [
Δ
σ
(
t
)
/
σ
(
t
)
]
max
⁡ is a negligible small amount compared with one; that is, [
Δ
σ
(
t
)
/
σ
(
t
)
]
max
⁡
≪
1, where [
Δ
σ
(
t
)
/
σ
(
t
)
]
max
⁡ is the maximum of Δ
σ
(
t
)
/
σ
(
t
) and Δ
σ
(
t
)
=
σ
(
t
+
1
)
-
σ
(
t
), t
=
1,2
,
…
,
n; then, the standard deviation function σ
(
t
) can be depicted by
(3)
σ
(
t
)
=
c
E
|
x
t
+
1
-
μ
(
t
+
1
)
-
x
t
+
μ
(
t
)
|
,
where c is a positive real number. See Appendix A for theorem proof.Generally speaking,[
Δ
σ
(
t
)
/
σ
(
t
)
]
max
⁡ can be considered as a negligible small amount when two orders smaller than one. It can be inferred, from the above theorem, that the standard deviation σ
(
t
) has the same function form with the mean function of series |
x
t
+
1
-
μ
(
t
+
1
)
-
x
t
+
μ
(
t
)
|, t
=
1,2
,
…
,
n
-
1. Consequently, in the process of mean and standard deviation function form determination, we need to conduct the following steps: (1) obtain the trend estimator μ
⌢
(
t
) and derive series |
x
t
+
1
-
μ
⌢
(
t
+
1
)
-
x
t
+
μ
⌢
(
t
)
|, t
=
1,2
,
…
,
n
-
1; (2) determine the trend of series |
x
t
+
1
-
μ
⌢
(
t
+
1
)
-
x
t
+
μ
⌢
(
t
)
|, t
=
1,2
,
…
,
n
-
1, which we take as the function form of standard deviation function σ
(
t
); (3) take the result σ
⌢
(
t
) as the estimate of standard deviation function when the theorem condition can be satisfied.Otherwise, when the theorem conditions cannot be met; that is, neitherΔ
σ
(
t
)
/
σ
(
t
) is a constant nor [
Δ
σ
(
t
)
/
σ
(
t
)
]
max
⁡ is a negligible small amount compared with one, we have to change the determination strategy. In this case, its standard deviation function σ
(
t
) has the same function form with the trend item of series |
x
t
-
μ
(
t
)
|, t
=
1,2
,
…
,
n. That is, for correlation coefficient stationary series x
t, t
=
1,2
,
…
,
n, its standard deviation function can be determined by
(4)
σ
(
t
)
=
c
1
E
|
x
t
-
μ
(
t
)
|
,
where c
1 is a positive real number. See Appendices for proof.The function form determination of mean and variance focuses on accessing the trend items ofx
t, |
x
t
+
1
-
μ
⌢
(
t
+
1
)
-
x
t
+
μ
⌢
(
t
)
| or |
x
t
-
μ
⌢
(
t
)
|, t
=
1,2
,
…
,
n. We consider that the trend function contains nonperiodic part and periodic part. In this paper, we propose a rolling window method called “moving multiple-point average method” for the determination of sequence trend item. In order to determine the nonperiodic part in the trend item, the proposed method movingly fits on the whole sample data length n with the multiple-point average method. Meanwhile we adopt the sample periodogram method to obtain the periodic part. To better address this issue, the following steps can be performed.(1) Determine the periodic part of trend item with sample periodogram method.First, we suppose the existence of frequencyω
1
,
ω
2
,
…
,
ω
M, and then we can express the periodic part of series x
t, t
=
1,2
,
…
,
n, with periodogram method [33] as
(5)
x
t
=
x
-
+
∑
j
=
1
M
{
α
j
cos
⁡
[
ω
j
(
t
-
1
)
]
+
β
j
sin
⁡
[
ω
j
(
t
-
1
)
]
}
,
t
=
1,2
,
…
,
n
,
where x
- is constant mean of series x
t, t
=
1,2
,
…
,
n, which can be obtained by x
-
=
∑
t
=
1
n
x
t
/
n; M is the existing frequency number, which equals n
/
2 when sample size n is an even number or equals (
n
-
1
)
/
2 when n is an odd number; ω
j
=
2
j
π
/
n, j
=
1,2
,
…
,
M, is an existing frequency; α
j and β
j are cosine and sine coefficients corresponding to frequency ω
j, j
=
1,2
,
…
,
M.When sample sizen is an odd number, coefficients α
j and β
j, j
=
1,2
,
…
,
M, can be calculated by
(6)
α
j
=
2
n
∑
t
=
1
n
x
t
cos
⁡
⁡
[
ω
j
(
t
-
1
)
]
,
j
=
1,2
,
…
,
M
,
β
j
=
2
n
∑
t
=
1
n
x
t
sin
[
ω
j
(
t
-
1
)
]
,
j
=
1,2
,
…
,
M
.
When sample size n is an even number, coefficients α
j and β
j, j
=
1,2
,
…
,
M
-
1, can also be worked out by (6) and
(7)
α
M
=
1
n
∑
t
=
1
n
(
-
1
)
t
-
1
x
t
,
β
M
=
0
.Then, we introduce a parameterA
j
=
α
j
2
+
β
j
2 depicting the amplitude of frequency ω
j, j
=
1,2
,
…
,
M. When one or more A
j is  significantly greater than the other ones A
1
,
A
2
,
…
,
A
j
-
1
,
A
j
+
1
,
…
,
A
M, we can affirm that a periodic item with frequency ω
j exists. And then the existing periodic item with frequency ω
j can be expressed as α
j
cos
⁡
[
ω
j
(
t
-
1
)
]
+
β
j
sin
[
ω
j
(
t
-
1
)
]. For the circumstance of multiple frequencies, the periodic item is sum of the periodic items corresponding to each crest value
(8)
∑
j
{
α
j
cos
⁡
[
ω
j
(
t
-
1
)
]
+
β
j
sin
[
ω
j
(
t
-
1
)
]
}
.(2) Calculate nonperiodic part of trend item with moving multiple-point average method.With the former results of step (1), nonperiodic trend part can be obtained by(9)
x
t
-
∑
j
{
α
j
cos
⁡
[
ω
j
(
t
-
1
)
]
+
β
j
sin
⁡
[
ω
j
(
t
-
1
)
]
}
,
t
=
1,2
,
…
,
n
.
Then we select the point number of each averaging segment m (subsequence length) and moving time interval Δ based on the volatility of the obtained series from (9). Generally speaking, m is in the range of [
n
/
50
,
n
/
2
], where n is the sample length and Δ is in the range of [
m
/
10
,
m
/
2
]. In order to get an accurate periodic part of trend item, we should note that averaging point number m must be not less than ent
(
2
π
/
ω
min
⁡
), whereas ω
min
⁡ is the smallest frequency in the determined periodic function; that is, (8).Then, we implement movingm-point average on the whole sample data length n based on the moving time interval Δ and obtain a group of mean value (
t
i
¯
,
x
i
¯
) by
(10)
t
i
¯
=
1
m
∑
t
=
Δ
i
-
1
Δ
i
+
m
-
2
t
,
i
=
1,2
,
…
,
l
,
(11)
x
i
¯
=
1
m
∑
t
=
Δ
i
-
1
Δ
i
+
m
-
2
x
t
,
i
=
1,2
,
…
,
l
,
where l
=
ent
(
(
n
-
m
)
/
Δ
)
+
1 indicates the total averaging times.Consequently, fit the obtained group of mean value(
t
i
¯
,
x
i
¯
), i
=
1,2
,
…
,
l, to get the regression function f
(
t
), which is the nonperiodic trend function part of series x
t, t
=
1,2
,
…
,
n.(3) Redetermine the periodic part of trend item.Based on the nonperiodic trend functionf
(
t
) obtained above, the following series can be calculated:
(12)
x
t
-
f
(
t
)
,
t
=
1,2
,
…
,
n
.
Process the obtained series from (12) and then the periodic part function expressed by (8) can be redetermined with the periodogram method.(4) Repeat step (2) and step (3) until each parameter results in periodic part function∑
j
{
α
j
cos
⁡
[
ω
j
(
t
-
1
)
]
+
β
j
sin
[
ω
j
(
t
-
1
)
]
} and nonperiodic part function f
(
t
) becomes numerical stabilized. Then we can obtain the final trend item expression as
(13)
μ
⌢
(
t
)
=
f
(
t
)
+
∑
j
{
α
j
cos
⁡
[
ω
j
(
t
-
1
)
]
+
β
j
sin
⁡
[
ω
j
(
t
-
1
)
]
}
.The mean functionμ
(
t
) can be directly determined by implementing the above steps from (1) to (4) on series x
t, t
=
1,2
,
…
,
n. Based on the theorem given in Section 3, the standard deviation function σ
(
t
) can be obtained by conducting the same steps from (1) to (4) on series |
x
t
+
1
-
μ
(
t
+
1
)
-
x
t
+
μ
(
t
)
|
=
|
x
t
+
1
-
μ
⌢
(
t
+
1
)
-
x
t
+
μ
⌢
(
t
)
|, t
=
1,2
,
…
,
n
-
1, if the theorem condition can be satisfied. Otherwise, the standard deviation function σ
(
t
) can be obtained by conducting the same steps from (1) to (4) on series on sequence |
x
t
-
μ
(
t
)
|
=
|
x
t
-
μ
⌢
(
t
)
|, t
=
1,2
,
…
,
n.In a word, the function forms of meanμ
(
t
) and standard deviation σ
(
t
) can be determined. In order to facilitate the parameter estimation process, we depict the time-varying functions of mean and standard deviation by
(14)
μ
(
a
,
t
)
=
∑
i
=
0
r
a
i
ϕ
i
(
t
)
,
σ
(
b
,
t
)
=
∑
j
=
0
s
b
j
ψ
j
(
t
)
,
where ϕ
0
(
t
)
=
ψ
0
(
t
)
=
1, ϕ
i
(
t
), i
=
1,2
,
…
,
r, and ψ
j
(
t
),j
=
1,2
,
…
,
s, are functions that can be known through the above determination process and a
=
(
a
0
,
a
1
,
…
,
a
r
)
T and b
=
(
b
0
,
b
1
,
…
,
b
s
)
T are general sets of unknown parameters to be calculated.


## 4. Model Construction and Testing

### 4.1. CCARMA Model Parameter Estimation

Letx
t be a CCARMA series with mean μ
(
a
,
t
) and standard deviation σ
(
b
,
t
) and then transformed sequence y
t
=
[
x
t
-
μ
(
t
)
]
/
σ
(
t
) is an ARMA series according to concept 2 in Section 2.1. The relationship between x
t and y
t can be rewritten as x
t
=
μ
(
t
)
+
σ
(
t
)
y
t. Joint probability density functions (PDF) f
X
n
,
X
n
-
1
,
…
,
X
1
(
x
n
,
x
n
-
1
,
…
,
x
1
) can be derived by PDF f
Y
n
,
Y
n
-
1
,
…
,
Y
1
(
y
n
,
y
n
-
1
,
…
,
y
1
) through
(15)
f
X
n
,
X
n
-
1
,
…
,
X
1
(
x
n
,
x
n
-
1
,
…
,
x
1
)
=
f
Y
n
,
Y
n
-
1
,
⋯
,
Y
1
(
y
n
,
y
n
-
1
,
…
,
y
1
)
σ
(
t
n
)
σ
(
t
n
-
1
)
⋯
σ
(
t
1
)
.
See Appendices for proof.According to time series analysis theory [33], a common approximation of the likelihood function for ARMA process conditions on initial values of both y’s and ε’s. Based on the recommendation given by Box and Jenkins [34], we set ε’s to zero for k
=
0
,
-
1
,
…
,
-
q
+
1, and y’s to their actual values for k
=
0
,
-
1
,
…
,
-
p
+
1. Then the sequence {
ε
1
,
ε
2
,
…
,
ε
n
} can be calculated from {
y
1
,
y
2
,
…
,
y
n
}, and the conditional log likelihood is then
(16)
ln
⁡
L
=
-
n
-
λ
2
ln
⁡
(
2
π
)
-
(
n
-
λ
)
ln
⁡
σ
ε
-
1
2
σ
ε
2
∑
k
=
λ
+
1
n
ε
k
2
.Considering the GaussianCCARMA
(
p
,
q
) process x
t depicted by (1), suppose that we have a sample of n observations x
t, t
=
1,2
,
…
,
n. Maximum likelihood estimation with conditional likelihood function is utilized to estimate the vector of population parameters θ
=
(
a
0
,
a
1
,
…
,
a
r
,
b
0
,
b
1
,
…
,
b
s
,
φ
1
,
φ
2
,
…
,
φ
p
,
θ
1
,
θ
2
,
…
,
θ
q
,
σ
ε
)
T. Accordingly, a common approximation of the likelihood function for CCARMA process conditions on initial values of both x’s and ε’s. We set ε’s to zero for k
=
0
,
-
1
,
…
,
-
q
+
1, and x’s to their actual values for k
=
0
,
-
1
,
…
,
-
p
+
1. Then the sequence {
ε
1
,
ε
2
,
⋯
,
ε
n
} can be calculated from {
x
1
,
x
2
,
…
,
x
n
} by iterating on
(17)
ε
k
=
x
k
-
μ
(
a
,
t
k
)
σ
(
b
,
t
k
)
-
∑
i
=
1
p
φ
i
x
k
-
i
-
μ
(
a
,
t
k
-
i
)
σ
(
b
,
t
k
-
i
)
+
∑
i
=
1
q
θ
i
ε
k
-
i
for k
=
1,2
,
…
,
n. Based on the joint PDF relationship depicted by (15), the conditional log likelihood is then
(18)
ln
⁡
L
(
θ
)
=
-
n
-
λ
2
ln
⁡
(
2
π
)
-
(
n
-
λ
)
ln
⁡
σ
ε
-
1
2
σ
ε
2
∑
k
=
λ
+
1
n
ε
k
2
-
∑
k
=
1
n
ln
⁡
σ
(
b
,
t
k
)
,
where λ equals the maximum one of p and q; that is, λ
=
max
⁡
(
p
,
q
). Model parameters can be determined by solving the following equations:
(19)
∂
ln
⁡
L
(
θ
)
∂
a
l
=
-
1
σ
ε
2
∑
k
=
λ
+
1
n
ε
k
∂
ε
k
∂
a
l
=
0
,
l
=
0,1
,
…
,
r
,
(20)
∂
ln
⁡
L
(
θ
)
∂
b
l
=
-
1
σ
ε
2
∑
k
=
λ
+
1
n
ε
k
∂
ε
k
∂
b
l
-
∑
k
=
1
n
1
σ
(
b
,
t
k
)
∂
σ
(
b
,
t
k
)
∂
b
l
=
0
,
l
=
0,1
,
…
,
s
,
(21)
∂
ln
⁡
L
(
θ
)
∂
φ
l
=
-
1
σ
ε
2
∑
k
=
λ
+
1
n
ε
k
∂
ε
k
∂
φ
l
=
0
,
l
=
0,1
,
…
,
p
,
(22)
∂
ln
⁡
L
(
θ
)
∂
θ
l
=
-
1
σ
ε
2
∑
k
=
λ
+
1
n
ε
k
∂
ε
k
∂
θ
l
=
0
,
l
=
0,1
,
…
,
q
,
(23)
∂
ln
⁡
L
(
θ
)
∂
σ
ε
=
-
n
-
λ
σ
ε
+
1
σ
ε
3
∑
k
=
λ
+
1
n
ε
k
2
=
0
.Based on the initial values forε’s given above, we can conclude that ∂
ε
k
/
∂
a
l, ∂
ε
k
/
∂
b
l, ∂
ε
k
/
∂
φ
l, and ∂
ε
k
/
∂
θ
l equal zero for k
=
0
,
-
1
,
…
,
-
q
+
1. Then these values for k
=
1,2
,
…
,
n can be calculated by iterating on
(24)
∂
ε
k
∂
a
l
=
-
1
I
(
b
,
t
k
)
∂
μ
(
a
,
t
k
)
∂
a
l
+
∑
i
=
1
p
φ
i
I
(
b
,
t
k
-
i
)
∂
μ
(
a
,
t
k
-
i
)
∂
a
l
+
∑
i
=
1
q
θ
i
∂
ε
k
-
i
∂
a
l
=
0
,
∂
ε
k
∂
b
l
=
-
x
k
-
∂
μ
(
a
,
t
k
)
I
2
(
b
,
t
k
)
I
(
b
,
t
k
)
∂
a
l
+
∑
i
=
1
p
η
i
[
x
k
-
i
-
μ
(
a
,
t
k
-
i
)
]
I
2
(
b
,
t
k
-
i
)
∂
I
(
b
,
t
k
-
i
)
∂
b
l
+
∑
i
=
1
q
θ
i
∂
ε
k
-
i
∂
b
l
=
0
,
∂
ε
k
∂
η
l
=
-
x
k
-
l
-
∂
μ
(
a
,
t
k
-
l
)
I
(
b
,
t
k
-
l
)
+
∑
i
=
1
q
θ
i
∂
ε
k
-
i
∂
η
l
=
0
,
∂
ε
k
∂
θ
l
=
ε
k
-
l
+
∑
i
=
1
q
θ
i
∂
ε
k
-
i
∂
θ
l
=
0
.
Then (23) can be rewritten as
(25)
σ
ε
2
=
1
n
-
λ
∑
i
=
λ
+
1
n
ε
i
2
=
0
.Parametersa
0
,
a
1
,
…
,
a
r, b
0
,
b
1
,
…
,
b
s, φ
1
,
φ
2
,
…
,
φ
p, and θ
1
,
θ
2
,
…
,
θ
q can be obtained by solving (19) to (22), and parameter σ
ε can be calculated through (25).


### 4.2. CCARMA Model Testing

Here, we focus on the model testing method, which can demonstrate whether the proposed model is appropriate to describe the observation data. For the nonstationary seriesx
t, t
=
1,2
,
…
,
n, based on the parameter estimation process given above, we can obtain its mean function μ
(
t
) and standard deviation function σ
(
t
). Then we can implement the correlation coefficient stationary test of series x
t through the covariance stationary test of series y
t, where y
t
=
[
x
t
-
μ
(
t
)
]
/
σ
(
t
). Several tests have been proposed and applied to examine the covariance stationarity in literature. In this paper, we take the postsample prediction testing method presented by Pagan and Schwert [35] as an example to illustrate the testing procedure. Postsample prediction test for covariance stationarity is a nonparametric method, and it is facilitative to implement and familiar by scholars.First, obtain the sequencey
t, t
=
1,2
,
…
,
n, through the transformation y
t
=
[
x
t
-
μ
(
t
)
]
/
σ
(
t
), based on the determined mean function μ
(
t
) and standard deviation function σ
(
t
).Second, split seriesy
t, t
=
1,2
,
…
,
n, averagely into two parts, and calculate the sample variance σ
⌢
(
1
)
2
=
E
(
(
2
/
n
)
∑
i
=
1
n
/
2
y
i
2
) and σ
⌢
(
2
)
2
=
E
(
(
2
/
n
)
∑
i
=
1
+
n
/
2
n
y
i
2
) for each. Then the test statistic τ
⌢
=
σ
⌢
(
2
)
2
-
σ
⌢
(
1
)
2 follows
(26)
2
n
τ
⌢
~
N
[
0,2
(
R
0
+
2
∑
i
=
1
∞
R
i
)
]
.
If y
t
2 is a covariance stationary process with autocovariances R
i, and let ν
=
R
0
+
2
∑
i
=
1
∞
R
i; then it can be estimated by
(27)
ν
⌢
=
R
⌢
0
+
2
∑
i
=
1
8
R
⌢
i
(
1
-
i
9
)
,
where R
⌢
i is the estimated serial correlation coefficients of y
t
2 calculated over the whole sample.Finally, define null hypothesis thatH
0: y
t, t
=
1,2
,
…
,
n, is a covariance stationarity series versus alternative hypothesis that H
a: y
t, t
=
1,2
,
…
,
n, is not a covariance stationarity series. Construct test statistic B, and the rejection region can be expressed by
(28)
B
=
|
n
2
τ
⌢
2
ν
|
>
u
1
-
α
/
2
,
where u
p is the 100
uth percentile of the standard normal distribution and α is the selected significance level indicating the probability of type I error.


## 4.1. CCARMA Model Parameter Estimation

Letx
t be a CCARMA series with mean μ
(
a
,
t
) and standard deviation σ
(
b
,
t
) and then transformed sequence y
t
=
[
x
t
-
μ
(
t
)
]
/
σ
(
t
) is an ARMA series according to concept 2 in Section 2.1. The relationship between x
t and y
t can be rewritten as x
t
=
μ
(
t
)
+
σ
(
t
)
y
t. Joint probability density functions (PDF) f
X
n
,
X
n
-
1
,
…
,
X
1
(
x
n
,
x
n
-
1
,
…
,
x
1
) can be derived by PDF f
Y
n
,
Y
n
-
1
,
…
,
Y
1
(
y
n
,
y
n
-
1
,
…
,
y
1
) through
(15)
f
X
n
,
X
n
-
1
,
…
,
X
1
(
x
n
,
x
n
-
1
,
…
,
x
1
)
=
f
Y
n
,
Y
n
-
1
,
⋯
,
Y
1
(
y
n
,
y
n
-
1
,
…
,
y
1
)
σ
(
t
n
)
σ
(
t
n
-
1
)
⋯
σ
(
t
1
)
.
See Appendices for proof.According to time series analysis theory [33], a common approximation of the likelihood function for ARMA process conditions on initial values of both y’s and ε’s. Based on the recommendation given by Box and Jenkins [34], we set ε’s to zero for k
=
0
,
-
1
,
…
,
-
q
+
1, and y’s to their actual values for k
=
0
,
-
1
,
…
,
-
p
+
1. Then the sequence {
ε
1
,
ε
2
,
…
,
ε
n
} can be calculated from {
y
1
,
y
2
,
…
,
y
n
}, and the conditional log likelihood is then
(16)
ln
⁡
L
=
-
n
-
λ
2
ln
⁡
(
2
π
)
-
(
n
-
λ
)
ln
⁡
σ
ε
-
1
2
σ
ε
2
∑
k
=
λ
+
1
n
ε
k
2
.Considering the GaussianCCARMA
(
p
,
q
) process x
t depicted by (1), suppose that we have a sample of n observations x
t, t
=
1,2
,
…
,
n. Maximum likelihood estimation with conditional likelihood function is utilized to estimate the vector of population parameters θ
=
(
a
0
,
a
1
,
…
,
a
r
,
b
0
,
b
1
,
…
,
b
s
,
φ
1
,
φ
2
,
…
,
φ
p
,
θ
1
,
θ
2
,
…
,
θ
q
,
σ
ε
)
T. Accordingly, a common approximation of the likelihood function for CCARMA process conditions on initial values of both x’s and ε’s. We set ε’s to zero for k
=
0
,
-
1
,
…
,
-
q
+
1, and x’s to their actual values for k
=
0
,
-
1
,
…
,
-
p
+
1. Then the sequence {
ε
1
,
ε
2
,
⋯
,
ε
n
} can be calculated from {
x
1
,
x
2
,
…
,
x
n
} by iterating on
(17)
ε
k
=
x
k
-
μ
(
a
,
t
k
)
σ
(
b
,
t
k
)
-
∑
i
=
1
p
φ
i
x
k
-
i
-
μ
(
a
,
t
k
-
i
)
σ
(
b
,
t
k
-
i
)
+
∑
i
=
1
q
θ
i
ε
k
-
i
for k
=
1,2
,
…
,
n. Based on the joint PDF relationship depicted by (15), the conditional log likelihood is then
(18)
ln
⁡
L
(
θ
)
=
-
n
-
λ
2
ln
⁡
(
2
π
)
-
(
n
-
λ
)
ln
⁡
σ
ε
-
1
2
σ
ε
2
∑
k
=
λ
+
1
n
ε
k
2
-
∑
k
=
1
n
ln
⁡
σ
(
b
,
t
k
)
,
where λ equals the maximum one of p and q; that is, λ
=
max
⁡
(
p
,
q
). Model parameters can be determined by solving the following equations:
(19)
∂
ln
⁡
L
(
θ
)
∂
a
l
=
-
1
σ
ε
2
∑
k
=
λ
+
1
n
ε
k
∂
ε
k
∂
a
l
=
0
,
l
=
0,1
,
…
,
r
,
(20)
∂
ln
⁡
L
(
θ
)
∂
b
l
=
-
1
σ
ε
2
∑
k
=
λ
+
1
n
ε
k
∂
ε
k
∂
b
l
-
∑
k
=
1
n
1
σ
(
b
,
t
k
)
∂
σ
(
b
,
t
k
)
∂
b
l
=
0
,
l
=
0,1
,
…
,
s
,
(21)
∂
ln
⁡
L
(
θ
)
∂
φ
l
=
-
1
σ
ε
2
∑
k
=
λ
+
1
n
ε
k
∂
ε
k
∂
φ
l
=
0
,
l
=
0,1
,
…
,
p
,
(22)
∂
ln
⁡
L
(
θ
)
∂
θ
l
=
-
1
σ
ε
2
∑
k
=
λ
+
1
n
ε
k
∂
ε
k
∂
θ
l
=
0
,
l
=
0,1
,
…
,
q
,
(23)
∂
ln
⁡
L
(
θ
)
∂
σ
ε
=
-
n
-
λ
σ
ε
+
1
σ
ε
3
∑
k
=
λ
+
1
n
ε
k
2
=
0
.Based on the initial values forε’s given above, we can conclude that ∂
ε
k
/
∂
a
l, ∂
ε
k
/
∂
b
l, ∂
ε
k
/
∂
φ
l, and ∂
ε
k
/
∂
θ
l equal zero for k
=
0
,
-
1
,
…
,
-
q
+
1. Then these values for k
=
1,2
,
…
,
n can be calculated by iterating on
(24)
∂
ε
k
∂
a
l
=
-
1
I
(
b
,
t
k
)
∂
μ
(
a
,
t
k
)
∂
a
l
+
∑
i
=
1
p
φ
i
I
(
b
,
t
k
-
i
)
∂
μ
(
a
,
t
k
-
i
)
∂
a
l
+
∑
i
=
1
q
θ
i
∂
ε
k
-
i
∂
a
l
=
0
,
∂
ε
k
∂
b
l
=
-
x
k
-
∂
μ
(
a
,
t
k
)
I
2
(
b
,
t
k
)
I
(
b
,
t
k
)
∂
a
l
+
∑
i
=
1
p
η
i
[
x
k
-
i
-
μ
(
a
,
t
k
-
i
)
]
I
2
(
b
,
t
k
-
i
)
∂
I
(
b
,
t
k
-
i
)
∂
b
l
+
∑
i
=
1
q
θ
i
∂
ε
k
-
i
∂
b
l
=
0
,
∂
ε
k
∂
η
l
=
-
x
k
-
l
-
∂
μ
(
a
,
t
k
-
l
)
I
(
b
,
t
k
-
l
)
+
∑
i
=
1
q
θ
i
∂
ε
k
-
i
∂
η
l
=
0
,
∂
ε
k
∂
θ
l
=
ε
k
-
l
+
∑
i
=
1
q
θ
i
∂
ε
k
-
i
∂
θ
l
=
0
.
Then (23) can be rewritten as
(25)
σ
ε
2
=
1
n
-
λ
∑
i
=
λ
+
1
n
ε
i
2
=
0
.Parametersa
0
,
a
1
,
…
,
a
r, b
0
,
b
1
,
…
,
b
s, φ
1
,
φ
2
,
…
,
φ
p, and θ
1
,
θ
2
,
…
,
θ
q can be obtained by solving (19) to (22), and parameter σ
ε can be calculated through (25).


## 4.2. CCARMA Model Testing

Here, we focus on the model testing method, which can demonstrate whether the proposed model is appropriate to describe the observation data. For the nonstationary seriesx
t, t
=
1,2
,
…
,
n, based on the parameter estimation process given above, we can obtain its mean function μ
(
t
) and standard deviation function σ
(
t
). Then we can implement the correlation coefficient stationary test of series x
t through the covariance stationary test of series y
t, where y
t
=
[
x
t
-
μ
(
t
)
]
/
σ
(
t
). Several tests have been proposed and applied to examine the covariance stationarity in literature. In this paper, we take the postsample prediction testing method presented by Pagan and Schwert [35] as an example to illustrate the testing procedure. Postsample prediction test for covariance stationarity is a nonparametric method, and it is facilitative to implement and familiar by scholars.First, obtain the sequencey
t, t
=
1,2
,
…
,
n, through the transformation y
t
=
[
x
t
-
μ
(
t
)
]
/
σ
(
t
), based on the determined mean function μ
(
t
) and standard deviation function σ
(
t
).Second, split seriesy
t, t
=
1,2
,
…
,
n, averagely into two parts, and calculate the sample variance σ
⌢
(
1
)
2
=
E
(
(
2
/
n
)
∑
i
=
1
n
/
2
y
i
2
) and σ
⌢
(
2
)
2
=
E
(
(
2
/
n
)
∑
i
=
1
+
n
/
2
n
y
i
2
) for each. Then the test statistic τ
⌢
=
σ
⌢
(
2
)
2
-
σ
⌢
(
1
)
2 follows
(26)
2
n
τ
⌢
~
N
[
0,2
(
R
0
+
2
∑
i
=
1
∞
R
i
)
]
.
If y
t
2 is a covariance stationary process with autocovariances R
i, and let ν
=
R
0
+
2
∑
i
=
1
∞
R
i; then it can be estimated by
(27)
ν
⌢
=
R
⌢
0
+
2
∑
i
=
1
8
R
⌢
i
(
1
-
i
9
)
,
where R
⌢
i is the estimated serial correlation coefficients of y
t
2 calculated over the whole sample.Finally, define null hypothesis thatH
0: y
t, t
=
1,2
,
…
,
n, is a covariance stationarity series versus alternative hypothesis that H
a: y
t, t
=
1,2
,
…
,
n, is not a covariance stationarity series. Construct test statistic B, and the rejection region can be expressed by
(28)
B
=
|
n
2
τ
⌢
2
ν
|
>
u
1
-
α
/
2
,
where u
p is the 100
uth percentile of the standard normal distribution and α is the selected significance level indicating the probability of type I error.


## 5. Simulation Experiment

To assess the computational performance of the proposed method and determine whether the approach seems to give reasonable results, we study the effectiveness and the performance of presented methods in Sections3 and 4 from a Monte Carlo simulated example.The following zero meanCCAR
(
1
) model
(29)
x
t
σ
(
t
)
=
φ
x
t
-
1
σ
(
t
-
1
)
+
ε
t
,
ε
t
~
NID
[
0
,
σ
ε
2
]
is considered with four simulation function forms for the standard deviation, including linear function, quadratic function, periodic function, and the combination function of periodic function and linear function. The parameters assumed in each experiment are summarized in Table 1 for convenience. The sample size is n
=
200. The point number of averaging and moving time interval for each experiment is also listed in Table 1.Table 1
Model specifications by experiments.
Experiment
[
Δ
σ
(
t
)
/
σ
(
t
)
]
max
⁡
m
Δ
φ
c
1
c
2
c
3
(1)σ
(
t
)
=
c
1
+
c
2
t
0.0040
40
20
0.4
5
0.02
∖
0.6
5
0.02
∖
0.8
5
0.02
∖
(2)σ
(
t
)
=
c
1
+
c
2
t
+
c
3
t
3
0.0087
40
20
0.4
4
−0.01
0.0003
0.6
4
−0.01
0.0003
0.8
4
−0.01
0.0003
(3)σ
(
t
)
=
c
1
+
c
2
sin
⁡
⁡
(
0.03
π
t
)
0.0241
∖
∖
0.4
13
2
∖
0.6
13
2
∖
0.8
13
2
∖
(4)σ
(
t
)
=
c
1
+
c
2
t
+
c
3
sin
⁡
⁡
(
0.03
π
t
)
0.0154
∖
∖
0.4
13
0.02
2
0.6
13
0.02
2
0.8
13
0.02
2
Note:m and Δ are point number of each averaging and moving time interval; φ is the autocorrelation coefficient of simulation model; c
1, c
2, and c
3 are model parameters in standard deviation function; [
Δ
σ
(
t
)
/
σ
(
t
)
]
max
⁡ is the maximum of Δ
σ
(
t
)
/
σ
(
t
). Symbol “∖” indicates that the parameter does not exist.Series|
x
t
+
1
-
x
t
| and |
x
t
| can be taken to determine the standard deviation function σ
(
t
) by (3) (Theorem method) presented in Section 3, since the mean of series x
t equals zero; that is, μ
(
t
)
=
0. Simulation results are analyzed by the index of average percent relative error err
⁡
=
(
1
/
n
)
∑
t
=
1
n
|
σ
(
t
)
-
σ
⌢
(
t
)
|
×
100
/
σ
(
t
), which is presented in Table 2. The results are based on 20000 Monte Carlo simulations with innovations drawn from an IID Gaussian distribution.Table 2
Results of parameter and the average percent relative error for each experiment.
φ
0.4
0.6
0.8
0.4
0.6
0.8
Method
Equation(3)
Equation(4)
Equation(3)
Equation(4)
Equation(3)
Equation(4)
Equation(3)
Equation(4)
Equation(3)
Equation(4)
Equation(3)
Equation(4)
Experiment 1
Experiment 2
c
1
4.7635
4.3543
4.4568
4.9790
4.1984
6.6076
4.7635
4.3543
4.4568
4.9790
4.1984
6.6076
c
2
0.0190
0.0174
0.0178
0.0200
0.0168
0.0268
0.0190
0.0174
0.0178
0.0200
0.0168
0.0268
Per. err./%
7.71
13.48
11.57
8.85
16.15
33.71
7.71
13.48
11.57
8.85
16.15
33.71
Experiment 3
Experiment 4
c
1
12.3941
11.3212
11.6007
12.9740
10.9307
17.2514
12.4157
11.3302
11.5889
12.9288
10.9200
17.1642
c
2
1.9023
1.7321
1.7246
1.9175
1.6054
2.5234
0.0188
0.0173
0.0178
0.0203
0.0169
0.0276
c
3
∖
∖
∖
∖
∖
∖
1.8880
1.7320
1.7573
1.9839
1.6753
2.6349
Per. err./%
10.22
14.49
11.21
8.32
15.94
33.30
8.62
13.82
11.89
10.29
16.17
34.58
Note: Experiments 1 to 4 are defined in Table1. Symbol “∖” indicates that the parameter does not exist.Nonzero meanCCAR
(
p
) and CCARMA
(
p
,
q
) model with higher order and different parameter specifications are also considered to study the authentication of the theorem, and the results were very similar to those reported in this paper.


## 6. Empirical Results

In this section, we focus on the practical performance of the proposed approach. Experiments are presented for four different economic data sets presented in Section6.1. In Section 6.2, the mean and standard deviation functions are determined. Then, in Section 6.3, we apply the statistical test method discussed in Section 4. Finally, Section 6.4 is devoted to evaluating the forecasting performance of the correlation coefficient stationary method.


### 6.1. The Data Sets

Daily CIR. The daily returns to composite index for Shanghai from January 5, 1999, through September 30, 2003 (1131 observations).Monthly MPL. The monthly maximum power load for Guangxi from January 1990 through December 1999 (120 observations).Monthly M2. The monthly money supply for China from January 2000 through December 2009 (112 observations).Daily FX Rate EUR/USD. Euro to the United States dollar parity from January 1, 2005, through December 30, 2005 (260 observations).


### 6.2. Determination of Mean and Variance

In order to test the stability of the three data sets, we apply the determination method presented in Sections3 and 4 on the four data sets. All parameter results are summarized in Tables 3 and 4. Note that the standard deviation functions σ
(
t
) are obtained by (3) for Daily CIR, Monthly MPL, and Daily FX rate data sets and by (4) for Monthly M2, because the theorem condition can be satisfied for the former three sequences.Table 3
Determination results of mean and variance functions.
Trend functionμ
(
t
)
Daily CIR
0
Monthly MPL
1057.62
+
20.498
t
-
10.297
sin
⁡
⁡
(
π
t
6
)
+
15.926
cos
⁡
(
π
t
6
)
Monthly M2
119862.92
+
224.06
t
+
65.280
t
2
-
0.7726
t
3
+
4.077
×
1
0
-
3
t
4
Daily FX rate
1.2489
+
5.1943
×
1
0
-
3
t
-
1.1693
×
1
0
-
4
t
2
+
9.4416
×
1
0
-
7
t
3
-
3.3508
×
1
0
-
9
t
4
+
4.3846
×
1
0
-
12
t
5
+
0.01006
cos
⁡
(
3
π
t
65
)
+
0.01587
cos
⁡
(
3
π
t
130
)
Standard deviation functionσ
(
t
)
Daily CIR
0.5456
+
3.1940
×
1
0
-
2
t
-
1.9532
×
1
0
-
4
t
2
+
4.3798
×
1
0
-
7
t
3
-
4.1174
×
1
0
-
10
t
4
+
1.3751
×
1
0
-
13
t
5
Monthly MPL
116.07
-
7.9024
sin
⁡
⁡
(
π
t
6
)
-
4.0838
cos
⁡
(
π
t
6
)
Monthly M2
2579.50
-
350.24
t
+
17.622
t
2
-
0.2823
t
3
+
1.4410
×
1
0
-
3
t
4
Daily FX rate
5.4449
×
1
0
-
3Table 4
Model parameter results.
Data set
Modely
t
=
∑
i
=
1
p
φ
i
y
t
-
i
+
ε
t
p
φ
1
φ
2
Daily CIR
1
0.02861
∖
Monthly MPL
1
0.6353
∖
Monthly M2
2
1.1572
−0.1401
Daily FX rate
1
0.8616
∖


### 6.3. Testing for Correlation Coefficient Stationarity

As a second step, we apply the correlation coefficient stationarity test presented in Section4 on the four data sets. All results are summarized in Table 5. We first derive sequence y
t by transformation y
t
=
[
x
t
-
μ
(
t
)
]
/
σ
(
t
) with the results listed in Table 3. Let significance level α
=
0.05 and the rejection region is B
>
u
1
-
α
/
2
=
1.96.Table 5
Testing results of correlation coefficient stationarity.
Data set
DailyCIR
MonthlyMPL
MonthlyM2
Daily FXrate
Testing statisticB
0.6420
0.3340
0.9236
0.9741From the results, we can conclude that the four data sets are correlation coefficient stationary time series, and the application of the proposed method is reasonable.


### 6.4. Prediction

Based on the CCAR models constructed above, the forecast of original sequencex
t can be obtained by (1). In this section, we will consider autoregressive integrated moving average (ARIMA) model, variable differential (VD) model, generalized autoregressive conditional heteroscedasticity (GARCH) model, Grey prediction model GM
(
1,1
), and modified GM
(
1,1
) model which are considered as standard remarks to study its forecasting accuracy. The future 5 daily Shanghai composite index data, the future 24 monthly maximum power load data for Guangxi from January 2000 through December 2001, the future 8 monthly money supply data for China from May 2009 through December 2009, and the future 10 daily FX rate data EUR/USD are forecasted by the established models. The prediction results are summarized from Tables 6, 7, 8, and 9, and they are analyzed by the index of percent relative error (Per. err.) ε
i
=
|
x
n
+
i
-
x
⌢
n
+
i
|
×
100
/
x
n
+
i, i
=
1,2
,
…
,
l, and mean percent relative error (Mper. err) ε
-
=
(
1
/
l
)
∑
i
=
1
l
ε
i, where l is the prediction step ahead.Table 6
Prediction results for future five daily Shanghai composite indexes.
Time
Real value
CCAR(1) model
AR model
GARCH model
Prediction value
Per. err./%
Prediction value
Per. err./%
Prediction value
Per. err./%
Oct-8th-03
1371.69
1367.50
0.30
1366.17
0.40
1366.27
0.40
Oct-9th-03
1369.17
1371.69
0.18
1365.83
0.24
1366.93
0.16
Oct-10th-03
1404.01
1369.17
2.48
1364.84
2.78
1365.68
2.73
Oct-13th-03
1399.66
1404.01
0.31
1364.86
2.48
1365.73
2.42
Oct-14th-03
1388.17
1399.66
0.83
1364.97
1.67
1365.80
1.61
Mper. Err/%
0.82
1.51
1.46Table 7
Prediction results of the Guangxi monthly maximum power load (%).
Time
Per. err.
Time
Per. err.
Time
Per. err.
CCAR(1) model
VD model
CCAR(1) model
VD model
CCAR(1) model
VD model
Jan-00
3.34
1.48
Sep-00
9.43
5.10
May-00
1.34
0.29
Feb-00
1.13
1.25
Oct-00
6.32
2.35
Jun-00
2.85
3.78
Mar-00
2.11
5.50
Nov-00
3.20
1.14
Jul-00
3.18
3.18
Apr-00
1.88
5.78
Dec-00
4.40
1.12
Aug-01
1.47
2.01
May-00
1.64
3.66
Jan-01
1.73
1.51
Sep-01
2.68
2.02
Jun-00
6.11
1.26
Feb-01
3.86
5.77
Oct-01
1.16
0.96
Jul-00
7.60
4.65
Mar-01
0.25
0.76
Nov-01
5.02
5.51
Aug-00
7.78
4.53
Apr-01
1.41
0.95
Dec-01
0.59
1.09
Mper. Err of CCAR(1) Model
Mper. Err of VD Mode
3.35
2.73Table 8
Prediction results of China money supply.
Time
Real value/10−1
Billion RMB
CCAR(2) model
ARIMA(6,2,0) model
Prediction value/10−1
Billion RMB
Per. err./%
Prediction value/10−1
Billion RMB
Per. err./%
May-09
548263.51
550638.76
0.43
548033.39
0.04
Jun-09
568916.20
561106.15
1.37
565101.99
0.67
Jul-09
573102.85
571893.07
0.21
579773.54
1.16
Aug-09
576698.95
583009.86
1.09
595508.25
3.26
Sep-09
585405.34
594467.16
1.55
611409.40
4.44
Oct-09
586643.29
606275.88
3.35
625911.63
6.69
Nov-09
594604.72
618447.18
4.01
636980.09
7.13
Dec-09
610224.52
630992.50
3.40
651301.43
6.73
Mper. Err/%
1.9274
3.7665Table 9
Prediction results of Daily FX rate.
Time
Real value
CCAR(1) model
GM(1,1)
MGM(1,1)
Prediction value
per. err./%
Prediction value
per. err./%
Prediction value
per. err./%
Dec-19th-05
1.2007
1.1986
0.1751
1.1701
2.5474
1.1517
4.0793
Dec-20th-05
1.1864
1.2006
1.1993
1.1695
1.4220
1.1507
3.0054
Dec-21th-05
1.1830
1.2026
1.6565
1.1689
1.1881
1.1498
2.8053
Dec-22th-05
1.1879
1.2045
1.3948
1.1684
1.6448
1.1489
3.2809
Dec-23th-05
1.1865
1.2062
1.6615
1.1678
1.5779
1.1481
3.2377
Dec-26th-05
1.1856
1.2078
1.8730
1.1672
1.5524
1.1473
3.2312
Dec-27th-05
1.1832
1.2092
2.1993
1.1666
1.4020
1.1465
3.0982
Dec-28th-05
1.1830
1.2104
2.3197
1.1660
1.4346
1.1458
3.1410
Dec-29th-05
1.1843
1.2115
2.2925
1.1654
1.5919
1.1452
3.3025
Dec-30th-05
1.1832
1.2122
2.4541
1.1649
1.5496
1.1446
3.2638
Mper. Err/%
1.7226
1.5911
3.2245From the results, we can conclude that the application of the proposed method for correlation coefficient stationary time series is reasonable and effective. Furthermore, the presented methodology can be considered good and shows a promise for future applications in nonstationary time series analysis and forecasting.


## 6.1. The Data Sets

Daily CIR. The daily returns to composite index for Shanghai from January 5, 1999, through September 30, 2003 (1131 observations).Monthly MPL. The monthly maximum power load for Guangxi from January 1990 through December 1999 (120 observations).Monthly M2. The monthly money supply for China from January 2000 through December 2009 (112 observations).Daily FX Rate EUR/USD. Euro to the United States dollar parity from January 1, 2005, through December 30, 2005 (260 observations).


## 6.2. Determination of Mean and Variance

In order to test the stability of the three data sets, we apply the determination method presented in Sections3 and 4 on the four data sets. All parameter results are summarized in Tables 3 and 4. Note that the standard deviation functions σ
(
t
) are obtained by (3) for Daily CIR, Monthly MPL, and Daily FX rate data sets and by (4) for Monthly M2, because the theorem condition can be satisfied for the former three sequences.Table 3
Determination results of mean and variance functions.
Trend functionμ
(
t
)
Daily CIR
0
Monthly MPL
1057.62
+
20.498
t
-
10.297
sin
⁡
⁡
(
π
t
6
)
+
15.926
cos
⁡
(
π
t
6
)
Monthly M2
119862.92
+
224.06
t
+
65.280
t
2
-
0.7726
t
3
+
4.077
×
1
0
-
3
t
4
Daily FX rate
1.2489
+
5.1943
×
1
0
-
3
t
-
1.1693
×
1
0
-
4
t
2
+
9.4416
×
1
0
-
7
t
3
-
3.3508
×
1
0
-
9
t
4
+
4.3846
×
1
0
-
12
t
5
+
0.01006
cos
⁡
(
3
π
t
65
)
+
0.01587
cos
⁡
(
3
π
t
130
)
Standard deviation functionσ
(
t
)
Daily CIR
0.5456
+
3.1940
×
1
0
-
2
t
-
1.9532
×
1
0
-
4
t
2
+
4.3798
×
1
0
-
7
t
3
-
4.1174
×
1
0
-
10
t
4
+
1.3751
×
1
0
-
13
t
5
Monthly MPL
116.07
-
7.9024
sin
⁡
⁡
(
π
t
6
)
-
4.0838
cos
⁡
(
π
t
6
)
Monthly M2
2579.50
-
350.24
t
+
17.622
t
2
-
0.2823
t
3
+
1.4410
×
1
0
-
3
t
4
Daily FX rate
5.4449
×
1
0
-
3Table 4
Model parameter results.
Data set
Modely
t
=
∑
i
=
1
p
φ
i
y
t
-
i
+
ε
t
p
φ
1
φ
2
Daily CIR
1
0.02861
∖
Monthly MPL
1
0.6353
∖
Monthly M2
2
1.1572
−0.1401
Daily FX rate
1
0.8616
∖


## 6.3. Testing for Correlation Coefficient Stationarity

As a second step, we apply the correlation coefficient stationarity test presented in Section4 on the four data sets. All results are summarized in Table 5. We first derive sequence y
t by transformation y
t
=
[
x
t
-
μ
(
t
)
]
/
σ
(
t
) with the results listed in Table 3. Let significance level α
=
0.05 and the rejection region is B
>
u
1
-
α
/
2
=
1.96.Table 5
Testing results of correlation coefficient stationarity.
Data set
DailyCIR
MonthlyMPL
MonthlyM2
Daily FXrate
Testing statisticB
0.6420
0.3340
0.9236
0.9741From the results, we can conclude that the four data sets are correlation coefficient stationary time series, and the application of the proposed method is reasonable.


## 6.4. Prediction

Based on the CCAR models constructed above, the forecast of original sequencex
t can be obtained by (1). In this section, we will consider autoregressive integrated moving average (ARIMA) model, variable differential (VD) model, generalized autoregressive conditional heteroscedasticity (GARCH) model, Grey prediction model GM
(
1,1
), and modified GM
(
1,1
) model which are considered as standard remarks to study its forecasting accuracy. The future 5 daily Shanghai composite index data, the future 24 monthly maximum power load data for Guangxi from January 2000 through December 2001, the future 8 monthly money supply data for China from May 2009 through December 2009, and the future 10 daily FX rate data EUR/USD are forecasted by the established models. The prediction results are summarized from Tables 6, 7, 8, and 9, and they are analyzed by the index of percent relative error (Per. err.) ε
i
=
|
x
n
+
i
-
x
⌢
n
+
i
|
×
100
/
x
n
+
i, i
=
1,2
,
…
,
l, and mean percent relative error (Mper. err) ε
-
=
(
1
/
l
)
∑
i
=
1
l
ε
i, where l is the prediction step ahead.Table 6
Prediction results for future five daily Shanghai composite indexes.
Time
Real value
CCAR(1) model
AR model
GARCH model
Prediction value
Per. err./%
Prediction value
Per. err./%
Prediction value
Per. err./%
Oct-8th-03
1371.69
1367.50
0.30
1366.17
0.40
1366.27
0.40
Oct-9th-03
1369.17
1371.69
0.18
1365.83
0.24
1366.93
0.16
Oct-10th-03
1404.01
1369.17
2.48
1364.84
2.78
1365.68
2.73
Oct-13th-03
1399.66
1404.01
0.31
1364.86
2.48
1365.73
2.42
Oct-14th-03
1388.17
1399.66
0.83
1364.97
1.67
1365.80
1.61
Mper. Err/%
0.82
1.51
1.46Table 7
Prediction results of the Guangxi monthly maximum power load (%).
Time
Per. err.
Time
Per. err.
Time
Per. err.
CCAR(1) model
VD model
CCAR(1) model
VD model
CCAR(1) model
VD model
Jan-00
3.34
1.48
Sep-00
9.43
5.10
May-00
1.34
0.29
Feb-00
1.13
1.25
Oct-00
6.32
2.35
Jun-00
2.85
3.78
Mar-00
2.11
5.50
Nov-00
3.20
1.14
Jul-00
3.18
3.18
Apr-00
1.88
5.78
Dec-00
4.40
1.12
Aug-01
1.47
2.01
May-00
1.64
3.66
Jan-01
1.73
1.51
Sep-01
2.68
2.02
Jun-00
6.11
1.26
Feb-01
3.86
5.77
Oct-01
1.16
0.96
Jul-00
7.60
4.65
Mar-01
0.25
0.76
Nov-01
5.02
5.51
Aug-00
7.78
4.53
Apr-01
1.41
0.95
Dec-01
0.59
1.09
Mper. Err of CCAR(1) Model
Mper. Err of VD Mode
3.35
2.73Table 8
Prediction results of China money supply.
Time
Real value/10−1
Billion RMB
CCAR(2) model
ARIMA(6,2,0) model
Prediction value/10−1
Billion RMB
Per. err./%
Prediction value/10−1
Billion RMB
Per. err./%
May-09
548263.51
550638.76
0.43
548033.39
0.04
Jun-09
568916.20
561106.15
1.37
565101.99
0.67
Jul-09
573102.85
571893.07
0.21
579773.54
1.16
Aug-09
576698.95
583009.86
1.09
595508.25
3.26
Sep-09
585405.34
594467.16
1.55
611409.40
4.44
Oct-09
586643.29
606275.88
3.35
625911.63
6.69
Nov-09
594604.72
618447.18
4.01
636980.09
7.13
Dec-09
610224.52
630992.50
3.40
651301.43
6.73
Mper. Err/%
1.9274
3.7665Table 9
Prediction results of Daily FX rate.
Time
Real value
CCAR(1) model
GM(1,1)
MGM(1,1)
Prediction value
per. err./%
Prediction value
per. err./%
Prediction value
per. err./%
Dec-19th-05
1.2007
1.1986
0.1751
1.1701
2.5474
1.1517
4.0793
Dec-20th-05
1.1864
1.2006
1.1993
1.1695
1.4220
1.1507
3.0054
Dec-21th-05
1.1830
1.2026
1.6565
1.1689
1.1881
1.1498
2.8053
Dec-22th-05
1.1879
1.2045
1.3948
1.1684
1.6448
1.1489
3.2809
Dec-23th-05
1.1865
1.2062
1.6615
1.1678
1.5779
1.1481
3.2377
Dec-26th-05
1.1856
1.2078
1.8730
1.1672
1.5524
1.1473
3.2312
Dec-27th-05
1.1832
1.2092
2.1993
1.1666
1.4020
1.1465
3.0982
Dec-28th-05
1.1830
1.2104
2.3197
1.1660
1.4346
1.1458
3.1410
Dec-29th-05
1.1843
1.2115
2.2925
1.1654
1.5919
1.1452
3.3025
Dec-30th-05
1.1832
1.2122
2.4541
1.1649
1.5496
1.1446
3.2638
Mper. Err/%
1.7226
1.5911
3.2245From the results, we can conclude that the application of the proposed method for correlation coefficient stationary time series is reasonable and effective. Furthermore, the presented methodology can be considered good and shows a promise for future applications in nonstationary time series analysis and forecasting.


## 7. Discussion

In this paper, we discussed the category of correlation coefficient stationary series, a nonstationary time series with time-varying mean and variance. We proposed a moving determination method for its time-varying mean function and standard deviation function. We also discussed the correlation coefficient stationary test method.The determination principle of function form and order of meanμ
(
t
) and standard deviation σ
(
t
) cannot be separated from primary sequence analysis. It is worth noting in prediction problem that the function models of μ
(
t
) and σ
(
t
) should also trade off between sequence volatility and accuracy requirements.For mean function determination, the proposed moving method can be used to establish the mean function of all kinds of nonstationary time series and can help improve the accuracy, especially in situations of small sample size. For standard deviation function fixing, the presented theorem is only suitable for the correlation coefficient stationary sequence.In addition, the mean and standard deviation functions, depicted by (14) are general expressions. Take the Monthly MPL data set in Section 6 for example; we can assume that ϕ
0
(
t
)
=
ψ
0
(
t
)
=
1, ϕ
1
(
t
)
=
t, ϕ
2
(
t
)
=
ψ
1
(
t
)
=
sin
(
π
t
/
6
), ϕ
3
(
t
)
=
ψ
2
(
t
)
=
cos
⁡
(
π
t
/
6
), and a
=
(
a
0
,
a
1
,
a
2
,
a
3
)
′ and b
=
(
b
0
,
b
1
,
b
2
)
′.



---

*Source: 148432-2014-01-22.xml*