A Second Order Cumulant Spectrum Test That a Stochastic Process is Strictly Stationary and a Step Toward a Test for Graph Signal Strict Stationarity
AA Second Order Cumulant Spectrum BasedTest for Strict Stationarity
Douglas Patterson, Department of Finance, Virginia Tech,Melvin Hinich † , Denisa Roberts AmazonJanuary 23, 2018
Abstract
This article develops a statistical test for the null hypothesis of strict stationarityof a discrete time stochastic process. When the null hypothesis is true, the secondorder cumulant spectrum is zero at all the discrete Fourier frequency pairs presentin the principal domain of the cumulant spectrum. The test uses a frame (window)averaged sample estimate of the second order cumulant spectrum to build a teststatistic that has an asymptotic complex standard normal distribution. We derivethe test statistic, study the size and power properties of the test, and demonstrate itsimplementation with intraday stock market return data. The test has conservativesize properties and good power to detect varying variance and unit root in the presenceof varying variance.
Keywords:
Discrete Time, Second Order Cumulant Spectrum, Stochastic Process, TimeVarying Variance, Unit Root, High Frequency Stock Returns Douglas Patterson is Professor of Finance in the Department of Finance, Pamplin College of Business,Virginia Tech, Blacksburg, VA 24061 (email: [email protected]). Deceased September 6, 2010. Denisa Roberts is Machine Learning Scientist at Amazon, NY, NY (email: [email protected]).This work was completed before joining Amazon. a r X i v : . [ q -f i n . S T ] J a n . Introduction. A stochastic process is called strictly stationary if the joint distribution of x ( t + τ ), x ( t + τ ), . . . , x ( t n + τ ) is independent of τ for all time points { t , t , . . . , t n } . In otherwords, the time series has the same joint distribution across all time points t . The hypothe-sis of stationarity plays a significant role in time series analysis and is frequently used as anassumption in both academic research and practical applications. A unit root is one type ofviolation of stationarity in the strict sense. Kwiatkowski et al. (1992) developed a test forunit root (or level stationarity) commonly mentioned in the stationarity related literature,hereafter referred to as KPSS. Because of its popularity, we use KPSS as a starting pointof comparison for our work.Weak stationarity, commonly used as an assumption in financial time series forecastingmodels, assumes that the first two moments of the time series are constant through time(under suitable mixing conditions that the series is summable, and the first two momentsexist Brillinger (2001)). Hence a process is weakly stationary (or stationary in the widesense), if its expected value E { x ( t m ) } is constant and its autocovariance c x ( t m ) dependsonly on the difference ( t n + m − t n ), c x ( t m ) = c x ( t n + m − t n ) (Brillinger (2001)). A weaklystationary Gaussian process is also stationary in the strict sense Brillinger (2001). AGaussian process is a sequence of independent, identically normally distributed randomvariables { u t } .There are two paradigms for analyzing time series: the time domain and the frequencydomain. In the time domain, one considers the observed data { x ( t ) } directly and typicallymakes conjectures about its moments. In the frequency domain, one decomposes the timeseries into underlying frequencies and makes conjectures about spectra and cumulant spec-tra. The frequency domain offers its own advantages to decompose time series into trends2nd work with higher order moments more easily. In this article we develop a statisticaltest of the hypothesis that an observed stochastic process { x ( t ) } is stationary in the strictsense. For convenience, we will refer to our test as PHR. PHR uses sample estimates ofthe second order cumulant spectrum to test the null hypothesis of strict stationarity. The PHR second order cumulant spectrum based test of strict stationarity is different fromwidely used tests for structural change of a linear time series model. For example, Farleyand Hinich (1970) proposed a control chart to detect small shifts in the mean of a stochastictime series modeled as a white noise process plus a mean parameter that can be constantor shifting. Their proposed chart draws no conclusions on the stationarity of the process.Estimation and inference about the change in the mean of a stationary random processwas amply examined in literature. Many articles studies the estimation and inference aboutthe change in the mean of a random process. Hawkins (1977) proposed a likelihood ratiobased procedure to detect a change in the mean while assuming constant variance. Ina different approach, Dickinson et al. (2014) developed CUSUM and EWMA charts thatdetect changes in probability distributions via their scale parameter. Also, tests of a unitroot with or without drift against an alternative hypothesis of trend stationarity, can befound, for example, in Dickey and Fuller (1979) and Phillips and Perron (1988).The space of nonstrictly stationary processes is vast and, therefore, in this manuscriptwe focus the evaluation of the power of the PHR test on a few alternative hypotheses. Wedo not study PHR’s power against the alternative of a process with periodicities in themean and random variation in the periodic structure, as in Hinich and Wild (2001). In thespace of tests for strict stationarity in the time domain, Lima and Neri (2013) proposed atest for strict stationarity and compared its power to the KPSS test of level stationarity.3imilarly to Lima and Neri (2013) we compare the power of the PHR test to the power ofthe KPSS test. We calculate the power for the alternative hypotheses of unit root, varyingvariance, and their combination. Our PHR test excels in detecting stationarity violationsdue to second moment.We present the properties of the second order cumulant spectrum in Section 2, the testdevelopment in Section 3, empirical size and power evaluation via Monte Carlo simulationsand comparison with the KPSS test under various alternative hypotheses in Section 4. Thetest is demonstrated in 5 using high frequency (trade by trade) stock returns for the Fordcompany .
2. The Second Order Cumulant Spectrum of a Strictly Stationary Process.
In the frequency domain, the second order cumulant spectrum of a stochastic process { x ( t n ) } with t n = nτK ( f , f ) = ∞ (cid:88) n = −∞ ∞ (cid:88) n = −∞ E [ x ( t n ) x ( t n )] exp [ − π ( f t n + f t n )] , (1)where f i are discrete Fourier frequencies.The cumulant spectrum K is defined for frequency values in the square {− f ≤ f ≤ f , − f ≤ f ≤ f } , where f = 1 / (2 τ ) is called the Nyquist frequency. Fourier frequenciesdisplay the following properties and symmetries: • Frequencies are between zero and one, normalized by Nyquist frequency. The Nyquistfrequency f is equal to the inverse of twice the sampling rate of a discrete timeprocess. The time series must be sampled at specific time intervals. We assumefor simplicity that the sampling interval length denoted τ is one. When τ = 1 theNyquist frequency is f = 0 .
5. 4
The discrete Fourier transform of a sample of length T at frequency f , X ( f ) = (cid:80) T − t =0 x ( t ) exp {− i πf t } is equal to the conjugate of the Fourier transform at − f , X ( f ) = X ∗ ( − f ). • The Fourier transform is periodic at intervals of the sampling frequency for discretetime processes. Higher frequencies fold down into the (0 , f ) interval.Because of these symmetries, we focus the study of the second order cumulant spectrumin the principal domain only. The principal domain contains a subset of frequencies thatare not aliased or indistinguishable from each other (Bloomfield (2004)), also called funda-mental frequencies. The principal domain forms the triangle { < f ≤ f , − f < f ≤ f } .In general, to cross from the time domain into the frequency domain, we use the discreteFourier transform of the time series realization at a series of frequencies in the principaldomain. From Hinich (1994), if the time series { x ( t ) } is strictly stationary then thesecond order cumulant spectrum for a pair of frequencies in the principal domain ( f , f ), K ( f , f )is equal to K ( f , f ) = E [ X ( f ) X ( f )] = S ( f ) + O (1) if f + f = 0 mod 1 O (1) otherwise, (2)where X ( f i ) denote the discrete Fourier transform of the time series at the frequenciesincluded in the principal domain. The line f = − f is not in the principal domain (Hinich(1994) ). From Equation 2, K ( f , − f ) = δ ( f ) S ( f ), where δ ( f ) is the Dirac deltafunction, and S ( f ) = (cid:80) ∞ m = −∞ c x ( t m ) exp ( − i πf t m ) is the spectrum at frequency f and c x ( t m ) is the autocovariance of the process.In summary, if a time series is strictly stationary, then its second order cumulant spec-trum 5 ( f , f ) = 0 (3)for all frequency pairs ( f , f ) in the principal domain Hinich (1994).
3. Derivation of The Test Statistic for the Second Order Cumulant Spectrum Testof Strict Stationarity.
We now use the second order cumulant spectrum to construct the test for the null hypoth-esis that a stochastic process { x ( t n ) } is strictly stationary. The test uses frame averagedcumulant spectra at fundamental frequency pairs and has an asymptotic complex normaldistribution. To get frame averaged estimates for the cumulant spectrum K ( f , f ) and forthe spectrum S ( f ), where ( f , f ) are frequencies in the principal domain, we partition thetime series into frames of equal length, calculate the cumulant spectrum and the spectrumfor each frame and then average over all the frames. A time series of length T is partitionedinto P = [ T /L ] complete nonoverlapping frames of length L . We omit the last frame if ithas less than L observations. Then the frame averaged cumulant spectrum estimate isˆ K ( f , f ) = 1 P P (cid:88) p =1 K p ( f , f ) , (4)We can also express the principal domain for K in terms of integers ( k , k ), with f = k/L , and 0 < k ≤ L/ − k ≤ k ≤ k . The steps to get the frame averaged estimatesof the cumulant spectrum, ˆ K , and of the spectrum, ˆ S , are calculated following the steps: • Use a reasonable length time series, ideally of no less than T = 150 observations. • Split the time series into P = [ T /L ] complete consecutive nonoverlapping frames orwindows of L elements, where L is sufficiently large. Let the p frame of length L be6 x p (1) , x p (2) , ..., x p ( L ) } . • Calculate the discrete Fourier transform (DFT) for the p frame of length L as X p ( k ) = L − (cid:88) t =0 x [ t + ( pL )] exp {− i πf k [ t + ( pL )] } , (5)where f k = kL . • Compute the spectrum at each frequency in the principal domain of K for each frame p S p ( k ) = 1 L | X p ( k ) | = 1 L X p ( k ) X p ( − k ) (6) • Compute the second order cumulant spectrum for each frame pK p ( k , k ) = 1 L X p ( k ) X p ( k ) (7)for pairs of frequencies defined by 0 < k ≤ L/ − k < k ≤ k . • Get the frame averaged ˆ S and ˆ K over the P frames for each k and ( k , k ) pairrespectively.From the theory in Hinich (1994) and in Brillinger (2001), if L and P are sufficientlylarge, the expected value of the frame averaged spectrum is equal to its theoretical valueup to O ( L − ) , E (cid:104) ˆ S ( f k ) (cid:105) = S ( f k ) + O (cid:0) L − (cid:1) . (8)7imilarly to Equation 8, under strict stationarity, from Hinich (1994) E (cid:104) ˆ K ( k , k ) (cid:105) = K ( k , k ) + O (cid:0) L − (cid:1) . (9)We stated in Equation 3 that the second order cumulant spectrum K ( k , k ) for strictlystationary time series is zero for pairs of frequencies in the frequency domain. Furthermore,the variance of ˆ K , considering P frames, is equal to P − S ( k ) S ( k ) as ( L, P → ∞ ) fromHinich (1994).Next, we define the normalized second order cumulant spectrum asˆΓ ( k , k ) = ˆ K ( k , k ) (cid:112) S ( f k ) S ( f k ) . (10)We make a further simplifying assumption used in Hinich and Rothman (1998) that thetime series has been prewhitened, hence the theoretical spectrum S is assumed constantacross frequencies and, without loss of generality, equal to one. With this assumption, thestandard error of ˆ K becomes 1 / √ P . Then the estimate for Γ, ˆΓ simplifies to √ P ˆ K .Next we construct the complex valued pivotal quantity Y asˆ Y ( k , k ) = √ P (cid:104) ˆΓ ( k , k ) − Γ ( k , k ) (cid:105) . (11)From the central limit theorem, the pivotal quantity Y has an asymptotic complexnormal distribution. Hence, the real value (cid:60) ( Y ) and the imaginary value (cid:61) ( Y ) are bivari-ate standard normal random variables asymptotically. Under the alternative hypothesisof nonstationarity, the theoretical distribution depends on the type of strict stationarityviolation.Notice that we calculate a Y value for each pair of ( k , k ) in the principal domain. Thenumber of calculated cumulants spectra equals the number of pairs ( k , k ) in the principal8omain, which depends on the frame length L . We assume a sufficiently large frame length L so that the asymptotic theory holds.We can test the normality of the stacked vector of real and imaginary values using theKolmogorov Smirnoff (KS) one sample test. Let F emp ( x ) denote the empirical cumulativedistribution of the vector of stacked ( (cid:60) ( Y ) , (cid:61) ( Y )). We calculate F emp ( x ) as the fractionof values that are less than a given value x in the unit interval. In our case the vector x represent the real and imaginary parts for the PHR test statistic and the cumulativedistribution of interest is that of a standard normal variable. The KS statistic is ˆ D = sup 4. Empirical Evaluation of Size and Power of the PHR Test. Comparison to theKPSS Test. Next we proceed to evaluate the empirical size and power of PHR to detect departuresfrom the hypothesis of strict stationarity. We compare PHR to the KPSS test for levelstationarity implemented in the Trapletti and Hornik (2015) tseries R package. Developedin Kwiatkowski et al. (1992), the KPSS test tests the null hypothesis of first order (levelor absence of unit root) stationarity using the test statistic KP SS = 1ˆ σ T T (cid:88) j =1 (cid:32) j (cid:88) t =1 ( x ( t ) − ¯ x ( t )) (cid:33) , (12)where x t is the time series, with t = 1 . . . T , where T is the length of the times series, ¯ x ( t )is the mean of the time series, and ˆ σ is a consistent estimator of the long run variance.Lima and Neri (2013) studied the size and power of KPSS to detect various departuresfrom the strict stationarity hypothesis. Cavaliere and Taylor (2005) evaluated KPSS per-formance to detect unit root in the presence of varying second moment. In the following9ections we replicate several of the data generating processes in Lima and Neri (2013) andCavaliere and Taylor (2005) and compare PHR with KPSS in terms of power to detectdepartures from strict stationarity.We first generated strictly stationary data and evaluated how closely the empirical testsize tracked the nominal test size of five percent and one percent respectively. PHR hasconservative size properties. We then evaluated the power of the test. The null hypothesisof strict stationarity can be violated in various ways, for example unit root, seasonal root,time varying second moment, and long memory. We expect that under large sample sizeconditions PHR detects a wide variety of violations. We limit ourselves to evaluating thepower of PHR to detect unit root, varying second moment, and the combination of thetwo. The power of the test to detect time varying second moment in almost all consideredcases as well as the power to detect unit root in certain cases surpasses the KPSS test.We built an R R Core Team (2015) package to implement PHR (available upon request).In the R package we use the base R implementation of the one sample KS test to test thatthe vector of real and imaginary parts of the ˆ Y test statistic in Equation 11 likely comesfrom the standard normal distribution. The KS test assumes a continuous distributionof the tested sample and a large sample size. The approximate p value of the KS test iscalculated using the technique described in Wang et al. (2003). The p value is approximatefor two reasons. First, the test assume that the sample is large enough for the asymptotictheory to hold, and second, the test assumes that sample values have no ties (but ties arepossible due to the discrete nature of the values).Each data generation process in the size and power study was replicated10,000 timesand PHR and KPSS size and power were calculated as the percentage of times each testrejected the null of strict stationarity under each data generation scenario. For reproducibleresults, the random number generator seed was set and reused for each evaluated scenario.10e performed all the simulations in the R open source software and we created all thetables using the stargazer Hlavac (2015) R package. Simulation code is available uponrequest.Our PHR is built on the assumptions that time series are free from outliers, havefinite moments, do not trend, are prewhitened and have zero means. We evaluated theperformance of both PHR and KPSS test to detect a unit root when data are generated froma normal distribution or a t distribution with three, five and ten degrees of freedom. We alsoevaluated the performance of the tests after first trimming (drop top and bottom 0.5, oneor two percent of the data respectively), demeaning (subtract the mean), and detrendingthe data. The Borchers (2015) pracma package detrended the data by estimating andremoving a least squares fit. Demeaning, detrending, and trimming are commonly used infast Fourier transform applications.Lastly, we generated data with time varying second moment, as well as a combinationof varying second moment and unit root, and evaluated the power of PHR and KPSS. Forthe sake of brevity, we only included partial and essential tables in this manuscript. Allsimulation results are available upon request. First, we generated a white noise process and evaluated PHR empirical size sensitivity tosample size and frame length L . Simulation results in Table 1 show that it is best to havesample sizes above 500, at least 50 frames and a frame length of at least 10. A rule ofthumb is to set L and P equal to the square root of sample size T .Next, we generated strictly stationary data in order to evaluate the empirical sizes ofPHR and KPSS tests. We generated samples from a t distribution with ∞ (normal), three,five, ten and 15 degrees of freedom. Furthermore, we included an autocorrelation structure11able 1: Empirical Test Size for Various Sample Sizes T and Frame Lengths L (correspond-ing number of frames P=[T/L]). Nominal Test Size of 0.05. L T.250 T.500 T.1000 T.500010 0 . 062 0 . 059 0 . 062 0 . . 062 0 . 058 0 . 060 0 . . 195 0 . 084 0 . 063 0 . . 716 0 . 183 0 . . of +0 . − . The stationary data generation process is x ( t ) = u ( t ) . (13)The innovations u ( t ) follow the process u ( t ) = ρu ( t − 1) + v ( t ) , (14)where v ( t ) are independently and identically distributed with one of the t distributions.When ρ is zero, u ( t ) in Equation 14 process reduces to just the v ( t ) component. If ρ isdifferent from zero the time series has short term memory.12 .1.2 Simulation Results for Empirical Test Size. Table 2 shows size of test results for a white noise process with zero autocorrelation in theerror term. We expect the empirical size of the test to be close to the nominal size of thetest. Each entry in the table represents the percentage of times PHR and KPSS rejectedthe null hypothesis of strict stationarity when the null was in fact true, in 10,000 runs. Wealso compared the size of PHR and KPSS after demeaning, trimming and detrending thedata.When we do not preprocess the data, PHR showed slightly higher empirical size thanKPSS, too frequently rejecting the null hypothesis. For a larger sample size of T = 5000the empirical size of PHR is closer to the nominal size of the test. When we trim the dataat different rates, we noticed that PHR’s empirical test size decreased while KPSS sizeremained unchanged. PHR is sensitive to outliers. On the other hand KPSS is extremelysensitive to detrending, achieving zero empirical test size in all cases where we includeddetrending. Neither PHR nor KPSS reacted excessively to demeaning. PHR’s empiricalsize varies with the sample size T , while KPSS’s empirical size is less dependent on samplesize.The empirical sizes of PHR and KPSS are consistent across t distributions with three,five, ten and 15 degrees of freedom and zero autocorrelation. We only included here thetable for five degrees of freedom but results are similar across degrees of freedom. PHR issensitive to fat tails and fares better with appropriate trimming. Table 3 illustrates that aone percent trimming in the tails tames the empirical Type I error rate for PHR. Results(not included) for two percent trimming were similar. The KPSS is less sensitive to fattails but reacts strongly to detrending, as observed in previous tables.Next, we evaluated PHR and KPSS empirical size sensitivity to a ρ = 0 . ρ = 0) and Nominal Size α = 0 . itrim idetrend idemean KPSS.250 PHR.250 KPSS.500 PHR.500 KPSS.1000 PHR.1000 KPSS.5000 PHR.50000 0 0 0 . 049 0 . 062 0 . 052 0 . 084 0 . 050 0 . 063 0 . 052 0 . . 050 0 . 051 0 . 051 0 . 070 0 . 051 0 . 044 0 . 051 0 . . 062 0 0 . 085 0 0 . 064 0 0 . . 049 0 0 . 069 0 0 . 046 0 0 . . 049 0 . 062 0 . 052 0 . 084 0 . 050 0 . 063 0 . 052 0 . . 050 0 . 051 0 . 051 0 . 070 0 . 051 0 . 044 0 . 051 0 . . 062 0 0 . 085 0 0 . 064 0 0 . . 049 0 0 . 069 0 0 . 046 0 0 . Table 3: Empirical Size of Test Comparisons for t5 Data with Zero Autocorrelation inError Term( ρ = 0) and Nominal Size α = 0 . itrim idetrend idemean KPSS.250 PHR.250 KPSS.500 PHR.500 KPSS.1000 PHR.1000 KPSS.5000 PHR.50000 0 0 0 . 046 0 . 152 0 . 048 0 . 266 0 . 048 0 . 247 0 . 048 0 . . 048 0 . 066 0 . 048 0 . 106 0 . 047 0 . 068 0 . 048 0 . . 152 0 0 . 268 0 0 . 246 0 0 . . 001 0 . 067 0 . 001 0 . 105 0 . 001 0 . 066 0 . 001 0 . . 046 0 . 152 0 . 048 0 . 266 0 . 048 0 . 247 0 . 048 0 . . 048 0 . 066 0 . 048 0 . 106 0 . 047 0 . 068 0 . 048 0 . . 152 0 0 . 268 0 0 . 246 0 0 . . 001 0 . 067 0 . 001 0 . 105 0 . 001 0 . 066 0 . 001 0 . t distribution with three, five, tenand 15 degrees of freedom we notice the necessity for trimming. We also notice the strongKPSS reaction to detrending. Notice in Table 4 how when we have large sample size of T = 5000, as required by the asymptotic properties of PHR, PHR’s empirical size is veryclose to nominal size and stable over demeaning, detrending and trimming scenarios whileKPSS’s empirical size varies.Table 4: Empirical Size of Test Comparisons for Normal Data with Autocorrelation inError Term ρ = 0 . α = 0 . itrim idetrend idemean KPSS.250 PHR.250 KPSS.500 PHR.500 KPSS.1000 PHR.1000 KPSS.5000 PHR.50000 0 0 0 . 048 0 . 062 0 . 050 0 . 084 0 . 047 0 . 061 0 . 047 0 . . 045 0 . 055 0 . 051 0 . 067 0 . 048 0 . 044 0 . 049 0 . . 061 0 0 . 084 0 0 . 061 0 0 . . 055 0 0 . 068 0 0 . 045 0 0 . . 048 0 . 062 0 . 050 0 . 084 0 . 047 0 . 061 0 . 047 0 . . 045 0 . 055 0 . 051 0 . 067 0 . 048 0 . 044 0 . 049 0 . . 061 0 0 . 084 0 0 . 061 0 0 . . 055 0 0 . 068 0 0 . 045 0 0 . When the data generation process included negative − . Next we evaluate the power of PHR to detect unit root and compare its power to the powerof KPSS test. 15 .2.1 Data Generation Processes With Unit Root. We replicated the unit root data generation process from Lima and Neri (2013). Thegenerated data included a unit root component y ( t ) and a stationary component u ( t ). Thetuning parameter λ moderates the strength of the unit root process. Then the time series x ( t ) is x ( t ) = λy ( t ) + u ( t ) , (15)with y ( t ) = (cid:80) tj =1 v ( j ).The innovations u ( t ) follow the process in Equation 14. The y ( t ) process is a randomwalk with v ( t ) i.i.d from the same distribution as u ( t ) but independently from u ( t ). Wecompared the power of the PHR to the power of KPSS test to detect a unit root in processeswith relatively low λ of 0 . . 3. A larger λ led to 100% power for both tests, as observedin additional simulation studies not included here. For a process with high lambda (strong presence of unit root), the estimated cumulantspectrum values depart strongly from zero, the spectrum estimate explodes at frequenciesnear zero and the distribution of the resulting test statistic has much heavier tails than astandard normal. Under the null hypothesis of strict stationarity, the test statistic followsan asymptotic normal distribution.The power of PHR is greater for stronger unit root presence (larger λ value), havingdifficulty with smaller λ values. We also evaluated the impact of trimming, demeaning anddetrending on the power of the test.Table 5 illustrates the good performance of PHR in large sample sizes for normal data.16able 5: Power of Test Comparisons for Normal Data with a Unit Root ( λ = 0 . ρ = 0), and Nominal Size of Test α = 0 . itrim idetrend idemean KPSS.5000 PHR.5000 KPSS.250 PHR.250 KPSS.1000 PHR.10000 0 0 1 0 . 997 0 . 963 0 . 080 0 . 997 0 . . 996 0 . 956 0 . 072 0 . 996 0 . . 998 0 . 997 0 . 620 0 . 078 0 . 936 0 . . 998 0 . 997 0 . 597 0 . 072 0 . 928 0 . . 997 0 . 963 0 . 080 0 . 997 0 . . 996 0 . 956 0 . 072 0 . 996 0 . . 998 0 . 997 0 . 620 0 . 078 0 . 936 0 . . 998 0 . 997 0 . 597 0 . 072 0 . 928 0 . The power is very close to one for both PHR and KPSS in sample sizes of 5000. The powerof PHR deteriorates at lower sample sizes while KPSS holds its ground. Power deterioratesfor both KPSS and PHR with detrending and trimming.Table 6: Power of Test Comparisons for t5 Data with Unit Root ( λ = 0 . ρ = 0), and Nominal Size of Test α = 0 . itrim idetrend idemean KPSS.250 PHR.250 KPSS.1000 PHR.10000 0 0 0 . 968 0 . 148 0 . 998 0 . . 965 0 . 106 0 . 997 0 . . 628 0 . 148 0 . 939 0 . . 623 0 . 098 0 . 935 0 . . 968 0 . 148 0 . 998 0 . . 965 0 . 106 0 . 997 0 . . 628 0 . 148 0 . 939 0 . . 623 0 . 098 0 . 935 0 . In Table 6, we generated data from a t distribution with five degrees of freedom, no17orrelation in error term and a unit root with λ = 0 . 3. PHR performed similarly to thenormal case, displaying the already established sensitivity to sample size and fat tails. InTable 6 we see slightly higher power to detect the unit root process in smaller sample sizesas a reaction to fat tails. It is also interesting to notice how the power of KPSS to detectunit root decreases significantly with small sample sizes such as T = 250 when we applyhigher trimming, and markedly decreases with detrending.PHR and KPSS performed similarly when we included autocorrelation in the error term,as well as when λ decreased to 0 . Another alternative to the null hypothesis of stationarity is a varying second moment.The authors of Cavaliere and Taylor (2005) evaluated the power of the KPSS test andtwo other tests to detect a unit root in the presence of varying variance. We employedtheir data generating processes in a similar but not identical fashion and evaluated thepower of PHR to detect the alternative hypothesis of time varying second moment forthree different patterns of nonconstant variance. The null hypothesis is the hypothesis ofstrict stationarity of the time series, when the second order cumulant spectrum is equal tozero. The alternative hypothesis is nonstationarity due to varying second moment. Hence,our set up is different from Cavaliere and Taylor (2005), and similar to Lima and Neri(2013).Lima and Neri (2013) also evaluated several tests for strict stationarity and their poweragainst the alternative of unconditional heteroskedasticity. Lima and Neri (2013) statedthat tests for level stationarity (or unit root) such as the KPSS test (and its robust versioniKPSS) have low power against the alternative of time varying variance. Furthermore they18dded that KPSS is designed based on the fluctuation of the data around the sample meanand it has low power to detect how the scale of the distribution varies over time. A processwith time varying scale has varying second moment and hence is not strictly stationary.Therefore we are interested in the power of PHR to detect a varying second moment. Limaand Neri (2013) only discussed one pattern of heteroskedasticity or varying second momentwhile we considered the three patterns suggested in Cavaliere and Taylor (2005). Next, we describe the data generating processes with time varying variance patterns andthe simulation results for the power of PHR and KPSS test. The first data generatingprocess (DGP1) in Equation 16 is for a level stationary (no unit root) process with variance σ t . Under the null hypothesis of time constant variance, σ t = 1. Under the alternativehypothesis of varying variance, σ t follows one of the three processes in Equation 18. DGP1is x ( t ) = σ t u ( t ) , (16)with t = 1 . . . T and u ( t ) a white noise.The second data generating process (DGP2) in Equation 17 generates a process withunit root and variance σ t . Similarly to DGP1, σ t = 1 under the null hypothesis of constantvariance. Under the alternative hypothesis of varying variance, σ t follows one of the threeprocesses in Equation 18. DGP2 is x ( t ) = x ( t − 1) + σ t u ( t ) , (17)where u ( t ) are white noise and x (0) = 0. When we evaluated the power of PHR andKPSS tests to detect violations in DGP2, we effectively study two confounding sources19f nonstationarity, both unit root and varying second moment. In this case, the tests donot indicate which one of the sources of nonstationarity caused the rejection of the nullhypothesis.For clarity, we only allow σ t to vary according to the three patterns described below,while the initial paper by Cavaliere and Taylor (2005) allowed the same patterns in theirregular component u ( t ) as well. Hence, our simulation results will naturally differ fromtheirs. Equation 18 defines the three variance patterns as: single and multiple breaks invariance, smooth transition variance breaks and linear trending variance. Notice how thevariance is constant up to a certain point in time, m , after which it evolves according toone of the three functions. • Single Break σ t = σ + ( σ − σ ) , t ≥ mσ , t < m where σ = 1 and σ = cσ . • Smooth Transition σ t = σ + (cid:0) σ − σ (cid:1) W ( t ) , (18)where W ( t ) = (1 + exp ( − γ ( t − m ))) − , σ = 1 and σ = cσ . • Piecewise Linear Trend σ t = σ + ( σ − σ ) ( t − m ) (1 − m ) − , t ≥ mσ , t < m where σ = 1 and σ = cσ . 20he tuning parameters c and m dictate the scale and time of the change in varianceunder the alternative hypothesis of varying variance. We considered c values of { . , } ,corresponding to a decreased or increased variance, respectively. For the single breakand smooth transition cases, m takes the values { . , . , . } . For a sample size T of1000 observations, m = 0 . γ = 10 parameter in Equation 18 gives the speed or smoothness of transitionfrom constant variance for the smooth transition variance pattern.For the piecewise linear trend in Equation 19 m takes the values { , . , . } . In thiscase, when m > mT and then it changesto an increasing (when c = 4) or decreasing (when c = 0 . 25) linear trend. Higher m dictatesa later linear trend with a higher slope. When m = 0 the variance process follows a lineartrend starting at t = 1.Figure 1 displays the patterns in σ for each of single break, smooth transition andpiecewise linear trend. In each of the four plots, variance increases, corresponding to theshift parameter c = 4 and the time series sample is 1000 observations long. For the singlebreak and smooth transition cases, m = 0 . m = 0 . m = 0.Figure 2 illustrates both DGP1 and DGP2 processes with time varying variance accord-ing to the three patterns. We present next Monte Carlo simulation results to compare the power of PHR and KPSStest to detect time varying variance. Table 7 contains the power of the two tests to detectpiecewise linear trend in second moment (increase if c = 4 and decrease if c = 0 . 25) that21 200 400 600 800 1000 . . . . . . . m = 0.5 time index s i g m a Single BreakSmooth Transition . . . . . . . Piecewise Linear Trend time index s i g m a m = 0.5m = 0 Figure 1: The Three Functions for a Variance Increase and T = 1000. Left: Single Breakand Smooth Transition. Right: Piecewise Linear Trend.starts at different points in the sample (early if m = 0, half way if m = 0 . m = 0 . T = { , , } , and at five percent nominalsize of the test. PHR outperformed KPSS. Both PHR and KPSS disappointed for a fewcombinations of parameters m and c . These cases are for DGP1 process (no unit root) and c = 4 (increase in variance). Even smaller sample sizes fared well for PHR power. For largesample sizes, the two tests agreed.For T = 1000, PHR achieved 100% power in almost all cases, except when detecting avariance decrease that started at a later stage. The KPSS test had close to 99% power forDGP2, which is a process that includes both a unit root and time varying variance. Forthe DGP1 process, when we have time varying variance and no unit root, KPSS had powerclose to 99%. The power patterns for both PHR and KPSS are almost identical for the T = 5000 sample size. PHR power was slightly less but still above 95% when the samplesize was T = 250. However, when the sample size was T = 250 with varying variance22 50 100 150 200 250 − − − Level Stationary and Single Break in Variance time points dgp1 − Unit Root and Single Break time points dgp2 − − − Level Stationary and Smooth Transition Variance time points dgp1 − Unit Root and Smooth Transition Variance time points dgp2 − − Level Stationary and Piecewise Linear Variance time points dgp1 − − − Unit Root and Piecewise Linear Variance time points dgp2 Figure 2: Time Series with Time Varying Variance23able 7: Power of Tests to Detect Piecewise Linear Trend at Nominal Test Size α = 0 . dgp1 2 m c PHR.250 KPSS.250 PHR.1000 KPSS.1000 PHR.5000 KPSS.50001 0 4 0 . 090 0 . 054 0 . 742 0 . 056 1 0 . . 995 0 . 982 1 0 . 998 1 11 0 0.25 0 . 094 0 . 050 0 . 762 0 . 051 1 0 . . 996 0 . 982 1 0 . 997 1 11 0.5 4 0 . 273 0 . 098 1 0 . 099 1 0 . . 987 0 . 966 1 0 . 995 1 11 0.5 0.25 0 . 078 0 . 036 0 . 378 0 . 036 0 . 956 0 . . 997 0 . 984 1 0 . 998 1 11 0.9 4 0 . 447 0 . 114 1 0 . 120 1 0 . . 990 0 . 974 1 0 . 996 1 11 0.9 0.25 0 . 066 0 . 042 0 . 084 0 . 042 0 . 105 0 . . 998 0 . 984 1 0 . 998 1 1 and no unit root, PHR had lower power. We know that the construction of PHR relies onasymptotical properties and large enough sample sizes and frame lengths so a diminishedpower is not surprising for small sample sizes. Since KPSS was designed to detect violationsdue to unit root, perhaps its underperformance detecting varying variance is not surprising.However, PHR had greater power to detect both the combination of unit root and varyingvariance violations as well as only varying variance violation.Simulation results for the cases with a single break in variance and with a smoothtransition in variance are consistent with the power results in Table 7 and we skipped thetables for the sake of brevity but have them available upon request. 5. Application. High Frequency Stock Market Returns. When conducting statistical analyses of stock market rates of return, finance academicstypically assume that the return series are normally distributed and weakly stationary. If24 time series is stationary, the k time step ahead forecasts of the time series converge tothe time series mean, and the variance of the forecasting error converges to the time seriesvariance, as the forecasting horizon lengthens. A test that detects evidence against weakstationarity also detects evidence against strict stationarity. On the other hand, if we fail toreject weak stationarity, we can still reject strict stationarity, as strict stationarity requiresthat higher order moments are also time invariant. Next we illustrate the test power todetect nonstationarity in high frequency returns for one New York Stock Exchange (NYSE)traded stock, the Ford Motor Company. We obtained our data from Wharton ResearchData Services.We considered every Ford stock trade for the years 2001, 2002, and 2003. With tradesoccurring at random times during market trading hours (9:30AM to 4:00PM), we sampledprices at equal length intervals. We then transformed the prices by taking the first differenceof the natural logarithm in order to produce the stock rates of return. When conducting astatistical analysis of stock market rates of return, finance academics typically assume thatthe data are normally distributed and weakly stationary. In order to guard against outliersunduly affecting the test we trimmed the lowest and largest one percent of the observations.We had over 19,000 observations per year, after trimming. We also detrended the data inorder to remove any potential subtle trend in the trimmed returns. Table 8 shows theFord stock descriptive statistics, before any data processing. The means and standarddeviations are typical for high frequency return data. Positive kurtosis is also a commoncharacteristic. The series has small negative skewness in all years. Figure 3 shows one yearof return data with apparent time varying variance.For each evaluated year of data, PHR rejected the null hypothesis of stationarity at fivepercent significance level before and after demeaning, detrending and trimming. The KPSStest of stationarity did not detect this pattern before trimming but did after trimming at25 − . − . − . . . High Frequency Return Data for Ford Stock For 2002 time points − m i nu t e r e t u r n s Figure 3: High Frequency Ford Returns 2002one percent level.When we adjusted the data using a rolling window estimate of variance like in Limaand Neri (2013), the time varying variance violation was smoothed and PHR concurredwith KPSS in failing to reject stationarity. We estimated the smoothing correction factorusing a rolling window of 120 days.Table 8: Descriptive Statistics for Ford Stock 5-Minute Returns Year N Mean St. Dev. Min Max Skew. Kurt.2001 19,344 0.00000 0.010 − − − − − − − . Conclusion. In this article we used the frame averaged estimates of the second order cumulant spectrumcalculated for pairs of fundamental Fourier frequencies in the principal domain to build astatistical test of the null hypothesis that a stochastic process is strictly stationary. Weassume that the process has zero mean, is free of outliers, has no trend and has beenprewhitened. One example of a stationary process is white noise. We based the testdevelopment on the property that the second order cumulant spectrum of a time series inthe frequency domain is zero under strict stationarity. An application of the test to fiveminute sampled high frequency Ford stock returns showed nonstationarity.We evaluated the size and power properties of PHR and KPSS via Monte Carlo sim-ulations. We also evaluated the sensitivity of PHR to sample sizes, fat tails in the timeseries distribution, and autocorrelation in the error term. Some limitations of PHR aresensitivity to outliers, lower power to detect subtle unit root processes in smaller samplesizes, and higher Type I error rate in certain data scenarios. PHR had good power todetect second moment variations, and outperformed KPSS test in all considered scenar-ios. PHR also showed good power to detect a unit root in the presence of a time varyingsecond moment. In future research it would be interesting to study the power of PHR todetect other types of stationarity violations, such as long memory violations as describedin Ashley and Patterson (2010). PHR can be included in a battery of tests used to verifythe strict stationarity assumptions of a time series. References Ashley, R. A. and Patterson, D. M. (2010). Apparent long memory in time series as anartifact of a time-varying mean: Considering alternatives to the fractionally integrated27odel. Macroeconomic Dynamics , 14(S1):59–87.Bloomfield, P. (2004). Fourier Analysis of Time Series: an Introduction . John Wiley &Sons.Borchers, H. W. (2015). pracma: Practical Numerical Math Functions . R package version1.8.8.Brillinger, D. R. (2001). Time Series: Data Analysis and Theory , volume 36. Siam.Cavaliere, G. and Taylor, A. R. (2005). Stationarity tests under time-varying second mo-ments. Econometric Theory , 21(06):1112–1129.Dickey, D. A. and Fuller, W. A. (1979). Distribution of the estimators for autoregres-sive time series with a unit root. Journal of the American Statistical Association ,74(366a):427–431.Dickinson, R. M., Roberts, D. A. O., Driscoll, A. R., Woodall, W. H., and Vining, G. G.(2014). Cusum charts for monitoring the characteristic life of censored Weibull lifetimes. Journal of Quality Technology , 46(4):340.Farley, J. U. and Hinich, M. J. (1970). Detecting small mean shifts in time series. Man-agement Science , 17(3):189–199.Hawkins, D. M. (1977). Testing a sequence of observations for a shift in location. Journalof the American Statistical Association , 72(357):180–186.Hinich, M. A. and Wild, P. (2001). Testing time-series stationarity against an alternativewhose mean is periodic. Macroeconomic Dynamics , 5(03):380–412.28inich, M. J. (1994). Higher order cumulants and cumulant spectra. Circuits, Systems andSignal Processing , 13(4):391–402.Hinich, M. J. and Rothman, P. (1998). Frequency-domain test of time reversibility. Macroe-conomic Dynamics , 2(01):72–88.Hlavac, M. (2015). stargazer: Well-Formatted Regression and Summary Statistics Tables .Harvard University, Cambridge, USA. R package version 5.2.Kwiatkowski, D., Phillips, P. C., Schmidt, P., and Shin, Y. (1992). Testing the null hypoth-esis of stationarity against the alternative of a unit root: How sure are we that economictime series have a unit root? Journal of Econometrics , 54(1):159–178.Lima, L. R. and Neri, B. (2013). A test for strict stationarity. Uncertainty Analysis inEconometrics with Applications , pages 17–30.Phillips, P. C. and Perron, P. (1988). Testing for a unit root in time series regression. Biometrika , pages 335–346.Press, W., Flannery, B., Teukolsky, S., and Vetterling, W. (1992). Kolmogorov-Smirnovtest. Numerical Recipes in FORTRAN: The Art of Scientific Computing , pages 617–620.R Core Team (2015). R: A Language and Environment for Statistical Computing . RFoundation for Statistical Computing, Vienna, Austria.Trapletti, A. and Hornik, K. (2015). tseries: Time Series Analysis and ComputationalFinance . R package version 0.10-34.Wang, J., Tsang, W. W., and Marsaglia, G. (2003). Evaluating Kolmogorov’s distribution.