Spectral Inference under Complex Temporal Dynamics

TL;DR

This paper introduces a unified framework for inferring evolutionary Fourier power spectra in complex, possibly nonlinear, locally stationary processes, providing statistically valid confidence regions and versatile tools for time-frequency analysis.

Contribution

It develops a comprehensive theory and bootstrap-based methodology for constructing confidence regions for spectral densities in non-stationary processes, applicable to various inference tasks.

Findings

Constructs asymptotically correct simultaneous confidence regions for spectral densities.

Proposes a bootstrap method for practical implementation of confidence regions.

Enables visual and statistical evaluation of time-varying spectral patterns.

Abstract

We develop unified theory and methodology for the inference of evolutionary Fourier power spectra for a general class of locally stationary and possibly nonlinear processes. In particular, simultaneous confidence regions (SCR) with asymptotically correct coverage rates are constructed for the evolutionary spectral densities on a nearly optimally dense grid of the joint time-frequency domain. A simulation based bootstrap method is proposed to implement the SCR. The SCR enables researchers and practitioners to visually evaluate the magnitude and pattern of the evolutionary power spectra with asymptotically accurate statistical guarantee. The SCR also serves as a unified tool for a wide range of statistical inference problems in time-frequency analysis ranging from tests for white noise, stationarity and time-frequency separability to the validation for non-stationary linear models.

Tables

Table 1. Table 1. Simulated Uncoverage Probabilities for Example   7.1

		$N=400$		$N=800$
$n$	$B_n$	$\alpha =0.05$	$\alpha =0.1$	$\alpha =0.05$	$\alpha =0.1$
$72$	$36$	$0.03$	$0.06$	$0.03$	$0.06$
$72$	$32$	$0.04$	$0.07$	$0.04$	$0.08$
$72$	$28$	$0.04$	$0.08$	$0.05$	$0.10$
$54$	$36$	$0.03$	$0.06$	$0.03$	$0.06$
$54$	$32$	$0.04$	$0.08$	$0.04$	$0.08$
$54$	$28$	$0.04$	$0.09$	$0.05$	$0.09$
$36$	$32$	$0.05$	$0.11$	$0.05$	$0.11$
$36$	$28$	$0.07$	$0.13$	$0.07$	$0.14$

Table 2. Table 2. Simulated Uncoverage Probabilities for Example   7.2

		$N=400$		$N=800$
$n$	$B_n$	$\alpha =0.05$	$\alpha =0.1$	$\alpha =0.05$	$\alpha =0.1$
$72$	$36$	$0.03$	$0.06$	$0.02$	$0.05$
$72$	$32$	$0.04$	$0.08$	$0.04$	$0.08$
$72$	$28$	$0.05$	$0.10$	$0.05$	$0.11$
$54$	$36$	$0.03$	$0.06$	$0.03$	$0.07$
$54$	$32$	$0.05$	$0.09$	$0.04$	$0.09$
$54$	$28$	$0.05$	$0.10$	$0.04$	$0.10$
$36$	$32$	$0.06$	$0.12$	$0.06$	$0.12$
$36$	$28$	$0.08$	$0.15$	$0.08$	$0.14$

Table 3. Table 3. Simulated Uncoverage Probabilities for Example   7.3

		$N=400$		$N=800$
$n$	$B_n$	$\alpha =0.05$	$\alpha =0.1$	$\alpha =0.05$	$\alpha =0.1$
$72$	$36$	$0.04$	$0.09$	$0.04$	$0.09$
$72$	$32$	$0.05$	$0.11$	$0.05$	$0.12$
$72$	$28$	$0.06$	$0.12$	$0.06$	$0.13$
$54$	$36$	$0.06$	$0.10$	$0.05$	$0.10$
$54$	$32$	$0.06$	$0.11$	$0.06$	$0.11$
$54$	$28$	$0.06$	$0.12$	$0.07$	$0.13$
$36$	$32$	$0.07$	$0.14$	$0.08$	$0.14$
$36$	$28$	$0.07$	$0.14$	$0.08$	$0.15$

Table 4. Table 4. Simulated Uncoverage Probabilities for Example   7.4

		$N=400$		$N=800$
$n$	$B_n$	$\alpha =0.05$	$\alpha =0.1$	$\alpha =0.05$	$\alpha =0.1$
$72$	$36$	$0.05$	$0.11$	$0.06$	$0.12$
$72$	$32$	$0.06$	$0.12$	$0.06$	$0.13$
$72$	$28$	$0.07$	$0.14$	$0.08$	$0.15$
$54$	$36$	$0.05$	$0.10$	$0.06$	$0.12$
$54$	$32$	$0.05$	$0.11$	$0.06$	$0.12$
$54$	$28$	$0.06$	$0.12$	$0.07$	$0.13$
$36$	$32$	$0.08$	$0.14$	$0.08$	$0.14$
$36$	$28$	$0.08$	$0.14$	$0.09$	$0.17$

Table 5. Table 5. Simulated Uncoverage Probabilities for Example   7.5

		$N=400$		$N=800$
$n$	$B_n$	$\alpha =0.05$	$\alpha =0.1$	$\alpha =0.05$	$\alpha =0.1$
$72$	$36$	$0.04$	$0.08$	$0.05$	$0.08$
$72$	$32$	$0.05$	$0.11$	$0.05$	$0.10$
$72$	$28$	$0.06$	$0.13$	$0.05$	$0.11$
$54$	$36$	$0.03$	$0.07$	$0.04$	$0.09$
$54$	$32$	$0.05$	$0.10$	$0.06$	$0.11$
$54$	$28$	$0.06$	$0.13$	$0.07$	$0.13$
$36$	$32$	$0.06$	$0.11$	$0.07$	$0.14$
$36$	$28$	$0.09$	$0.16$	$0.09$	$0.17$

Table 6. Table 6. Simulated Uncoverage Probabilities with Tuning Parameters Selected by the MV Method (Numbers in the Parentheses Represent the Average Width of the SCR)

	$N=400$				$N=800$
	$n$	$B_n$	$\alpha =0.05$	$\alpha =0.1$	$n$	$B_n$	$\alpha =0.05$	$\alpha =0.1$
Example 7.1 (Table 1)	$54$	$30$	$0.06$ ($0.51$)	$0.10$ ($0.45$)	$72$	$30$	$0.05$ ($0.47$)	$0.09$ ($0.43$)
Example 7.2 (Table 2)	$52$	$32$	$0.05$ ($0.37$)	$0.09$ ($0.33$)	$69$	$30$	$0.04$ ($0.35$)	$0.09$ ($0.30$)
Example 7.3 (Table 3)	$54$	$32$	$0.06$ ($0.52$)	$0.11$ ($0.47$)	$72$	$32$	$0.05$ ($0.49$)	$0.12$ ($0.43$)
Example 7.4 (Table 4)	$50$	$32$	$0.06$ ($0.53$)	$0.12$ ($0.45$)	$70$	$32$	$0.06$ ($0.48$)	$0.12$ ($0.43$)
Example 7.5 (Table 5)	$52$	$32$	$0.04$ ($0.52$)	$0.09$ ($0.46$)	$69$	$31$	$0.06$ ($0.44$)	$0.11$ ($0.40$)

Table 7. Table 7. Simulated Accuracy of Hypothesis Testing

Nominal Level		$\alpha =0.05$	$\alpha =0.1$		$\alpha =0.05$	$\alpha =0.1$
Example 7.6	$N=400$	$0.06$	$0.12$	$N=800$	$0.04$	$0.09$
Example 7.7	$N=800$	$0.05$	$0.10$	$N=1200$	$0.04$	$0.09$
Example 7.8	$N=400$	$0.04$	$0.09$	$N=800$	$0.06$	$0.12$

Table 8. Table 8. Real Data: p-values for Testing (a) Stationarity, (b) Time-Varying White Noise, (c) Time-Frequency Separability (Correlation Stationarity).

$H_0$	Stationarity	TV White Noise	Separability
Earthquake	${0.0011}^{\ast \ast }$	${0.012}^{\ast }$	${0.064}^{+}$
Explosion	${0.0005}^{\ast \ast \ast }$	${0.033}^{\ast }$	$0.61$
SP500	${0.0001}^{\ast \ast \ast }$	$0.99$	$0.99$
SP500 (Abs)	${0.0004}^{\ast \ast \ast }$	${0.037}^{\ast }$	${0.048}^{\ast }$

Table 9. Table 9. p-values for Validating Time-varying ARMA Models to Absolute SP500

Model	p-value	Model	p-value
tv-AR$\left(1\right)$	${0.0066}^{\ast \ast }$	tv-ARMA$(1,1)$	${0.019}^{\ast }$
tv-AR$\left(2\right)$	${0.0015}^{\ast \ast }$	tv-ARMA$(2,1)$	$0.79$
tv-AR$\left(3\right)$	${0.0015}^{\ast \ast }$	tv-ARMA$(3,1)$	$0.77$
tv-AR$\left(4\right)$	${0.0012}^{\ast \ast }$	tv-ARMA$(4,1)$	$0.78$
tv-AR$\left(5\right)$	${0.0012}^{\ast \ast }$	tv-ARMA$(5,1)$	$0.84$

Equations

$$X_i^{(N)}=G\left(i/N,{\mathcal{F}}_i\right),i=1,\dots ,N,$$

$$\Vert G(u,{\mathcal{F}}_i)-G(s,{\mathcal{F}}_i){\Vert }_q\le C|u-s|.$$

$$G(u,{\mathcal{F}}_i)=\sum _{j=0}^{\infty }{a}_j(u){ϵ}_{i-j},$$

$$G(u,{\mathcal{F}}_i)=a(u)[G(u,{\mathcal{F}}_{i-1}){]}^{+}+b(u)[-G(u,{\mathcal{F}}_{i-1}){]}^{+}+{ϵ}_i,0\le u\le 1,$$

$${\delta }_p(k):=\underset{0\le u\le 1}{\mathrm{supinfsup}}\Vert G(u,{\mathcal{F}}_k)-G(u,({\mathcal{F}}_{-1},{ϵ}_0^{\prime },{ϵ}_1,\dots ,{ϵ}_k)){\Vert }_p.$$

$$X_i=R(X_{i-1},{ϵ}_i),$$

$$X_i(u)=R(u,X_{i-1}(u),{ϵ}_i),0\le u\le 1.$$

$$\underset{u\in [0,1]}{\mathrm{supinfsup}}\underset{x\ne y}{\mathrm{supinfsup}}\frac{\Vert R(u,x,{ϵ}_0)-R(u,y,{ϵ}_0){\Vert }_{\alpha }}{|x-y|}<1.$$

$$r(u,k):=\text{Cov}\left(G(u,{\mathcal{F}}_0),G(u,{\mathcal{F}}_k)\right).$$

$$|r(u,k)-r(s,k)|/|u-s|\le C,$$

$$r(u+{\delta }_u,k)-r(u,k)=\mathcal{O}(n/N),\forall -n/N\le {\delta }_u\le n/N.$$

$$f(u,\theta ):=\frac{1}{2\pi }\sum _{k\in \mathbb{Z}}r(u,k)\exp (\sqrt{-1}k\theta ).$$

$$J_n(u,\theta ):=\sum _{i=1}^N\tau \left(\frac{i-\lfloor uN\rfloor }{n}\right)X_i\exp (\sqrt{-1}\theta i).$$

$$I_n(u,\theta ):=\frac{1}{2\pi n}|J_n(u,\theta ){|}^2.$$

$$\hat{r}(u,k):=\frac{1}{n}\sum _{i=1}^N\tau \left(\frac{i-\lfloor uN\rfloor }{n}\right)\tau \left(\frac{i+k-\lfloor uN\rfloor }{n}\right)X_i{X}_{i+k},$$

$$I_n(u,\theta )=\frac{1}{2\pi }\sum _{k=-n}^n\hat{r}(u,k)\exp (\sqrt{-1}\theta k).$$

$${\hat{f}}_n(u,\theta ):=\frac{1}{2\pi }\sum _{k=-B_n}^{B_n}\hat{r}(u,k)a(k/B_n)\exp (\sqrt{-1}k\theta ).$$

$$\begin{array}{rl}{Z}_{u,j}^{(n)}& =\frac{{\sum }_{k=1}^N\tau \left(\frac{k-\lfloor uN\rfloor }{n}\right)X_k\cos (k2\pi j/n)}{\sqrt{\pi nf(u,2\pi j/n)}},\\ Z_{u,j+m}^{(n)}& =\frac{{\sum }_{k=1}^N\tau \left(\frac{k-\lfloor uN\rfloor }{n}\right)X_k\sin (k2\pi j/n)}{\sqrt{\pi nf(u,2\pi j/n)}},j=1,\dots ,m,\end{array}$$

$$Z_{U,J}=(Z_{u_1,j_1}^{(n)},\dots ,Z_{u_1,j_p}^{(n)},\dots ,Z_{u_q,j_1}^{(n)},\dots ,Z_{u_q,j_p}^{(n)}{)}^T$$

$$\underset{J}{\mathrm{supinfsup}}\underset{c\in {\Omega }_{p,q}}{\mathrm{supinfsup}}\underset{x}{\mathrm{supinfsup}}|P(c^T{Z}_{U,J}\le x)-\Phi (x)|=o(1),$$

$$\underset{u}{\mathrm{supinfsup}}\underset{\theta \in [0,\pi ]}{\mathrm{maxinfsup}}\sqrt{n/B_n}|{\hat{f}}_n(u,\theta )-\mathbb{E}({\hat{f}}_n(u,\theta ))|={\mathcal{O}}_{\mathbb{P}}(\sqrt{\log n}).$$

$$\underset{u\to 0}{\mathrm{liminfsup}}{u}^{-2}[1-a(u)]=C,$$

$$\underset{u}{\mathrm{supinfsup}}\underset{\theta }{\mathrm{supinfsup}}\left[\mathbb{E}{\hat{f}}_n(u,\theta )-f(u,\theta )-\frac{C}{B_n^2}{f}^{\prime \prime }(u,\theta )\right]=o(1/B_n^2),$$

$$\sqrt{n/B_n}\{{\hat{f}}_n(u,\theta )-\mathbb{E}({\hat{f}}_n(u,\theta ))\}\Rightarrow \mathcal{N}(0,{\sigma }_u^2(\theta )),$$

$$\begin{array}{rl}& \mathbb{P}[\underset{(u,\theta )\in {\mathcal{G}}_N}{\mathrm{maxinfsup}}\frac{n}{B_n}\frac{|{\hat{f}}_n(u,\theta )-\mathbb{E}({\hat{f}}_n(u,\theta )){|}^2}{f^2(u,\theta ){\int }_{-1}^1{a}^2(t)\mathrm{d}t}\\ & -2\log B_n-2\log C_n+\log (\pi \log B_n+\pi \log C_n)\le x]\to e^{-e^{-x/2}}.\end{array}$$

$$\begin{array}{rl}\underset{\{u_i,{\theta }_j\}\in {\tilde{\mathcal{G}}}_N}{\mathrm{maxinfsup}}\underset{\{u:|u-u_i|\le {\delta }_{u,n},\theta :|\theta -{\theta }_j|\le {\delta }_{\theta ,n}\}}{\mathrm{supinfsup}}& \sqrt{n/B_n}\left|\left[{\hat{f}}_n(u,\theta )-\mathbb{E}({\hat{f}}_n(u,\theta ))\right]-\left[{\hat{f}}_n(u_i,{\theta }_j)-\mathbb{E}({\hat{f}}_n(u_i,{\theta }_j))\right]\right|\\ & =o_{\mathbb{P}}(\sqrt{\log n}).\end{array}$$

$$\begin{array}{rl}& \sqrt{n/B_n}\left|\underset{(u,\theta )\in (0,1)\times [0,\pi )}{\mathrm{supinfsup}}\left|{\hat{f}}_n(u,\theta )-\mathbb{E}({\hat{f}}_n(u,\theta ))\right|-\underset{\{u_i,{\theta }_j\}\in {\tilde{\mathcal{G}}}_N}{\mathrm{maxinfsup}}\left|{\hat{f}}_n(u_i,{\theta }_j)-\mathbb{E}({\hat{f}}_n(u_i,{\theta }_j))\right|\right|=o_{\mathbb{P}}(\sqrt{\log n}).\end{array}$$

$$H_0:\forall \theta ,f(u,\theta )=g(u),u\in [0,1]$$

$$g_0(u):=\arg \underset{\tilde{g}}{\mathrm{mininfsup}}\frac{1}{\pi }{\int }_0^{\pi }|f(u,\theta )-\tilde{g}(u){|}^2\mathrm{d}\theta .$$

$$\hat{g}(u):=\frac{1}{\pi }{\int }_0^{\pi }{\hat{f}}_n(u,\theta )\mathrm{d}\theta \approx \frac{1}{\pi }{\int }_0^{\pi }f(u,\theta )\mathrm{d}\theta =g_0(u).$$

Full text

1. 1 Introduction
2. 2 Problem Formulation
3. 3 Fourier Transforms
4. 4 Consistency and Asymptotic Normality
5. 5 Maximum Deviations
6. 6 Bootstrap and Tuning Parameter Selection
7. 7 Simulations and Data Analysis
8. 8 Proofs of Main Results
9. A Supplemental Material

1. Introduction

It is well known that the frequency content of many real-world stochastic processes evolves over time. Motivated by the limitations of the traditional spectral methods in analyzing non-stationary signals, time-frequency analysis has become one of the major research areas in applied mathematics and signal processing [Coh95, Grö01, Dau90]. Based on various models or representations of the non-stationary signal and its time-varying spectra, time-frequency analysis aims at depicting temporal and spectral information simultaneously and jointly. Roughly speaking, there are three major classes of algorithms in time-frequency analysis: linear algorithms such as short time Fourier transforms (STFT) and wavelet transforms [All77, Mey92, Dau92]; bilinear time-frequency representations such as the Wigner–Ville distribution and more generally the Cohen’s class of bilinear time–frequency distributions [Coh95, HBB92] and nonlinear algorithms such as the empirical mode decomposition method [HSLW+98] and the synchrosqueezing transform [DLW11]. Though there exists a vast literature on defining and estimating the time-varying frequency content, statistical inference such as confidence region construction and hypothesis testing has been paid little attention to in time-frequency analysis.

It is clear that the subject and the goals of time-frequency analysis and non-stationary time series analysis are highly overlapped. Unfortunately it seems that the non-stationary spectral domain theory and methodology in the time series literature have been developed largely independently from time-frequency analysis. One major effort in non-stationary time series analysis lies in forming general classes of non-stationary time series models through their evolutionary spectral representation. Among others, [Pri65] proposed the notion of evolutionary spectra in a seminar paper. In another seminal work, [Dah97] defined a general and theoretically tractable class of locally stationary time series models based on their time-varying spectral representation. [NvK00] studied a class of locally stationary time series from an evolutionary wavelet spectrum perspective and investigated the estimation of the latter spectrum. A second line of research in the non-stationary spectral domain literature involves adaptive estimation of the evolutionary spectra. See for instance [Ada98] for a binary segmentation based method, [ORSM01] for an automatic estimation procedure based on the smooth localized complex exponential (SLEX) transform and [FN06] for a Haar–Fisz technique for the estimation of the evolutionary wavelet spectra. On the statistical inference side, there exists a small number of papers utilizing the notion of evolutionary spectra to test some properties, especially second order stationarity, of a time series. See for instance [Pap10, DPV11, DSR11, JSR15] for tests of stationarity based on properties of the Fourier periodogram or spectral density. See also [Nas13] for a test of stationarity based on the evolutionary wavelet spectra. On the other hand, however, to date there have been no results on the joint and simultaneous inference of the evolutionary spectrum itself for general classes of non-stationary and possibly nonlinear time series to the best of our knowledge.

The purpose of the paper is to develop a unified theory and methodology for the joint and simultaneous inference of the evolutionary spectral densities for a general class of locally stationary and possibly nonlinear processes. From a time-frequency analysis perspective, the purpose of the paper is to provide a unified and asymptotically correct method for the simultaneous statistical inference of the STFT-based evolutionary power spectra, one of the most classic and fundamental algorithms in time-frequency analysis. Let $\{X^{(N)}_{i}\}_{i=1}^{N}$ be the observed time series or signal. One major contribution of the paper is that we establish a maximum deviation theory for the STFT-based spectral density estimates over a nearly optimally dense grid ${\mathcal{G}}_{N}$ in the joint time-frequency domain. Here the optimality of the grid refers to the best balance between computational burden and (asymptotic) correctness in depicting the overall time-frequency stochastic variation of the estimates. We refer the readers to Section 5.1 for a detailed definition and discussion of the optimality. The theory is established for a very general class of possibly nonlinear locally stationary processes which admit a time-varying physical representation in the sense of [ZW09] and serves as a foundation for the joint and simultaneous time-frequency inference of evolutionary spectral densities. Specifically, we are able to prove that the spectral density estimates on ${\mathcal{G}}_{N}$ are asymptotically independent quadratic forms of $\{X_{i}^{(N)}\}_{i=1}^{N}$. Consequently, the maximum deviation of the spectral density estimates on ${\mathcal{G}}_{N}$ behaves asymptotically like a Gumbel law. The key technique used in the proofs is a joint time-frequency Gaussian approximation to a class of diverging dimensional quadratic forms of non-stationary time series, which may have wider applicability in evolutionary power spectrum analysis.

A second main contribution of the paper is that we propose a simulation based bootstrap method to implement simultaneous statistical inferences to a wide range of problems in time-frequency analysis. The motivation of the bootstrap is to alleviate the slow convergence of the maximum deviation to its Gumbel limit. The bootstrap simply generates independent normally distributed pseudo samples of length $N$ and approximate the distribution of the target maximum deviation with that of the normalized empirical maximum deviations of the spectral density estimates from the pseudo samples. The similar idea was used in, for example [WZ07, ZW10], for different problems. The bootstrap is proved to be asymptotically correct and performs reasonably well in the simulations. One important application of the bootstrap is to construct simultaneous confidence regions (SCR) for the evolutionary spectral density, which enables researchers and practitioners to visually evaluate the magnitude and pattern of the evolutionary power spectra with asymptotically accurate statistical guarantee. In particular, the SCR helps one to visually identify which variations in time and/or frequency are genuine and which variations are likely to be produced by random fluctuations. See Section 7.4 for two detailed applications in earthquake and explosion signal processing and finance. On the other hand, the SCR can be applied to a wide range of tests on the structure of the evolutionary spectra or the time series itself. Observe that typically under some specific structural assumptions, the time-varying spectra can be estimated with a faster convergence rate than those estimated by STFT without any prior information. Therefore a generic testing procedure is to estimate the evolutionary spectra under the null hypothesis and check whether the latter estimated spectra can be fully embedded into the SCR. This is a very general procedure and it is asymptotically correct as long as the evolutionary spectra estimated under the null hypothesis converges faster than the SCR. Furthermore, the test achieves asymptotically the power $1$ for local alternatives whose evolutionary spectra deviate from the null hypothesis with a rate larger than the order of the width of the SCR. Specific examples include tests for non-stationary white noise, weak stationarity and time-frequency separability as well as model validation for locally stationary ARMA models and so on. See Section 5.2 for a detailed discussion and Section 7.4 for detailed implementations of the tests in real data.

Finally, we would like to mention that, under the stationarity assumption, the inference of the spectral density is a classic topic in time series analysis. There is a vast literature on the topic and we will only list a very small number of representative works. Early works on this topic include [Par57, WN67, Bri69, And71, Ros84] among others where asymptotic properties of the spectral density estimates were established under various linearity, strong mixing and joint cumulant conditions. For recent developments see [LW10, PP12, WZ18] among others.

The rest of the paper is organized as follows. We first formulate the problem in Section 2. In Section 3, we study the STFT and show that the STFTs are asymptotically independent Gaussian random variables under very mild conditions. In Section 4, we study the asymptotic properties of the STFT-based spectral density estimates, including consistency and asymptotic normality. In Section 5, we establish a maximum deviation theory for the STFT-based spectral density estimates over a nearly optimally dense grid in the joint time-frequency domain. In Section 6, we discuss tuning parameter selection and propose a simulation-based bootstrap method to implement the simultaneous statistical inference. Simulations and real data analysis are given in Section 7. Proofs of the main results are deferred to Section 8 and many details of the proofs have been put in Appendix A.

2. Problem Formulation

We first define locally stationary time series and their instantaneous covariance and spectral density. Throughout the article, we assume the time series $\{X_{i}^{(N)}\}_{i=1}^{N}$ is centered, i.e. $\mathbb{E}[X_{i}^{(N)}]=0$. Furthermore, for a random variable $X$, define $\|X\|_{q}:=[\mathbb{E}|X|^{q}]^{1/q}$ and use $\|\cdot\|$ to denote $\|\cdot\|_{2}$ for simplicity.

Definition 2.1**.**

(Locally stationary time series [ZW09]) We say $\{X_{i}^{(N)}\}_{i=1}^{N}$ is a locally stationary time series if there exists a nonlinear filter $G$ such that

[TABLE]

where $\mathcal{F}_{i}=(\dots,\epsilon_{0},\dots,\epsilon_{i-1},\epsilon_{i})$ and $\epsilon_{i}$’s are i.i.d. random variables. Furthermore, the nonlinear filter $G$ satisfies the stochastic Lipschitz continuity condition, $\operatorname{SLC}(q)$, for some $q>0$; that is, there exists $C>0$ such that for all $i$ and $u,s\in(0,1)$, we have

[TABLE]

Remark 2.2*.*

For time series $X_{1},X_{2},\dots,X_{N}$, we rescale the time index as $t_{i}=i/N$, $i=1,\dots,N$. Then $\{t_{i}\}$ forms a dense grid in $[0,1]$. The rescaled time $u\in[0,1]$ is a natural extension of $\{t_{i}\}_{i=1}^{N}$ to be continuum. This rescaling provides an asymptotic device for studying locally stationary time series, which was first introduced by [Dah97]. In particular, the rescaling together with the stochastic Lipschitz continuity assumption ensure that for each $X_{i}$, there is a diverging number of data points in its neighborhood with similar distributional properties.

Example 2.3**.**

(Locally stationary linear time series) Let $\epsilon_{i}$ be i.i.d. random variables and

[TABLE]

where $a_{j}(u)\in\mathcal{C}^{1}[0,1]$ for $j=0,1,\dots$. This model was considered in [Dah97]. Verification of the SLC assumption is discussed in [ZW09, Propositions 2 and 3].

Example 2.4**.**

(Time varying threshold AR models) Let $\epsilon_{i}\in\mathcal{L}^{q},q>0$ be i.i.d. random variables with distribution function $F_{\epsilon}$ and density $f_{\epsilon}$. Consider the model

[TABLE]

where $a(\cdot),b(\cdot)\in\mathcal{C}^{1}[0,1]$. Then if $\operatorname*{\mathrm{sup}\vphantom{\mathrm{infsup}}}_{u}[|a(u)|+|b(u)|]<1$, the $\operatorname{SLC}(q)$ assumption holds. See also [ZW09, Section 4] for more discussions on checking the SLC assumption for locally stationary nonlinear time series.

For simplicity, we will use $X_{i}$ to denote $X_{i}^{(N)}$ in this paper. Without loss of generality, we assume $X_{i}=0$ for any $i>N$. We adopt the physical dependence measure [ZW09] to describe the dependence structure of the time series.

Definition 2.5**.**

(Physical dependence measure) Let $\{\epsilon_{i}^{\prime}\}$ be an i.i.d. copy of $\{\epsilon_{i}\}$. Consider the locally stationary time series $\{X_{i}\}_{i=1}^{N}$. Assume $\operatorname*{\mathrm{max}\vphantom{\mathrm{infsup}}}_{1\leq i\leq N}\|X_{i}\|_{p}<\infty$. For $k\geq 0$, define the $k$-th physical dependence measure by

[TABLE]

Next, we extend the geometric-moment contraction (GMC) condition [SW07] to the non-stationary setting.

Definition 2.6**.**

(Geometric-moment contraction) We say that the locally stationary time series $\{X_{i}\}_{i=1}^{N}$ is $\operatorname{GMC}(p)$ if for any $k$ we have $\delta_{p}(k)=\mathcal{O}(\rho^{k})$ for some $\rho\in(0,1)$.

Let $\mathcal{P}_{k}(X):=\mathbb{E}(X\,|\,\mathcal{F}_{k})-\mathbb{E}(X\,|\,\mathcal{F}_{k-1})$ and $\tilde{X}_{k}^{[\ell]}:=\mathbb{E}(X_{k}\,|\,\epsilon_{k-\ell+1},\dots,\epsilon_{k})$ be the $\ell$-dependent conditional expectations of $X_{k}$. From the $\operatorname{GMC}(2)$ condition and $\operatorname*{\mathrm{sup}\vphantom{\mathrm{infsup}}}_{k}\|X_{k}\|<\infty$, one can easily verify that $\operatorname*{\mathrm{sup}\vphantom{\mathrm{infsup}}}_{k}\sum_{j=-\infty}^{k}\|\mathcal{P}_{j}X_{k}\|<\infty$ and $\quad\operatorname*{\mathrm{lim}\vphantom{\mathrm{infsup}}}_{\ell\to\infty}\operatorname*{\mathrm{sup}\vphantom{\mathrm{infsup}}}_{k}\|X_{k}-\tilde{X}_{k}^{[\ell]}\|=0$. We refer to Remark A.1 and [SW07] for more discussions on the GMC condition.

Example 2.7**.**

(Non-stationary nonlinear time series) Many stationary nonlinear time series models are of the form

[TABLE]

where $\epsilon_{i}$ are i.i.d. and $R$ is a measurable function. A natural extension to a locally stationary setting is to incorporate the time index $u$ via

[TABLE]

[ZW09, Theorem 6 ] showed that one can have a non-stationary process $X_{i}=X_{i}^{(N)}=G(i/N,\mathcal{F}_{i})$ and the $\operatorname{GMC}(\alpha)$ condition holds, if $\operatorname*{\mathrm{sup}\vphantom{\mathrm{infsup}}}_{u}\|R(u,x_{0},\epsilon_{i})\|_{\alpha}<\infty$ for some $x_{0}$, and

[TABLE]

See [ZW09, Section 4.2] for more details.

Definition 2.8**.**

(Instantaneous covariance) Let $u\in[0,1]$. The instantaneous covariance at $u$ is defined by

[TABLE]

Remark 2.9*.*

The assumption of $\operatorname{SLC}(q)$ together with $\operatorname*{\mathrm{sup}\vphantom{\mathrm{infsup}}}_{i}\mathbb{E}|X_{i}|^{p}<\infty$, where $1/p+1/q=1$, implies the instantaneous covariance $r(u,k)$ is Lipschitz continuous. That is, for all $k$ and for all $u,s\in[0,1],u\neq s$, we have

[TABLE]

for some finite constant $C$. The proof is given in Section A.16. Therefore, uniformly on $u$, for any positive integer $n\leq N$, we have

[TABLE]

Particularly, if we choose $n=o(\sqrt{N})$ then $r(u+\delta_{u},k)-r(u,k)=o(1/n),\forall-n/N\leq\delta_{u}\leq n/N$.

Next, we define the evolutionary spectral density using the instantaneous covariance.

Definition 2.10**.**

(Instantaneous spectral density) Let $u\in[0,1]$. The spectral density at $u$ is defined by

[TABLE]

Remark 2.11*.*

In the definition of instantaneous spectral density, $u\in[0,1]$ represents the rescaled time (see Remark 2.2 for more discussions) and $\theta\in[0,2\pi)$ represents the frequency. Different from the usual spectral density for stationary process, the instantaneous spectral density is a two dimensional function of $u$ and $\theta$, which captures the spectral density variation in both time and frequency. The usual spectral density for stationary process is a one-dimensional function of $\theta$ and is static over time. The notion of instantaneous spectral density is useful for capturing the dynamics of the spectral evolution over time.

Remark 2.12*.*

Note that, for any fixed time point $u$, $r(u,k)$ is a non-negative definite function on the integers. Hence Bochner’s Theorem (or Herglotz Representation Theorem) implies that the covariance function $r(u,k)$ and the spectral density function $f(u,\theta)$ has a one-to-one correspondence at each rescaled time point $u$ under the GMC condition. Therefore, $r(u,k)$ defined in Definition 2.8 has a one-one-one correspondence to the spectral density $f(u,\theta)$ defined in Definition 2.10 for short range dependent locally stationary time series defined in our paper.

In this paper, we always assume $f_{*}:=\operatorname*{\mathrm{inf}\vphantom{\mathrm{infsup}}}_{u,\theta}f(u,\theta)>0$, which is a natural assumption in the time series literature (see e.g. [SW07, LW10]). Finally, we define the STFT, the local periodogram, and the STFT-based spectral density estimates.

Definition 2.13**.**

(Short-time Fourier transform) Let $\tau(\cdot)\leq\tau_{*}<\infty$ be a kernel with support $[-1/2,1/2]$ such that $\tau\in\mathcal{C}^{1}([-1/2,1/2])$ and $\int\tau^{2}(x)\mathrm{d}x=1$. Let $n$ be the number of data in a local window and $\theta\in[0,2\pi)$. Then the STFT is defined by

[TABLE]

Definition 2.14**.**

(Local periodogram)

[TABLE]

Remark 2.15*.*

Note that defining

[TABLE]

then we can write $I_{n}(u,\theta)$ as

[TABLE]

It is well known that $I_{n}(u,\theta)$ is an inconsistent estimator of $f(u,\theta)$ due to the fact that $\hat{r}(u,k)$ are inconsistent when $k$ is large. A natural and classic way to overcome this difficulty is to restrict the above summation to relatively small $k$’s only. This leads to the following.

Definition 2.16**.**

(STFT-based spectral density estimator) Let $a(\cdot)$ be an even, Lipschitz continuous kernel function with support $[-1,1]$ and $a(0)=1$; let $B_{n}$ be a sequence of positive integers with $B_{n}\to\infty$ and $B_{n}/n\to 0$. Then the STFT-based spectral density estimator is defined by

[TABLE]

Remark 2.17*.*

The modified $\hat{f}_{n}(u,\theta)$ in Eq. 17 is not always non-negative as it depends on the property of the kernel function $a(\cdot)$. According to [And91, pp.822], if the kernel function further satisfies $\frac{1}{2\pi}\int_{-\infty}^{\infty}a(x)\exp(-\sqrt{-1}\theta x)\mathrm{d}x\geq 0$ for any $\theta\in[0,2\pi)$, then the modified $\hat{f}_{n}(u,\theta)$ in Eq. 17 is always non-negative. For example, the Bartlett kernel, $a(x)=(1-|x|)\mathbf{1}_{\{|x|\leq 1\}}$, and the Parzen kernel, $a(x)=(1-6x^{2}+6|x|^{3})\mathbf{1}_{\{0\leq|x|\leq 1/2\}}+2(1-|x|)^{3}\mathbf{1}_{\{1/2<|x|\leq 1\}}$.

3. Fourier Transforms

In this section, we study the STFT and show that the STFTs are asymptotically independent and normally distributed under mild conditions. More specifically, when we consider frequencies $\{2\pi j/n:j=1,\dots,n\}$, we show that uniformly over a grid of $u$ and $j$, $\{J_{n}(u,2\pi j/n)\}$ are asymptotically independent and normally distributed random variables.

Denote the real and imaginary parts of $\{J_{n}(u,2\pi j/n)/\sqrt{\pi nf(u,2\pi j/n)}\}$ by

[TABLE]

where $m:=\lfloor(n-1)/2\rfloor$. Then, we have the following result.

Theorem 3.1**.**

Assume $\operatorname{GMC}(2)$, $\operatorname{SLC}(2)$, and $\operatorname*{\mathrm{sup}\vphantom{\mathrm{infsup}}}_{k}\mathbb{E}(X_{k}^{2})<\infty$. Let $\Omega_{p,q}=\{c\in\mathbb{R}^{pq}:|c|=1\}$, where $|\cdot|$ denotes Euclidean norm, and

[TABLE]

for $J=(j_{1},\dots,j_{p})$ satisfies $1\leq j_{1},\dots,j_{p}\leq 2m$ and $U=(u_{1},\dots,u_{q})$ satisfies $0<u_{1}<\dots<u_{q}<1$. Then for any fixed $p,q\in\mathbb{N}$, as $n\to\infty$, we have that

[TABLE]

where $\Phi(x)$ is the cumulative distribution function of the standard normal distribution.

Proof.

See Section 8.1. ∎

The above theorem shows that if we select any $p$ elements from the canonical frequencies $\{2\pi j/n,j=1,\dots,n\}$ and $q$ well-separated points from the re-scaled time, the STFTs are asymptotically independent on the latter time-frequency grid. Moreover, the vector formed by these STFTs is asymptotically jointly normally distributed.

4. Consistency and Asymptotic Normality

In this section, we study the asymptotic properties of the smoothed periodogram estimator $\hat{f}_{n}(u,\theta)$.

4.1. Consistency

The consistency result for the local spectral density estimate $\hat{f}_{n}(u,\theta)$ is as follows.

Theorem 4.1**.**

Assume $\operatorname{GMC}(2)$, $\operatorname{SLC}(2)$, and there exists $\delta\in(0,4]$ such that $\operatorname*{\mathrm{sup}\vphantom{\mathrm{infsup}}}_{i}\mathbb{E}(|X_{i}|^{4+\delta})<\infty$. Let $B_{n}\to\infty$, $B_{n}=\mathcal{O}(n^{\eta})$, $0<\eta<\delta/(4+\delta)$. Then

[TABLE]

Proof.

See Section 8.2. ∎

Later we will see from Theorem 5.3 that the order $\mathcal{O}_{\mathbb{P}}(\sqrt{\log n})$ on the right hand side of Eq. 20 is indeed optimal.

Remark 4.2*.*

Assume $\operatorname*{\mathrm{sup}\vphantom{\mathrm{infsup}}}_{i}\mathbb{E}|X_{i}|^{p}<\infty$ with $p>4$ and $\operatorname{SLC}(q)$ with $1/p+1/q=1$. If we further assume the kernel $\tau(\cdot)$ is an even function and $r(u,k)$ is twice continuously differentiable with respect to $u$, then under $\operatorname{GMC}(2)$, whenever $n=o(N^{2/3})$, $B_{n}=o(\operatorname*{\mathrm{min}\vphantom{\mathrm{infsup}}}\{n,N^{1/3}\})$, and $\operatorname*{\mathrm{sup}\vphantom{\mathrm{infsup}}}_{u}\sum_{k\in\mathbb{Z}}k^{2}|r(u,k)|<\infty$, if $a(\cdot)$ is locally quadratic at [math], i.e.

[TABLE]

where $C$ is a nonzero constant, then we have

[TABLE]

where $f^{\prime\prime}(u,\theta):=-\frac{1}{2\pi}\sum_{k\in\mathbb{Z}}k^{2}r(u,k)\exp(\sqrt{-1}k\theta)$. The proof is given in Section A.13. Therefore, the consistency of $\hat{f}_{n}(u,\theta)$ is implied by combining Theorem 4.1 and Eq. 22.

4.2. Asymptotic Normality

Developing an asymptotic distribution for the local spectral density estimate is an important problem in spectral analysis of non-stationary time series. This allows one to perform statistical inference such as constructing point-wise confidence intervals and performing point-wise hypothesis testing. In the following, we derive a central limit theorem for $\hat{f}_{n}(u,\theta)$.

Theorem 4.3**.**

Assume $\operatorname{GMC}(2)$, $\operatorname{SLC}(2)$, and $\operatorname*{\mathrm{sup}\vphantom{\mathrm{infsup}}}_{i}\mathbb{E}(|X_{i}|^{4+\delta})<\infty$ for some $\delta>0$, $B_{n}\to\infty$ and $B_{n}=o(n/(\log n)^{2+8/\delta})$. Then

[TABLE]

where $\Rightarrow$ denotes weak convergence, $\sigma^{2}_{u}(\theta)=[1+\eta(2\theta)]f^{2}(u,\theta)\int_{-1}^{1}a^{2}(t)\mathrm{d}t$ and $\eta(\theta)=1$ if $\theta=2k\pi$ for some integer $k$ and $\eta(\theta)=0$ otherwise.

Proof.

See Section 8.3. ∎

5. Maximum Deviations

The asymptotic normality for $\hat{f}_{n}(u,\theta)$ derived in the last section cannot be used to construct simultaneous confidence regions (SCR) over $u$ and $\theta$. For simultaneous spectral inference under complex temporal dynamics, one needs to know the asymptotic behavior of the maximum deviation of $\hat{f}_{n}(u,\theta)$ from $f(u,\theta)$ on the joint time-frequency domain, which is an extremely difficult problem. In this section, we establish a maximum deviation theory for the STFT-based spectral density estimates over a dense grid in the joint time-frequency domain. Such results serve as a theoretical foundation for the joint time-frequency inference of the evolutionary spectral densities.

• •

Condition (a): Define $\mathcal{U}=\{u_{1},\dots,u_{C_{n}}\}$ where $C_{n}=|\mathcal{U}|$ and $\frac{n}{2N}<u_{i}<1-\frac{n}{2N},i=1,\dots,C_{n}$. For any $u_{i_{1}},u_{i_{2}}\in\mathcal{U}$ with $i_{1}\neq i_{2}$, we assume that $|u_{i_{1}}-u_{i_{2}}|\geq\frac{n}{N}(1-1/(\log B_{n})^{2})$.

• •

Condition (b): Assume $\operatorname*{\mathrm{sup}\vphantom{\mathrm{infsup}}}_{k}\mathbb{E}|X_{k}|^{p}<\infty$ where $p>4$, and $\operatorname{SLC}(q)$ where $1/p+1/q=1$. Let $\alpha$ be a constant such that $\frac{3}{4(p-1)}<\alpha<\frac{1}{4}$. Then assume $C_{n}=o[\operatorname*{\mathrm{min}\vphantom{\mathrm{infsup}}}\{(nB_{n})^{2\alpha(p-1)-1},B_{n}^{1+2\alpha(p-2)}n^{-2-2\gamma}\}]$ for some $\gamma>0$.

• •

Condition (c): Assume that $a(\cdot)$ is an even and bounded function with bounded support $[-1,1]$, $\operatorname*{\mathrm{lim}\vphantom{\mathrm{infsup}}}_{x\to 0}a(x)=a(0)=1$, $\int_{-1}^{1}a^{2}(x)\mathrm{d}x<\infty$, and $\sum_{j\in\mathbb{Z}}\operatorname*{\mathrm{sup}\vphantom{\mathrm{infsup}}}_{|s-j|\leq 1}|a(jx)-a(sx)|=\mathcal{O}(1)$ as $x\to 0$.

• •

Condition (d): There exists $0<\delta_{1}<\delta_{2}<1$ and $c_{1},c_{2}>0$ such that for all large $n$, $c_{1}n^{\delta_{1}}\leq B_{n}\leq c_{2}n^{\delta_{2}}$.

Note that Conditions (c) and (d) are very mild. Condition (a) implies that the time interval between any two time points on the grid $\mathcal{U}$ cannot be too close. Condition (b) implies that the total number of the selected time points is not too large.

Remark 5.1*.*

Condition (a) implies that $C_{n}\leq\frac{N}{n}(1-\frac{n}{N})(1-\frac{1}{(\log B_{n})^{2}})=\mathcal{O}(N/n)$. Although we do not assume $\{u_{i}\}$ to be equally spaced, we suggest in practice choosing $\{u_{i}\}$ equally spaced and $C_{n}=\frac{N}{n}(1-\frac{n}{N})(1-\frac{1}{(\log B_{n})^{2}})$ to avoid the tricky problem on how to choose the $u_{i}$’s and the $C_{n}$.

Definition 5.2**.**

(Dense Grid $\mathcal{G}_{N}$) Let $\mathcal{G}_{N}$ be a collection of time-frequency pairs such that $(u,\theta)\in\mathcal{G}_{N}$ if $u\in\mathcal{U}$ and $\theta\in\{\frac{i\pi}{B_{n}},i=0,\dots,B_{n}\}$.

The following theorem states that the maximum deviation of the spectral density estimates behaves asymptotically like a Gumbel distribution.

Theorem 5.3**.**

Under $\operatorname{GMC}(2)$ and Conditions (a)–(d), we have that, for any $x\in\mathbb{R}$,

[TABLE]

Proof.

See Section 8.4. ∎

Theorem 5.3 states that the spectral density estimates $\hat{f}_{n}(u,\theta)$ on a dense grid $\mathcal{G}_{N}$ consisting of $C_{n}\times B_{n}$ total number of pairs of $(u,\theta)$ are asymptotically independent quadratic forms of $\{X_{i}\}_{i=1}^{N}$. Furthermore, the maximum deviation of the spectral density estimates on $\mathcal{G}_{N}$ converges to a Gumbel law. This result can be used to construct SCR for the evolutionary spectral densities. Note that Theorem 5.3 is established for a very general class of possibly nonlinear locally stationary processes for the joint and simultaneous time-frequency inference of the evolutionary spectral densities.

5.1. Near optimality of the grid selection

Note that there is a trade-off on how dense the grid should be chosen. On the one hand, we hope the grid is dense enough to asymptotically correctly depict the whole time-frequency stochastic variation of the estimates. On the other hand, making the grid too dense is a waste of computational resources since it does not reveal any extra useful information on the overall variability of the estimates. In the following, we define the notion of asymptotically uniform variation matching of a sequence of dense grids. The purpose of the latter notion is to mathematically determine how dense a sequence of grids should be such that it will adequately capture the overall stochastic variation of the spectral density estimates on the joint time-frequency domain.

Definition 5.4**.**

(Asymptotically uniform variation matching of grids) Consider a given sequence of bandwidths $(n,B_{n})$, and let $\{\tilde{\mathcal{G}}_{N}\}$ be a sequence of grids of time-frequency pairs $\{(u_{i},\theta_{j})\}$ with time and frequencies equally spaced i.e. $|u_{i+1}-u_{i}|=\delta_{\theta,n}$ and $|\theta_{j+1}-\theta_{j}|=\delta_{u,n}$, respectively. Then the sequence $\{\tilde{\mathcal{G}}_{N}\}$ is said to be asymptotically uniform variation matching if

[TABLE]

Note that we have previously shown in Theorem 4.1 that the uniform stochastic variation of $\sqrt{n/B_{n}}\hat{f}_{n}(u,\theta)$ on $(u,\theta)\in(0,1)\times[0,\pi)$ has the order $\mathcal{O}_{\mathbb{P}}(\sqrt{\log n})$. In combination with Theorem 5.3, we can see the order $\mathcal{O}_{\mathbb{P}}(\sqrt{\log n})$ cannot be improved. Therefore, by a simple chaining argument, we can show if a sequence of grids $\{\tilde{\mathcal{G}}_{N}\}$ is an asymptotically uniform variation matching, then

[TABLE]

Hence, the uniform stochastic variation of $\hat{f}_{n}(u,\theta)$ on $(u,\theta)\in\tilde{\mathcal{G}}_{N}$ is asymptotically equal to the uniform stochastic variation of $\hat{f}_{n}(u,\theta)$ on $(u,\theta)\in(0,1)\times[0,\pi)$. In other words, $\operatorname*{\mathrm{max}\vphantom{\mathrm{infsup}}}_{\{u_{i},\theta_{j}\}\in\tilde{\mathcal{G}}_{N}}\left|\hat{f}_{n}(u_{i},\theta_{j})-\mathbb{E}(\hat{f}_{n}(u_{i},\theta_{j}))\right|$ and $\operatorname*{\mathrm{sup}\vphantom{\mathrm{infsup}}}_{(u,\theta)\in(0,1)\times[0,\pi)}\left|\hat{f}_{n}(u,\theta)-\mathbb{E}(\hat{f}_{n}(u,\theta))\right|$ have the same limiting distribution.

However, a grid that is asymptotically uniform variation matching may be unnecessarily dense which causes a waste of computational resources without depicting any additional useful information. The optimal grid should balance between computational burden and asymptotic correctness in depicting the overall time-frequency stochastic variation of the estimates. Furthermore, if the grid is too dense, the limiting distribution is different from our main result and is unknown to the best of our knowledge. Therefore, we hope to choose a sequence of grids as sparse as possible provided it is (nearly) asymptotically uniform variation matching.

Next, we show the sequence of grids used in Theorem 5.3 is indeed nearly optimal in this sense. Recall that in Theorem 5.3, the interval between adjacent frequencies is of order $\delta_{\theta,n}=\Omega(1/B_{n})$ and the averaged interval between two adjacent time indices is of order $\delta_{u,n}=\Omega(n/N)$, where we define $a_{n}=\Omega(b_{n})$ if $1/a_{n}=\mathcal{O}(1/b_{n})$. In the following, we show that if we choose a sequence of slightly denser grids with $\delta_{\theta,n}=\mathcal{O}\left(\frac{1}{B_{n}(\log n)^{\alpha}}\right)$ and $\delta_{u,n}=\mathcal{O}\left(\frac{n}{N(\log n)^{\alpha}}\right)$ where $\alpha$ is any fixed positive constant, then the latter sequence of grids is asymptotically uniform variation matching. Since $\alpha$ can be chosen arbitrarily close to zero, the dense grids in Theorem 5.3 are nearly optimal.

Theorem 5.5**.**

Under the assumptions of Theorem 5.3, a sequence of grids with equally spaced time and frequency intervals $\delta_{u,n}$ and $\delta_{\theta,n}$ is asymptotically uniform variation matching if $\delta_{u,n}=\mathcal{O}\left(\frac{n}{N(\log n)^{\alpha}}\right)$ and $\delta_{\theta,n}=\mathcal{O}\left(\frac{1}{B_{n}(\log n)^{\alpha}}\right)$ for some $\alpha>0$.

Proof.

See Section A.14. ∎

5.2. Applications of the Simultaneous Confidence Regions

In this subsection, we illustrate several applications of the proposed SCR for joint time-frequency inference. These examples include testing time-varying white noise (Example 5.6), testing stationarity (Example 5.7), testing time-frequency separability or correlation stationarity (Example 5.8), and validating time-varying ARMA models (Example 5.9).

These examples demonstrate that our maximum deviation theory can serve as a foundation for the joint and simultaneous time-frequency inference. In particular, as far as we know, there is no existing methodology in the literature for testing time-frequency separability of locally stationary time series, nor model validation for time-varying ARMA models, although they are certainly very important problems. On the other hand, our proposed SCR serves as an asymptotically valid and visually friendly tool for the above purposes (see Examples 5.8 and 5.9).

In order to implement the tests, observe that typically under some specific structural assumptions, the time-varying spectra can be estimated with a faster convergence rate than those estimated by the STFT. Therefore, to test the structure of the evolutionary spectra under the null hypothesis, a generic procedure is to check whether the estimated spectra under the null hypothesis can be fully embedded into the SCR. Note that this very general procedure is asymptotically correct as long as the evolutionary spectra estimated under the null hypothesis converges faster than the SCR. The test achieves asymptotic power $1$ for local alternatives whose evolutionary spectra deviate from the null hypothesis with a rate larger than the order of the width of the SCR.

Example 5.6**.**

(Testing time-varying white noise) White noise is a collection of uncorrelated random variables with mean [math] and time-varying variance $\sigma^{2}(u)$. It can be verified that testing time-varying white noise is equivalent to testing the following null hypothesis:

[TABLE]

for some time-varying function $g(\cdot)$. Consider the following optimization problem:

[TABLE]

That is, we would like to find a function of $u$ which is closest to $f(u,\theta)$ in $L_{2}$ distance. Direct calculations show that $g_{0}(u)=\frac{1}{\pi}\int_{0}^{\pi}f(u,\theta)\mathrm{d}\theta$. Therefore, under the null hypothesis we can estimate the function $g$ in Eq. 27 by

[TABLE]

It can be shown that under the null hypothesis the convergence rate of $\hat{g}(u)$ uniformly over $u$ is $\mathcal{O}_{\mathbb{P}}(\sqrt{\log n}/\sqrt{n})$, which is faster than the rate of SCR which is $\mathcal{O}_{\mathbb{P}}(\sqrt{\log n}/\sqrt{n/B_{n}})$. Therefore, we can apply the proposed SCR to test time-varying white noise.

Example 5.7**.**

(Testing stationarity) Under the null hypothesis that the time series is stationary, it is equivalent to testing

[TABLE]

for some function $h(\cdot)$. Consider the following optimization problem:

[TABLE]

That is, we would like to find a function of $\theta$ which is closest to $f(u,\theta)$ in $L_{2}$ distance. Direct calculations show that $h_{0}(\theta)=\int_{0}^{1}f(u,\theta)\mathrm{d}u$. Therefore, under the null hypothesis, we can estimate the function $h$ in Eq. 30 by

[TABLE]

It can be shown that the convergence rate of $\hat{h}(\theta)$ uniformly over $\theta$ is $\mathcal{O}_{\mathbb{P}}(\sqrt{\log n}/\sqrt{N/B_{n}})$, which is faster than the rate $\mathcal{O}_{\mathbb{P}}(\sqrt{\log n}/\sqrt{n/B_{n}})$ of the SCR. Therefore, we can apply the proposed SCR to test stationarity.

Example 5.8**.**

(Testing time-frequency separability or correlation stationarity) We call a non-stationary time series time-frequency separable if $f(u,\theta)=g(u)h(\theta)$ for some functions $g(\cdot)$ and $h(\cdot)$. If a non-stationary time series is time-frequency separable, the frequency curves across different times are parallel to each other. Similarly, the time curves across different frequencies are parallel to each other as well. Therefore, the property of time-frequency separability enables one to model the temporal and spectral behaviors of the time-frequency function separately. Furthermore, it can be verified that testing time-frequency separability is equivalent to testing correlation stationarity for locally stationary time series, i.e. $\textrm{corr}(X_{i},X_{i+k})=l(k)$, for some function $l(\cdot)$. Without loss of generality, we can formulate the null hypothesis as

[TABLE]

for some constant $C_{0}$ and $\int_{0}^{1}g(u)\mathrm{d}u=1$ and $\int_{0}^{\pi}h(\theta)=1$. Under the null hypothesis, we can estimate $C_{0}$, $g(u)$ and $h(\theta)$ by

[TABLE]

and we can estimate $f(u,\theta)$ by $\hat{C}_{0}\hat{g}(u)\hat{h}(\theta)$. It can be shown that the convergence rates of $\hat{C}_{0}$, $\hat{g}(u)$, and $\hat{h}(\theta)$ are $\mathcal{O}_{\mathbb{P}}(1/\sqrt{N})$, $\mathcal{O}_{\mathbb{P}}(\sqrt{\log n}/\sqrt{n})$, and $\mathcal{O}_{\mathbb{P}}(\sqrt{\log n}/\sqrt{N/B_{n}})$, respectively. All of them are faster than the convergence rate of the SCR which is $\mathcal{O}_{\mathbb{P}}(\sqrt{\log n}/\sqrt{n/B_{n}})$. Therefore, we can apply the proposed SCR to test the null hypothesis.

Example 5.9**.**

(Validating time-varying ARMA models) Consider the null hypothesis that the time series follows the following time-varying ARMA model

[TABLE]

where $a_{0}(u)=1$, $a_{i}(\cdot),b_{i}(\cdot)\in\mathcal{C}^{1}[0,1]$, and $\epsilon_{i}$ are uncorrelated random variables with mean [math] and variance $1$. Under the null hypothesis, $\{X_{i}\}$ is a locally stationary time series with spectral density

[TABLE]

The spectral density can be fitted using the generalized Whittle’s method [Dah97], where $a_{i}(t/N)$ and $b_{i}(t/N)$ are estimated by minimizing a generalized Whittle function and $p$ and $q$ are selected, for example, by AIC. Note that under the null hypothesis, the spectral density estimated using Whittle’s method has a convergence rate $\mathcal{O}_{\mathbb{P}}(\sqrt{\log n}/\sqrt{n})$ which is faster than the rate $\mathcal{O}_{\mathbb{P}}(\sqrt{\log n}/\sqrt{n/B_{n}})$ by the STFT-based methods without prior information. Therefore, to test the fitted non-parametric time-varying ARMA model, we can plot the non-parametric spectral density using the estimated time-varying parameters $a_{i}(\cdot)$ and $b_{i}(\cdot)$. Under the null hypothesis, the non-parametric spectral density should fall within our SCR with the prescribed probability asymptotically.

The benefits of spectral domain approach to various hypothesis testing problems depend on the specific application. For example, for tests of stationarity, the test based on evolutionary spectral density is technically easier than the corresponding tests in the time domain. The main reason is that the time domain test needs to consider time-invariance of $r(u,k)$ for a diverging number of $k$ and hence is a high-dimensional problem. On the other hand, the spectral domain test of stationarity only needs to check that $f(u,\theta)$ does not depend on $u$. Similar arguments apply to the test of white noise. For another example, we proposed a frequency domain method for the problem of model validation of non-stationary linear models. However, technically it is difficult to approach this problem from the time domain. Furthermore, for many time series signals in engineering applications, the most important information is embedded in the frequency domain. Therefore, in engineering and signal processing applications, frequency domain methods are typically more favourable and are widely used. Therefore frequency-domain-based tests are preferable in many such applications.

6. Bootstrap and Tuning Parameter Selection

In Section 6.1, we propose a simulation based bootstrap method to implement simultaneous statistical inferences. The motivation of the bootstrap procedure is to alleviate the slow convergence of the maximum deviation to its Gumbel limit in Theorem 5.3. We discuss methods for tuning parameter selection in Section 6.2.

6.1. The Bootstrap Procedure

Although Theorem 5.3 shows that SCR can be constructed using the Gumbel distribution, the convergence rate in Theorem 5.3 is too slow to be useful in moderate samples. We propose a bootstrap procedure to alleviate the slow convergence of the maximum deviations. One important application of the bootstrap is to construct SCR in moderate sample cases.

Let $\{\epsilon_{1},\dots,\epsilon_{N}\}$ be i.i.d. $\mathcal{N}(0,1)$ random variables. Defining

[TABLE]

and

[TABLE]

it can be easily verified that the following analogy of Theorem 5.3 holds.

[TABLE]

Therefore, we propose to construct the SCR for $\{\hat{f}_{n}(u,\theta)\}$ using the empirical distribution of $\hat{f}^{\epsilon}_{n}(u,\theta)$. More specifically, we generate $\{\epsilon_{i}\}_{i=1}^{N}$ independently for $N_{\textrm{MC}}$ times. Let $\bar{f}^{\epsilon}_{n}(u,\theta)$ be the sample mean of $\{\hat{f}^{\epsilon}_{n,m}(u,\theta),m=1,\dots,N_{\textrm{MC}}\}$ from the $N_{\textrm{MC}}$ Monte Carlo experiments. Then we compute the empirical distribution of

[TABLE]

to approximate the distribution of

[TABLE]

which can be employed to construct the SCR. For example, for a given $\alpha\in(0,1)$, we estimate the $(1-\alpha)$-th quantile $\gamma^{2}_{1-\alpha}$ from the bootstrapped distribution using $\hat{f}^{\epsilon}_{n}(u,\theta)$, which also approximately satisfies

[TABLE]

Therefore, the constructed SCR is

[TABLE]

Note that in small sample cases, the lower bound for the confidence region can be [math] if the estimated $\gamma_{1-\alpha}$ is larger than $1$. This happens when $N$ is not large enough and large $B_{n}$ and $C_{n}$ are selected. For large sample sizes, the estimated $\gamma_{1-\alpha}^{2}$ is typically much smaller than $1$. In that case, we can further use the following approximation

[TABLE]

Then the SCR can be constructed as

[TABLE]

Overall, the practical implementation is given as follows

1. (1)

Select $B_{n}$ and $n$ using the tuning parameter selection method described in Section 6.2; 2. (2)

Compute the critical value using bootstrap described in Section 6.1; 3. (3)

Compute the spectral density estimates by Eq. 17; 4. (4)

Compute the SCR defined in Section 6.1 using the spectral density estimates and the critical value obtained by the bootstrap.

Note that the validity of the proposed bootstrap procedure is asymptotically justified by Theorem 5.3. On the other hand, theoretical justification for the superiority of the bootstrap procedure for moderate samples is extremely difficult, as it requires deriving higher order asymptotics of the maximum deviation of the time-varying spectral densities. We will investigate this problem in some future work.

6.2. Tuning parameter selection

Choosing $B_{n}$ and $n$ in practice is a non-trivial problem. In our Monte Carlo experiments and real data analysis, we find that the minimum volatility (MV) method [PRW99, Zho13] performs reasonably well. Specifically, the MV method uses the fact that the estimator $\hat{f}_{n}(u,\theta)$ becomes stable when the block size $n$ and the bandwidth $B_{n}$ are in an appropriate range. More specifically, we first set a proper interval for $n$ as $[n_{l},n_{r}]$. In our simulations and data analysis, we choose $n_{l}=2N^{\eta}$ and $n_{r}=3N^{\eta}$ if $N\leq 1000$, $n_{l}=2.5N^{\eta}$ and $n_{r}=4N^{\eta}$ if $1000<N\leq 2000$ and $n_{l}=3N^{\eta}$ and $n_{r}=5N^{\eta}$ if $N>2000$, where $\eta=0.48$. Although the rule for setting $n_{l}$ and $n_{r}$ is ad-hoc, it works well in our simulations and data analysis. In practice, one can also either choose $n_{l}$ and $n_{r}$ based on prior knowledge of the data, or select them by visually evaluating the fitted evolutionary spectral densities. A reasonable value of $n$ should not produce too rough or too smooth estimates of the spectral density. We remark that $n_{l}$ and $n_{u}$ are only upper and lower bounds of the candidate bandwidths. Simulations show that the simulated coverage probabilities are typically not sensitive to the choices of $n_{l}$ and $n_{r}$.

In order to use the MV method, we first form a two-dimensional grid of all candidate pairs of $(n,B_{n})$ such that $n\in[n_{l},n_{r}]$ and $B_{n}<n/\log(n)$. Then, for each candidate pair $(n,B_{n})$, we estimate $\hat{f}_{n}(u,\theta)$ using the candidate pair for a fixed time-frequency grid of $(u,\theta)$. Next, we compute the average variance of the spectral density estimates $\hat{f}_{n}(u,\theta)$ over the neighborhood of each candidate pair on the two-dimensional grid of all candidate pairs of $(n,B_{n})$. Finally, we choose the pair of $(n,B_{n})$ which gives the lowest average variance. We refer to [PRW99, Zho13] for more detailed discussions of the MV method.

Note that cross validation is another popular method for choosing bandwidths [DR19]. However, it is a difficult task to implement cross validation in the context of time-varying spectral density estimation. Finally, it is well-known that choosing theoretically optimal bandwidths is an extremely difficult problem. We hope to investigate this problem in some future work.

7. Simulations and Data Analysis

In this section, we study the performance of the proposed SCR via simulations and real data analysis. In Section 7.1, the accuracy of the proposed bootstrap procedure is studied; The accuracy of tuning parameter selection is considered in Section 7.2; The accuracy and power for hypothesis testing is studied in Section 7.3; Finally, we perform real data analysis in Section 7.4. Throughout this section, the kernel $\tau(\cdot)$ is chosen to be a re-scaled Epanechnikov kernel such that $\int\tau^{2}(x)\mathrm{d}x=1$, and the kernel $a(\cdot)$ is a re-scaled tri-cube kernel such that $a(0)=1$. The two kernel functions are defined as follows.

[TABLE]

In all the simulations, we ran Monte Carlo experiments for $N_{\textrm{MC}}=10000$. The results for SCR and hypothesis testing are obtained by averaging over $1000$ independent datasets.

7.1. Accuracy of Bootstrap

In this subsection, we study the accuracy of the proposed bootstrap procedure for moderate finite samples (e.g. $N=400$ or $N=800$). We consider different examples of locally stationary time series models described in the following Examples 7.1, 7.2, 7.3, 7.4 and 7.5.

Example 7.1**.**

(Time-varying AR model) We have

[TABLE]

where $\{\epsilon_{i}\}$ are i.i.d. $\mathcal{N}(0,1)$. In this example, we choose $a(u)=0.3\cos(2\pi u)$. Then the model is locally stationary in the sense that the AR$(1)$ coefficient $a(u)=0.3\cos(2\pi u)$ changes smoothly on the interval $[0,1]$. The simulated uncoverage probabilities of the SCR are shown in Table 1.

Example 7.2**.**

(Time-varying ARCH model) Consider the following time-varying ARCH$(1)$ model:

[TABLE]

where $\{\epsilon_{i}\}$ are i.i.d. standard normally distributed random variables, $a_{0}(u)>0,a_{1}(u)>0$ and $a_{0}(u)+a_{1}(u)<1$. Note that $\{X_{i}\}$ is a white noise sequence. In this example, we choose $a_{0}(u)=0.6$ and $a_{1}(u)=0.3\sin(\pi u)$. The simulated uncoverage probabilities of the SCR are shown in Table 2.

Example 7.3**.**

(Time-varying Markov switching model) Suppose $\{S_{i}\}$ is a Markov chain on state space $\{0,1\}$ with transition matrix $P$. Consider the following time-varying Markov switching model

[TABLE]

where $\{\epsilon_{i}\}$ are i.i.d. standard normally distributed random variables, $|b_{1}|<1$, and $|b_{2}|<1$. In this example, we choose $P=\begin{bmatrix}0.9&0.1\\ 0.5&0.5\end{bmatrix}$, $b_{1}(u)=0.4\cos(2\pi u)$, and $b_{2}(u)=0.1\sin(2\pi u)$. The simulated uncoverage probabilities of the SCR are shown in Table 3.

Example 7.4**.**

(Time-varying threshold AR model) Suppose $\{\epsilon_{i}\}$ are i.i.d. standard normally distributed random variables and consider the following threshold AR model

[TABLE]

where $\operatorname*{\mathrm{sup}\vphantom{\mathrm{infsup}}}_{u\in[0,1]}[|a(u)|+|b(u)|]<1$. In this example, we choose $a(u)=0.3\cos(2\pi u)$ and $b(u)=0.3\sin(2\pi u)$. The simulated uncoverage probabilities of the SCR are shown in Table 4.

Example 7.5**.**

(Time-varying bilinear process) Let $\{\epsilon_{i}\}$ be i.i.d. standard normally distributed random variables and consider the following model

[TABLE]

where $b^{2}(u)+c^{2}(u)<1$. In this example, we choose $b(u)=0.3\cos(2\pi u)$ and $c(u)=0.1\sin(2\pi u)$. The simulated uncoverage probabilities of the SCR are shown in Table 5.

According to the results in Tables 1, 2, 3, 4 and 5, one can see that the proposed bootstrap works well when $B_{n}$ and $n$ are chosen in a relatively wide range. In the next subsection, we discuss the MV method for selecting $B_{n}$ and $n$ in practice.

7.2. Accuracy of Tuning Parameter Selection

We apply the MV method described in Section 6.2 to select the tuning parameters for Examples 7.1, 7.2, 7.3, 7.4 and 7.5. For all examples, $N=400$ and $N=800$ are considered. The bootstrap accuracy is shown in Table 6. Furthermore, according to Eq. 45, we also included the average width of the SCR over $(u,\theta)\in\mathcal{G}_{N}$ in Table 6. From Table 6, we can see that the coverage probabilities of the SCR with bandwidths selected by the MV method are accurate. Furthermore, the average width of the SCR decreases as $N$ increases.

7.3. Accuracy and Power of Hypothesis Testing

In this subsection, we study the accuracy and power of hypothesis testing using the proposed SCR. We consider Example 7.6 for testing stationarity and Example 7.7 for testing time-varying white noise. Furthermore, we also consider another example of non-parametric ARMA model validation, which is given in Example 7.8.

Example 7.6**.**

(Time-varying ARCH model) Consider the following model

[TABLE]

where $a_{0}(u)=0.3$ and $a_{1}(u)=0.2+\delta u$. Observe that when $\delta=0$, the model is stationary. When $\delta=0$, the accuracy of the hypothesis testing for stationarity is studied for two cases, one with $N=400$ and the other with $N=800$, where $n$ and $B_{n}$ are selected by the MV method. We have shown the simulated Type I error rates of the SCR in Table 7. Next, we study the power of the hypothesis testing for stationarity using the proposed SCR by increasing $\delta$. We study both $0.05$ and $0.1$ level tests. The simulated powers of the SCR for $N=800$ and $N=600$ are shown in Fig. 2. One can see that, for both $N=800$ and $N=600$, the simulated Type I error rates of the SCR (when $\delta=0$) are accurate. Furthermore, the simulated power of the SCR increases with $N$.

Example 7.7**.**

(Time-varying MA model) Consider the following model:

[TABLE]

where we let $a_{0}(u)=0.7+0.9\cos(2\pi u)$ and $a_{1}(u)=\delta a_{0}(u)$. Clearly, when $\delta=0$, the model generates a time-varying white noise. When $\delta=0$, we study the accuracy of the hypothesis testing for time-varying white noise using the proposed SCR. The accuracy by the SCR is shown in Table 7, one with $N=800$ and the other with $N=1200$. The tuning parameters $n$ and $B_{n}$ are selected by the MV method. We then test time-varying white noise using our proposed SCR by increasing $\delta$ for $N=800$ and $N=600$. The simulated powers of the SCR are shown in Fig. 2. According to Fig. 2, one can see that the simulated coverage probabilities for $N=600$ are slightly below the nominal level. This is because the structure of the time series in this example is complicated. A sample size with $N=600$ is not large enough for the local stationarity of the time series to be fully captured statistically. On the other hand, for sample size $N=800$, the simulated Type I error rates are accurate and the powers are significantly higher than the case of $N=600$.

Example 7.8**.**

(Validating time-varying AR model) Consider the following time-varying AR model

[TABLE]

where $a_{0}(u)=1$, $a_{j}(\cdot)$ and $\sigma(\cdot)$ are smooth functions, $\epsilon_{i}$ are i.i.d. with mean [math] and variance $1$. Then $\{X_{i}\}$ is a locally stationary time series with spectral density

[TABLE]

In this example, we generate time series with $p=1$, $a_{1}(u)=0.3+0.2u$, $\sigma(u)=1+0.3u+0.2u^{2}$, and length $N=400$ or $N=800$.

For each generated time series, we fit a time-varying AR model with $p=1$ by minimizing the local Whittle likelihood [Dah97]. We can then test if the spectral density of the fitted non-parametric time-varying AR model falls into the proposed SCR. The simulated coverage probabilities of the SCR are shown in Table 7, where $n$ and $B_{n}$ are selected by the MV method. We can see that, under the null hypothesis, the non-parametric time-varying AR model is validated since the simulated Type I error rates match quite well with the nominal levels of the proposed SCR test.

7.4. Real Data Analysis

In this subsection, we present some real data analysis . We study an earthquake and explosion data set from seismology in Example 7.9 and then daily SP500 return from finance in Example 7.10. Observe that all time series are relatively long with $N>2000$. For tuning parameter selection, we use the MV method to search $(n,B_{n})$ within the region $B_{n}<n/\log(n)$. Hypothesis tests are performed, including testing stationarity, time-varying white noise, and time-frequency separability on all the data sets.

Example 7.9**.**

(Earthquakes and explosions [SS17]) In this example, we study an earthquake signal and an explosion signal from a seismic recording station [SS17]. The recording instruments in Scandinavia are observing earthquakes and mining explosions with one of each shown in Fig. 4 and Fig. 4, respectively. The two time series (see Fig. 4 and Fig. 4) each has length $N=2048$ representing two phases or arrivals along the surface, denote by phase $P$: $\{X_{i}:i=1,\dots,1024\}$ and phase $S$: $\{X_{i}:i=1025,\dots,2048\}$. The general problem of interest is in distinguishing or discriminating between waveforms generated by earthquakes and those generated by explosions. The original data came from the technical report by [Bla93]. According to [Bla93], the original earthquake and explosion signals have been filtered with a $3$-pole highpass Butterworth filter with the corner frequency at $1$ Hz to improve the signal-to-noise ratio. Then the amplitudes of the waveforms have been rescaled so the maximum amplitude for each signal is equal. According to [Bla93, Figure 2a and 2b], the unit for time is $0.02$ second and the values of the earthquake and explosion data are rescaled to be no more than $1$.

From the time domain (see Figs. 4 and 4), one can observe that rough amplitude ratios of the first phase $P$ to the second phase $S$ are different for the two data sets, which tend to be smaller for earthquakes than for explosions. From the spectral density estimates and their confidence regions, the $S$ component for the earthquake (see Fig. 4) shows power at the low frequencies only, and the power remains strong for a long time. In contrast, the explosion (see Fig. 4) shows power at higher frequencies than the earthquake, and the power of the $P$ and $S$ waves does not last as long as in the case of the earthquake.

Moreover, we notice from the confidence region at selected times and frequencies that the spectral density of explosion has the similar shape at different times, as well as at different frequencies (see Figs. 8 and 8). However, the spectral density of earthquakes does not seem to have this property (see Figs. 6 and 6). This may suggest that the explosion data are correlation stationary or time-frequency separable. We further perform hypothesis tests on both data sets to confirm our observation (see Table 8). The p-values for testing stationarity and time-varying white noise for both earthquake and explosion are quite small, which implies that earthquake and explosion time series are not stationary and not time-varying white noise. However, the p-values for the hypothesis of time-frequency separability (i.e., correlation stationary) is $0.61$ for explosion, but $0.064$ for earthquake. This interesting result discovers a potential important difference between earthquake and explosion: at least from the analyzed data, explosion tends to be time-frequency separable (correlation stationary) but earthquake does not.

There are two main benefits from knowing that explosion time series are time-frequency separable but earthquake time series are not. First, this reveals an important structural property of the time-frequency behavior for explosion signals. Since time-frequency separability implies the time curves for different frequencies are parallel and the frequency curves for different times are parallel as well, this directly suggests a parsimonious model for explosion time series using two one-dimensional models. Second, for the classification of earthquake and explosion signals, time-frequency separability provides a non-linear feature of the explosion that could potentially serve the purpose. Since most commonly used features for classification are linear features, time-frequency separability is potentially important for feature extraction in order to improve the accuracy in classification tasks. However, since we only have analyzed one pair of earthquake and explosion signals, further studies with a large database of earthquake and explosion signals are needed to confirm this property for explosions which we leave to a future work.

Example 7.10**.**

(SP500 daily returns)

In this example, we analyze daily returns of SP500 from September 23rd, 1991 to August 17th, 2018. We plot the original time series, the spectral density estimates and their confidence regions in Fig. 10. Observing that the SCR in Fig. 10 appears to be quite flat over frequencies, it is reasonable to ask if the time series may be modeled as time-varying white noise. Actually, in the finance literature, it is commonly believed that stock daily returns behave like time-varying white noise. We further confirm this observation by performing hypothesis tests. The results (see Table 8) show that the SP500 time series is not stationary but it is likely to be a time-varying white noise since the p-value for testing time-varying white noise is $0.99$. Furthermore, the p-value for testing time-frequency separability is also quite large which is $0.99$.

Next, we turn our focus to the absolute value of SP500 daily returns. Volatility forecasting, i.e. forecasting future absolute values or squared values of the return, is a key problem in finance. The celebrated ARCH/GARCH models are equivalent to exponential smoothings of the absolute or squared returns. The optimal weights in the smoothing are determined fully by the evolutionary spectral density. Hence, to optimally forecast the evolutionary volatility, one way is to fit the absolute returns by an appropriate non-stationary linear model, then apply the fitted model to forecast the future volatility. To date, to our knowledge, there exists no methodology for validating non-stationary linear models. In the following, we demonstrate that the proposed SCR is a useful tool for validating non-stationary linear models for absolute SP500 daily returns.

We first remove the local mean of the original SP500 time series by kernel smoothing. The spectral density estimates and the SCRs are shown in Figs. 10, 12 and 12. We observe from the plots that the spectral density of the absolute SP500 returns behaves quite differently from the original SP500 time series. For example, unlike the case for the original SP500 time series, the SCR for the absolute SP500 in Fig. 12 is not flat over frequencies anymore. We perform the same hypothesis tests again to the absolute SP500 time series. The results (see Table 8) show that the p-value for testing time-varying white noise is $0.037$, which is much smaller than that of the original SP500 time series. Furthermore, the p-value for testing time-frequency separability is $0.048$ which is also much smaller than the one for the original SP500 data.

Finally, we fit time-varying non-stationary linear models for the absolute SP500 daily returns with mean removed by kernel smoothing. We first fit various time-varying AR or ARMA models

[TABLE]

to the absolute returns by minimizing the local Whittle likelihood [Dah97]. We then validate if the fitted spectral densities from the time-varying AR or ARMA models fall into the proposed SCR. The p-values for validating time-varying AR/ARMA models are shown in Table 9. One can see that, the p-values for the tv-AR models are quite small, which implies that no tv-AR models up to order $5$ is appropriate for fitting absolute SP500 daily returns. For tv-ARMA models, the p-value for the tv-ARMA$(1,1)$ model equals $0.019$. This suggests that this tv-ARMA model is not appropriate for fitting the absolute SP500 daily returns either. In contrast, the corresponding p-value for validating the tv-ARMA$(2,1)$ model is $0.79$. This interesting observation suggests that the tv-ARMA$(2,1)$ model may be appropriate to fit the absolute returns. We further plot the spectral densities of the fitted time-varying AR$(1)$, AR$(4)$, AR$(5)$, ARMA$(1,1)$, ARMA$(2,1)$, and ARMA$(3,1)$ models in Fig. 13. From Fig. 13, one can see that the fitted spectral densities by the tv-AR models are quite different from the STFT-based spectral density estimates. For tv-ARMA models, the spectral density estimates by the tv-ARMA$(1,1)$ model are not close to the STFT-based spectral density estimates either. Therefore, based on the proposed SCR, we conclude that the tv-ARMA$(2,1)$ model is an appropriate candidate for the analyzed data and can be used for short-term future volatility forecasting.

8. Proofs of Main Results

8.1. Proof of Theorem 3.1

We prove Theorem 3.1 in two steps. In the first step, we show in Section 8.1.1 that Theorem 3.1 is true for $q=1$. In this case, we let $\Omega_{p}=\{c\in\mathbb{R}^{p}:|c|=1\}$, $Z_{u,J}=(Z_{u,j_{1}}^{(n)},\dots,Z_{u,j_{p}}^{(n)})^{T}$ for $J=(j_{1},\dots,j_{p})$ satisfies $1\leq j_{1},\dots,j_{p}\leq 2m$ (recall that $m=\lfloor(n-1)/2\rfloor$). We prove for any fixed $p\in\mathbb{N}$, as $n\to\infty$, we have that

[TABLE]

In the second step of the proof, we show in Section 8.1.2 that for fixed $q\in\mathbb{N}$, for any given $0<u_{1}<\dots<u_{q}<1$, we have $\{(c^{(i)})^{T}Z_{u_{i},J},i=1,\dots,q\}$ are asymptotically independent uniformly over $\{c^{(i)}\in\mathbb{R}^{p}:|c^{(i)}|=1\}$ for $i=1,\dots,q$. Finally, Theorem 3.1 is proved by combining the two parts.

8.1.1. Proof of Eq. 59

We denote $2\pi j/n$ by $\theta_{j}$ in this proof. With out loss of generality, we restrict $J=\{j_{1},\dots,j_{p}\}\in\{1,\dots,m\}$. Let $c=(c_{1},\dots,c_{p})$, define $\mu_{u,k}:=\sum_{\ell=1}^{p}\frac{c_{\ell}\cos(k\theta_{j_{\ell}})}{\sqrt{\pi f(u,\theta_{j_{l}})}}$. Then $\mu_{u,k}\leq\sum_{\ell=1}^{p}\frac{|c_{\ell}|}{\sqrt{\pi f_{*}}}\leq\frac{p}{\sqrt{\pi f_{*}}}=:\mu_{*},\quad\forall c\in\Omega_{p},\quad\forall J$. Furthermore, defining

[TABLE]

and $\eta:=\left(\frac{\|T_{u,n}-\tilde{T}_{u,n}\|}{\sqrt{n}}\right)^{1/2}$, we have the following key lemmas.

Lemma 8.1**.**

Under the assumptions of Theorem 3.1, we have

[TABLE]

Proof.

See Section A.1. ∎

Lemma 8.2**.**

Under the assumptions of Theorem 3.1, we have

[TABLE]

Proof.

See Section A.2. ∎

Lemma 8.3**.**

Under the assumptions of Theorem 3.1, we have

[TABLE]

where $\delta_{n}\to 0$ as $n\to\infty$ uniformly over $J$, $c$ and $u$.

Proof.

See Section A.3. ∎

Using the above results, we can then prove Eq. 59 as follows. First, by Chebyshev inequality and $\eta=\left(\frac{\|T_{u,n}-\tilde{T}_{u,n}\|}{\sqrt{n}}\right)^{1/2}$, we have

[TABLE]

Next, according to Lemma 8.1, uniformly over $J$, $c$ and $u$, for any fixed $\ell$, as $n\to\infty$, we have that

[TABLE]

By Lemma 8.3, we have that

[TABLE]

Note that $\delta_{n}\to 0$ as $n\to\infty$. Also, by Lemma 8.2, uniformly over $J,c,u,n$, we have $\eta\to 0$ as $\ell\to\infty$. Finally, letting $n\to\infty$ then $\ell\to\infty$, we have that $\operatorname*{\mathrm{sup}\vphantom{\mathrm{infsup}}}_{x}\left|\mathbb{P}\left(\frac{T_{u,n}}{\sqrt{n}}\leq x\right)-\Phi\left(x\right)\right|\to 0$, uniformly over $J$, $c$, and $u$.

8.1.2. Proof of asymptotically independence of $\{(c^{(i)})^{T}Z_{u_{i},J},i=1,\dots,q\}$

We can write $T_{u_{i},n}$ and $\tilde{T}_{u_{i},n}$ defined in Eq. 60 as $T_{u_{i},n,c^{(i)}}$ and $\tilde{T}_{u_{i},n,c^{(i)}}$. Then by Lemma 8.2, it suffices to show that $\{\tilde{T}_{u_{i},n,c^{(i)}},i=1,\dots,q\}$ are asymptotically independent uniformly over $\{c^{(i)}\in\mathbb{R}^{p}:|c^{(i)}|=1\}$. Note that in the definition of $\tilde{T}_{u_{i},n,c^{(i)}}$, $\tilde{X}_{k}$ is $\ell$-dependent, therefore, $\tilde{T}_{u_{1},n,c^{(i)}}$ and $\tilde{T}_{u_{2},n,c^{(i)}}$ with $u_{2}>u_{1}$ are independent if $\lfloor(u_{2}-u_{1})N\rfloor>\ell+2n$. Since $0<u_{1}<\dots<u_{q}<1$ are fixed, $\operatorname*{\mathrm{min}\vphantom{\mathrm{infsup}}}_{i\neq j}|u_{i}-u_{j}|>0$ is bounded away from zero. Therefore, $\{\tilde{T}_{u_{i},n,c^{(i)}},i=1,\dots,q\}$ are independent if $\ell<\lfloor(\operatorname*{\mathrm{min}\vphantom{\mathrm{infsup}}}_{i\neq j}|u_{i}-u_{j}|)N\rfloor-2n$. Choosing $\ell=o(n)$ and $n=o(N)$, we have $\{\tilde{T}_{u_{i},n,c^{(i)}},i=1,\dots,q\}$ are asymptotically independent.

8.2. Proof of Theorem 4.1

Throughout the proof, we use $\|\cdot\|$ to denote $\|\cdot\|_{2}$ for simplicity. We define $X_{u,i,n}:=\tau\left(\frac{i-\lfloor n/2\rfloor}{n}\right)X_{\lfloor uN\rfloor+i-\lfloor n/2\rfloor}$. For simplicity we will omit the index $n$ and use $X_{u,i}$ for $X_{u,i,n}$. Define $Y_{u,i}:=Y_{u,i}(\theta)=\frac{1}{2\pi}\sum_{k=-B_{n}}^{B_{n}}X_{u,i}X_{u,i+k}a(k/B_{n})\cos(k\theta)$, $g_{n}(u,\theta):=\sum_{i=1}^{n}Y_{u,i}(\theta)$, and $h_{n}(u,\theta):=\frac{1}{\sqrt{nB_{n}}}g_{n}(u,\theta)-\sqrt{n/B_{n}}\hat{f}_{n}(u,\theta)$, we have that

[TABLE]

Next, denote $\tilde{X}_{k}^{[\ell]}$ as the $\ell$-dependent conditional expectation of $X_{k}$, $\tilde{X}_{u,i}^{[\ell]}$ as the $\ell$-dependent conditional expectation of $X_{u,i}$, and $\tilde{Y}_{u,i}$ as the correspondence of sum using $\tilde{X}_{u,i}^{[\ell]}$ instead of $X_{u,i}$, and $\tilde{g}_{n}$ as the correspondence of $g_{n}$ using $\tilde{Y}_{u,i}$ instead of $Y_{u,i}$. Note that under $\operatorname{GMC}(2)$ and $\operatorname*{\mathrm{sup}\vphantom{\mathrm{infsup}}}_{i}\mathbb{E}|X_{i}|^{4+\delta}<\infty$, we know $\operatorname{GMC}(4)$ holds. Then we have the following results.

Lemma 8.4**.**

Under the assumptions of Theorem 4.1, $\operatorname{GMC}(4)$ holds with $0<\rho<1$, then

[TABLE]

[TABLE]

[TABLE]

Proof.

See Section A.4. ∎

Next, we apply the block method to $\{\tilde{Y}_{u,i}(\theta)\}$. Define

[TABLE]

where $k_{n}:=\lfloor n/(p_{n}+q_{n})\rfloor$. Let $p_{n}=q_{n}=\lfloor n^{1-4\eta/\delta}(\log n)^{-8/\delta-4}\rfloor$ (i.e. same block length) and $\ell=\ell_{n}=\lfloor-9\log n/\log\rho\rfloor$ (Note $B_{n}=o(p_{n})$ since $\eta<\delta/(4+\delta)$). Then $U_{u,r}(\theta),r=1,\dots,k_{n}$ are independent (not identically distributed) block sums with block length $p_{n}$, and $V_{u,r}(\theta),r=1,\dots,k_{n}-1$ are independent block sums with block length $q_{n}$. Define $U_{u,r}^{\prime}(\theta):=U_{u,r}(\theta)\mathbf{1}(|U_{u,r}(\theta)|\leq d_{n})$ where $d_{n}=\lfloor\sqrt{nB_{n}}(\log n)^{-1/2}\rfloor$. Then we have the following results.

Lemma 8.5**.**

Under the assumptions of Theorem 4.1, we have that

[TABLE]

[TABLE]

[TABLE]

Furthermore, we have that

[TABLE]

where the $o(1)$ term holds uniformly over $\theta$, $r$ and $u$.

Proof.

See Section A.5. ∎

Lemma 8.6**.**

Let $U_{u,i}(\theta)$ be one of the block sums with block length $p_{n}$. Then we have that

[TABLE]

Proof.

See Section A.6. ∎

Using the previous results Eqs. 68, 70 and 73, we have that

[TABLE]

First, we can show that the third term of the right hand side of Eq. 77 is $o_{\mathbb{P}}(\sqrt{\log n})$. This is because by Eq. 72, it suffices to show $\frac{\sqrt{p_{n}\ell_{n}}B_{n}}{\sqrt{nB_{n}}}=o(\sqrt{\log n})$ and this can be easily verified using $p_{n}=n^{1-4\eta/\delta}(\log n)^{-8/\delta-4}$, $B_{n}=\mathcal{O}(n^{\eta})$ and $\delta\leq 4$.

Next, we show that the right hand side of the first two terms of Eq. 77 have a order of $\mathcal{O}_{\mathbb{P}}(\sqrt{\log n})$. Let $H_{u,n}(\theta)=\sum_{r=1}^{k_{n}}[U_{u,r}(\theta)-\mathbb{E}(U_{u,r}(\theta))]$ and $H_{u,n}^{\prime}(\theta)=\sum_{r=1}^{k_{n}}[U_{u,r}^{\prime}(\theta)-\mathbb{E}(U_{u,r}^{\prime}(\theta))]$. Let $\theta_{j}=\pi j/t_{n},j=0,\dots,t_{n}$ where $t_{n}=\lfloor B_{n}\log(B_{n})\rfloor$. Then, since both $H_{u,n}$ and $H_{u,n}^{\prime}$ have trigonometric polynomial forms, we can apply the following result from [WN67, Corollary 2.1].

Lemma 8.7**.**

Let $p(\lambda)=\sum_{v=-k}^{k}\alpha_{v}\exp(iv\lambda)$ be a trigonometric polynomial. Let $\lambda_{i}=\pi(i/rk),|i|\leq rk$. Then $\operatorname*{\mathrm{max}\vphantom{\mathrm{infsup}}}_{|\lambda|\leq\pi}|p(\lambda)|\leq\operatorname*{\mathrm{max}\vphantom{\mathrm{infsup}}}_{|i|\leq rk}|p(\lambda_{i})/(1-3\pi r^{-1})|$.

Proof.

See [WN67, Corollary 2.1]. ∎

By setting $k=B_{n}$ and $r=\log(B_{n})$ in Lemma 8.7, we get

[TABLE]

By Eqs. 74 and 75, there exists a constant $C_{1}$ such that

[TABLE]

Let $\alpha_{n}:=(C_{1}nB_{n}\log n)^{1/2}$, by the union upper bound,

[TABLE]

Then we apply Bernstein’s inequality (see Lemma A.3) to $\mathbb{P}(|H_{u,n}^{\prime}(\theta_{j})|\geq 4\alpha_{n})$. This leads to, uniformly over $u$ and $\theta_{j}$,

[TABLE]

Therefore, uniformly over $u$, we have that $\mathbb{P}(\operatorname*{\mathrm{max}\vphantom{\mathrm{infsup}}}_{0\leq j\leq t_{n}}|H_{u,n}^{\prime}(\theta_{j})|\geq 4\alpha_{n})=\mathcal{O}(t_{n})\mathcal{O}(1/n)=o(1)$. Let $U_{u,n}^{*}(\theta)=U_{u,n}(\theta)-U_{u,n}^{\prime}(\theta)$ and $H_{u,n}^{*}(\theta)=H_{u,n}(\theta)-H_{u,n}^{\prime}(\theta)$. By the union upper bound and Chebyshev’s inequality

[TABLE]

Using Lemma 8.6, $\operatorname*{\mathrm{sup}\vphantom{\mathrm{infsup}}}_{u}\operatorname*{\mathrm{max}\vphantom{\mathrm{infsup}}}_{i}\operatorname*{\mathrm{sup}\vphantom{\mathrm{infsup}}}_{\theta}\|U_{u,i}(\theta)\|_{2+\delta/2}=\mathcal{O}(\ell_{n}\sqrt{p_{n}B_{n}})$, and

[TABLE]

we have that

[TABLE]

Using $p_{n}=n^{1-4\eta/\delta}(\log n)^{-8/\delta-4}$ we have $p_{n}^{\delta/4}=(n^{\delta/4-\eta})(\log n)^{-2-\delta}$. Therefore,

[TABLE]

Finally, $B_{n}=\mathcal{O}(n^{\eta}),\delta\leq 4$ implies $B_{n}^{\delta/4}=\mathcal{O}(n^{\eta})$, so we have that $\sum_{j=0}^{t_{n}}\frac{\sum_{i=1}^{k_{n}}\textrm{var}(U_{u,i}^{*}(\theta_{j}))}{16\alpha_{n}^{2}}=o(1)$. Therefore, uniformly over $u$, we have $\operatorname*{\mathrm{max}\vphantom{\mathrm{infsup}}}_{\theta}|H_{u,n}^{\prime}(\theta)|=\mathcal{O}_{\mathbb{P}}(\alpha_{n})$ and $\operatorname*{\mathrm{max}\vphantom{\mathrm{infsup}}}_{\theta}|H_{u,n}^{*}(\theta)|=\mathcal{O}_{\mathbb{P}}(\alpha_{n})$. Then $\operatorname*{\mathrm{max}\vphantom{\mathrm{infsup}}}_{\theta}|H_{u,n}(\theta)|=\operatorname*{\mathrm{max}\vphantom{\mathrm{infsup}}}_{\theta}|H_{u,n}^{\prime}(\theta)+H_{u,n}^{*}(\theta)|=\mathcal{O}_{\mathbb{P}}(\alpha_{n})=\mathcal{O}_{\mathbb{P}}(\sqrt{nB_{n}\log n})$. So Eq. 77 has the order of $\mathcal{O}_{\mathbb{P}}(\sqrt{\log n})$.

8.3. Proof of Theorem 4.3

Throughout the proof, we use $\|\cdot\|$ to denote $\|\cdot\|_{2}$ for simplicity. We define $Y_{u,i}$, $g_{n}$, $h_{n}$, $\tilde{X}_{k}^{[\ell]}$, $\tilde{Y}_{u,i}$, $\tilde{g}_{n}$ the same as in Section 8.2. Therefore, Lemma 8.4 holds. Next, we apply the block method to $\{\tilde{Y}_{u,i}(\theta)\}$. Define

[TABLE]

where $k_{n}:=\lfloor n/(p_{n}+q_{n})\rfloor$. Let $\psi_{n}=n/(\log n)^{2+8/\delta}$, $p_{n}=\lfloor\psi_{n}^{2/3}B_{n}^{1/3}\rfloor$, and $q_{n}=\lfloor\psi_{n}^{1/3}B_{n}^{2/3}\rfloor$. Then we have $p_{n},q_{n}\to\infty$ and $q_{n}=o(p_{n})$. Since $\ell_{n}=\mathcal{O}(\log n)$, we have $2B_{n}+\ell_{n}=o(q_{n})$ and $k_{n}=\lfloor n/(p_{n}+q_{n})\rfloor\to\infty$. Note that $U_{u,r}(\theta),r=1,\dots,k_{n}$ are independent (not identically distributed) block sums with block length $p_{n}$, and $V_{u,r}(\theta),r=1,\dots,k_{n}$ are independent block sums with block length $q_{n}$. Now the proof of Lemma 8.5 still follows.

Defining $a_{n}/b_{n}\to 1$ by $a_{n}\sim b_{n}$, we have the following result.

Lemma 8.8**.**

Let the sequence $s_{n}\in\mathbb{N}$ satisfy $s_{n}\leq n$, $s_{n}=o(n)$ and $B_{n}=o(s_{n})$. Under $\textrm{GMC}(4)$ we have that

[TABLE]

where $\sigma^{2}_{u}(\theta)=[1+\eta(2\theta)]f^{2}(u,\theta)\int_{-1}^{1}a^{2}(t)\mathrm{d}t$ and $\eta(\theta)=1$ if $\theta=2k\pi$ for some integer $k$ and $\eta(\theta)=0$ otherwise.

Proof.

See Section A.7. ∎

According to Lemmas 8.8 and 8.4, for each block $U_{u,r},r=1,\dots,k_{n}$, we have that

[TABLE]

where $\mathcal{L}_{r}=\{j\in\mathbb{N}:(r-1)(p_{n}+q_{n})+1\leq j\leq r(p_{n}+q_{n})-q_{n}\}$. Similarly, we can also show that $\|V_{u,r}-\mathbb{E}(V_{u,r})\|\sim(q_{n}B_{n}\sigma^{2}_{u})^{1/2}+\mathcal{O}(q_{n}B_{n}\rho^{\ell_{n}/4})$. Then, since $q_{n}=o(p_{n})$, we have that

[TABLE]

which implies that $\frac{\sum_{r}(V_{u,r}-\mathbb{E}(V_{u,r}))}{\sqrt{nB_{n}}}\Rightarrow 0$. Also, by Eq. 68, we have that $\textrm{var}(h_{n}(u,\theta))=\mathcal{O}(B_{n}/n)=\mathcal{O}((\log n)^{-2-8/\delta})$, which implies that $h_{n}(u,\theta)-\mathbb{E}(h_{n}(u,\theta))\Rightarrow 0$. Therefore, by

[TABLE]

we only need to show that $\frac{\sum_{r}(U_{u,r}-\mathbb{E}(U_{u,r}))}{\sqrt{nB_{n}}}\Rightarrow\mathcal{N}(0,\sigma^{2}_{u})$. We can check the conditions of Lemma A.2 (the Berry–Esseen lemma) as follows.

[TABLE]

By Lemma 8.6, we know $\|U_{u,r}\|_{2+\delta/2}=\mathcal{O}(\ell_{n}\sqrt{p_{n}B_{n}})$, which implies

[TABLE]

Note that $k_{n}=\lfloor n/(p_{n}+q_{n})\rfloor\sim n\psi^{-2/3}B_{n}^{-1/3}\sim n^{1/3}(\log n)^{(4/3+16/3\delta)}B_{n}^{-1/3}$, which implies $k_{n}^{-1}=\mathcal{O}((\log n)^{-4/3-16/3\delta})$. Then $\ell_{n}k_{n}^{-\delta/4}=\mathcal{O}((\log n)(\log n)^{(-\delta/3-4/3)})=\mathcal{O}((\log n)^{(-\delta/3-1/3)})\to 0$. Therefore, the result holds by Lemma A.2.

8.4. Proof of Theorem 5.3

Define $D_{n}=C_{n}B_{n}$, $\theta_{i}=\frac{i\pi}{B_{n}},i=0,\dots,B_{n}$, and $\alpha_{n,k}=a(k/B_{n})\cos(k\theta)$. We use the previous definitions of $X_{u,k}$ and the $\ell$-dependent $\tilde{X}_{u,k}^{[\ell]}$ as in Section 8.2. Let $g_{n}(u,\theta):=[2\pi n\hat{f}_{n}(u,\theta)-\sum_{k=1}^{n}X_{u,k}^{2}]-\mathbb{E}[2\pi n\hat{f}_{n}(u,\theta)-\sum_{k=1}^{n}X_{u,k}^{2}]$, where $\ell=\lfloor n^{\gamma}\rfloor$ for fixed $\gamma>0$ which is close to zero. Note that

[TABLE]

Therefore, we have $g_{n}(u,\theta)=\sum_{1\leq k,k^{\prime}\leq n,k\neq k^{\prime}}\alpha_{n,k-k^{\prime}}[X_{u,k}X_{u,k^{\prime}}-\mathbb{E}(X_{u,k}X_{u,k^{\prime}})]$. Then let $\tilde{g}_{n}(u,\theta)$ be the corresponding version of $g_{n}(u,\theta)$ using $\ell$-dependent $\{\tilde{X}_{u,k}^{[\ell]}\}$ instead of $\{X_{u,k}\}$. Define $X_{u,k}^{\prime}=\tilde{X}_{u,k}^{[\ell]}\mathbf{1}_{\left|\tilde{X}_{u,k}^{[\ell]}\right|\leq(nB_{n})^{\alpha}}$ where $\alpha<\frac{1}{4}$. Next, let $\bar{X}_{u,k}:=X_{u,k}^{\prime}-\mathbb{E}X_{u,k}^{\prime}$ and define

[TABLE]

In the following, we show that $g_{n}(u,\theta)$ can be approximated by $\tilde{g}_{n}(u,\theta)$.

Lemma 8.9**.**

Under the assumptions of Theorem 5.3, we have $\operatorname*{\mathrm{max}\vphantom{\mathrm{infsup}}}_{u\in\mathcal{U}}\operatorname*{\mathrm{max}\vphantom{\mathrm{infsup}}}_{0\leq i\leq B_{n}}\mathbb{E}|g_{n}(u,\theta_{i})-\tilde{g}_{n}(u,\theta_{i})|=o(n^{1+\gamma}\rho^{\lfloor n^{\gamma}\rfloor})$ and $\operatorname*{\mathrm{max}\vphantom{\mathrm{infsup}}}_{u\in\mathcal{U}}\operatorname*{\mathrm{max}\vphantom{\mathrm{infsup}}}_{0\leq i\leq B_{n}}\frac{|g_{n}(u,\theta_{i})-\tilde{g}_{n}(u,\theta_{i})|}{\sqrt{nB_{n}}}=o_{\mathbb{P}}(1)$.

Proof.

See Section A.8. ∎

Next, we show that $\tilde{g}_{n}(u,\theta)$ can be approximated by $\bar{g}_{n}(u,\theta)$.

Lemma 8.10**.**

Under the assumptions of Theorem 5.3, we have that

[TABLE]

Proof.

See Section A.9. ∎

According to Lemma 8.9 and Lemma 8.10, together with $\operatorname*{\mathrm{max}\vphantom{\mathrm{infsup}}}_{i}|\tilde{g}_{n}(u,\theta_{i})-\bar{g}_{n}(u,\theta_{i})|\leq\operatorname*{\mathrm{max}\vphantom{\mathrm{infsup}}}_{\theta}|\tilde{g}_{n}(u,\theta)-\bar{g}_{n}(u,\theta)|$, we have that $\operatorname*{\mathrm{max}\vphantom{\mathrm{infsup}}}_{u\in\mathcal{U}}\operatorname*{\mathrm{max}\vphantom{\mathrm{infsup}}}_{0\leq i\leq B_{n}}\frac{|g_{n}(u,\theta)-\tilde{g}_{n}(u,\theta)|^{2}}{nB_{n}}=o_{\mathbb{P}}(1)$ and

[TABLE]

Since $\operatorname*{\mathrm{max}\vphantom{\mathrm{infsup}}}_{u}\operatorname*{\mathrm{max}\vphantom{\mathrm{infsup}}}_{i}|\mathbb{E}\tilde{g}_{n}(u,\theta_{i})-\mathbb{E}\bar{g}_{n}(u,\theta_{i})|\leq\mathbb{E}(\operatorname*{\mathrm{max}\vphantom{\mathrm{infsup}}}_{u}\operatorname*{\mathrm{max}\vphantom{\mathrm{infsup}}}_{i}|\tilde{g}_{n}(u,\theta_{i})-\bar{g}_{n}(u,\theta_{i})|)$, it suffices to show that, for $D_{n}=B_{n}C_{n}$, we have that

[TABLE]

Let $p_{n}=\lfloor B_{n}^{1+\beta}\rfloor$, $q_{n}=B_{n}+\ell$, $\ell=\lfloor n^{\gamma}\rfloor$ and $k_{n}=\lfloor n/(p_{n}+q_{n})\rfloor$, where $\gamma$ is small enough and $\beta>0$ is sufficiently close to zero. Split the interval $[1,n]$ into alternating big and small blocks $H_{j}$ and $I_{j}$ by

[TABLE]

Define $\bar{Y}_{u,k}:=\bar{X}_{u,k}\sum_{s=1}^{k-1}\alpha_{n,k-s}\bar{X}_{u,s}$. Then $\bar{g}_{n}=\sum_{k=1}^{n}(\bar{Y}_{u,k}-\mathbb{E}\bar{Y}_{u,k})$. For $1\leq j\leq k_{n}+1$, let

[TABLE]

Then $\bar{g}_{n}=\sum_{j=1}^{k_{n}}U_{j}+\sum_{j=1}^{k_{n}+1}V_{j}$. Next, define a truncated and normalized version of $U_{j}$ as

[TABLE]

In the following, we show that $\bar{g}_{n}(u,\theta_{i})-\mathbb{E}(\bar{g}_{n}(u,\theta_{i}))$ can be approximated by $\sum_{j=1}^{k_{n}}\bar{U}_{j}(u,\theta_{i})$.

Lemma 8.11**.**

Under the assumptions of Theorem 5.3, we have that

[TABLE]

Proof.

See Section A.10. ∎

Furthermore, we show in the following that $\sum_{j=1}^{k_{n}}\bar{U}_{j}(u,\theta_{i})$ can be ignored if $i\notin[(\log B_{n})^{2},B_{n}-(\log B_{n})^{2}]$.

Lemma 8.12**.**

Under the assumptions of Theorem 5.3, we have that

[TABLE]

Proof.

See Section A.11. ∎

Finally, we complete the proof of Eq. 24 by the following result.

Lemma 8.13**.**

Under the assumptions of Theorem 5.3, we have that

[TABLE]

Proof.

See Section A.12. ∎

Acknowledgement

The authors are grateful to the anonymous referees for their many helpful comments and suggestions which significantly improved the quality of the paper.

\AtNextBibliography

A Supplemental Material

Remark A.1*.*

Denote $X_{u,i}:=G(i/N,\mathcal{F}_{u,i})$ where $\mathcal{F}_{u,i}=(\dots,\epsilon_{\lfloor uN\rfloor},\epsilon_{\lfloor uN\rfloor+1},\dots,\epsilon_{\lfloor uN\rfloor+i})$. Let $\epsilon_{k}^{\prime}$ be an i.i.d. copy of $\epsilon_{k}$ and $X_{u,i}^{\prime}:=G(i/N,\mathcal{F}_{u,i}^{\prime})$ where $\mathcal{F}_{u,i}^{\prime}=(\dots,\epsilon_{0}^{\prime},\dots,\epsilon_{\lfloor uN\rfloor}^{\prime},\epsilon_{\lfloor uN\rfloor+1},\dots,\epsilon_{\lfloor uN\rfloor+i})$ is a coupled version of $\mathcal{F}_{u,i}$. Then under $\operatorname{GMC}(p)$, $p>0$, there exist $C>0$ and $0<\rho=\rho(p)<1$ that do not depend on $u$, such that for any $u$ and $i$, we have

[TABLE]

This is because, when $\operatorname{GMC}(p)$ holds, we have $\operatorname*{\mathrm{sup}\vphantom{\mathrm{infsup}}}_{u}\mathbb{E}(|X_{u,i}^{\prime}-X_{u,i}|^{p})\leq\sum_{k=i}^{\infty}\delta_{p}(k)\leq\mathcal{O}(\sum_{k=i}^{\infty}\rho^{k})=\mathcal{O}(\rho^{i})$.

Furthermore, it can be easily shown that if $\textrm{GMC}(2)$ holds, then $\operatorname*{\mathrm{sup}\vphantom{\mathrm{infsup}}}_{u}|r(u,k)|=\mathcal{O}(\rho^{k})$ for some $\rho\in(0,1)$. Also, if $\operatorname*{\mathrm{sup}\vphantom{\mathrm{infsup}}}_{i}\|X_{i}\|_{p}<\infty$ and $\operatorname{GMC}(\alpha)$ holds with any given $\alpha>0$, then $X_{i}$ is $\operatorname{GMC}(\alpha)$ with any $\alpha\in(0,p)$. In particular, if $\operatorname{GMC}(\alpha)$ holds with some $\alpha\geq 2$, then we must have $\operatorname*{\mathrm{sup}\vphantom{\mathrm{infsup}}}_{u}\sum_{k=-\infty}^{\infty}|r(u,k)|<\infty$ since $\operatorname*{\mathrm{sup}\vphantom{\mathrm{infsup}}}_{u}|r(u,k)|=\mathcal{O}(\rho^{k})=o(k^{-2})$. Also, if $\operatorname{GMC}(2)$ holds as well as $\operatorname*{\mathrm{sup}\vphantom{\mathrm{infsup}}}_{i}\mathbb{E}(|X_{i}|^{4+\delta})<\infty$ for some $\delta>0$, then $\operatorname{GMC}(4)$ holds.

Lemma A.2**.**

(Berry-Esseen) If $\{X_{i},i\geq 1\}$ are independent random variables with $\mathbb{E}(X_{i})=0$, $s_{n}^{2}=\sum_{i=1}^{n}\mathbb{E}(X_{i}^{2})>0$, $\sum_{i=1}^{n}\mathbb{E}|X_{i}|^{2+\delta}<\infty$, for some $\delta\in(0,1]$ and $S_{n}=\sum_{i=1}^{n}X_{i}$, there exists a universal constant $C_{\delta}$ such that

[TABLE]

Proof.

See [CT88, pp. 304]. ∎

A.1. Proof of Lemma 8.1

Define $d_{u,n}(h)=\frac{1}{n}\sum_{k=1+h}^{n}\mu_{u,k}\mu_{u,k-h}$ for $0\leq h\leq n-1$ and $d_{u,n}(h)=0$ if $h\geq n$. Since

[TABLE]

using

[TABLE]

we get that uniformly over $J$, $c$ and $u$, there exists $K_{0}$ such that

[TABLE]

Next, we can write $\|T_{u,n}\|^{2}/n$ as

[TABLE]

Furthermore, defining

[TABLE]

we have that

[TABLE]

By the assumptions that $\tau\in\mathcal{C}^{1}([-1/2,1/2])$, $\int\tau^{2}(x)\mathrm{d}x=1$, together with $\operatorname*{\mathrm{sup}\vphantom{\mathrm{infsup}}}_{u}|r(u,h)|=o(h^{-2})$, and $\sum_{h=1}^{\infty}|r(u,h)|<\infty$, we have $f_{n}(u,\theta)=f(u,\theta)+o(1)$, uniformly over $u$ and $\theta$. This implies that

[TABLE]

Therefore, uniformly over $J$ and $c$, we have that

[TABLE]

Finally, since $\operatorname*{\mathrm{sup}\vphantom{\mathrm{infsup}}}_{u}\sum_{h}|r(u,h)|<\infty$, we have $\operatorname*{\mathrm{sup}\vphantom{\mathrm{infsup}}}_{u}\sum_{h>n}|r(u,h)|\to 0$. Also, as $n\to\infty$,

[TABLE]

Therefore, $\left|\frac{\|T_{u,n}\|^{2}}{n}-1\right|\to 0$.

A.2. Proof of Lemma 8.2

Throughout the proof, we write $\tilde{X}_{k}^{[\ell]}$ as $\tilde{X}_{k}$ for short. Note that $\int\tau(x)\tau(x+h)\mathrm{d}x\leq\frac{1}{2}\int\left[\tau(x)^{2}+\tau(x+h)^{2}\right]\mathrm{d}x=1$. For simplicity of the proof, we can assume that there exists some finite $\tau_{*}$ such that

[TABLE]

Then we have that

[TABLE]

Since the upper bound does not depend on $u$, the convergence holds uniformly over $u$.

A.3. Proof of Lemma 8.3

In this proof, we omit subscript $u$ for simplicity and write $\tilde{X}_{k}^{[\ell]}$ as $\tilde{X}_{k}$ for short. Since $\operatorname*{\mathrm{sup}\vphantom{\mathrm{infsup}}}_{k}\mathbb{E}(X_{k}^{2})<\infty$, we have that

[TABLE]

By the property of conditional expectation, we have $\mathbb{E}(\tilde{X}_{k}^{2})<\mathbb{E}(X_{k}^{2})$. Therefore, defining

[TABLE]

we can get $\operatorname*{\mathrm{lim}\vphantom{\mathrm{infsup}}}_{n\to\infty}g_{n}(r)=0$ for all given $r>0$. Also $g_{n}$ is non-decreasing with $r$. Then there exists a sequence $\{r_{n}\}$ such that $r_{n}\!\uparrow\!\infty$ and $g_{n}(r_{n})\to 0$. Note that $r_{n}$ does not depend on $u$.

For simplicity, we use $\tilde{X}_{u,k}$ to denote $\tilde{X}_{\lfloor uN\rfloor+k-\lfloor n/2\rfloor}$. Let $Y_{u,k}=\tilde{X}_{u,k}\mathbf{1}(|\tilde{X}_{u,k}|\leq\sqrt{n}/r_{n})$ and $T_{u,n,Y}=\sum_{k=1}^{n}\mu_{u,k}Y_{u,k}$. Since $\mathbb{E}[\tilde{X}_{k}^{2}\mathbf{1}(|\tilde{X}_{k}|\geq\sqrt{n}/r)]=o(1/r_{n}^{2})$ by the definition of $r_{n}$, we have $\|Y_{u,k}-\tilde{X}_{u,k}\|=o(1/r_{n})$. Now since $Y_{u,k}-\tilde{X}_{u,k}$ is $\ell$-dependent, we divide each of $\{Y_{u,k}\}$ and $\{\tilde{X}_{u,k}\}$ into $\ell$ sub-sequences that each sub-sequences has $\lfloor n/\ell\rfloor$ independent elements. Then by the triangle inequality we can get

[TABLE]

Next, divide the sequence of $\{Y_{u,k}\}$ into pieces of length $p_{n}+\ell$ where $p_{n}=\lfloor r_{n}^{1/4}\rfloor$.

[TABLE]

where $B_{t}=\{a\in\mathbb{N}:1+(t-1)(p_{n}+\ell)\leq a\leq p_{n}+(t-1)(p_{n}+\ell)\}$. Note that for given $u$, $\{U_{u,t}\}$ are independent (but not identically distributed) for different $t$.

Define $V_{u,t}=\sum_{t=1}^{t_{n}}U_{u,t}$, then the difference between $V_{u,t}$ and $T_{u,n,Y}$ is the sum of those dropped $\ell$ terms in each piece. Since $\ell$ is fixed and there are $t_{n}$ blocks, we have $\|T_{u,n,Y}-V_{u,t}\|=\mathcal{O}(\sqrt{t_{n}})$.

Furthermore, since

[TABLE]

we have $P(|\tilde{X}_{k}|\geq\sqrt{n}/r))=o(1/n)$. Then, using

[TABLE]

we have $\mathbb{E}(Y_{k})=o(\frac{1}{\sqrt{n}r_{n}})$, which implies $|\mathbb{E}(V_{n})|=O(n)|\mathbb{E}(Y_{k})|=o(\sqrt{n}/r_{n})$.

Next, defining $W=(V_{n}-\mathbb{E}(V_{n}))/\sqrt{n}$ and $\Delta=\tilde{T}_{n}/\sqrt{n}-W$, we get

[TABLE]

Next, we apply Lemma A.2 to $\{U_{t}-\mathbb{E}(U_{t}),t=1,\dots,t_{n}\}$. Recall that $V_{n}=\sum_{t=1}^{t_{n}}U_{t}$ and $W=(V_{n}-\mathbb{E}(V_{n}))/\sqrt{n}$, then

[TABLE]

Next, we get upper bounds of $\mathbb{E}|U_{t}|^{3}$ and $\|V_{n}-\mathbb{E}(V_{n})\|^{-3}$. First, by Hölder’s inequality $\sum_{a\in B_{t}}|Y_{a}|\leq(\sum_{a\in B_{t}}|Y_{a}|^{3})^{1/3}(\sum_{a\in B_{t}}1)^{2/3}$, we have that

[TABLE]

For sequences $a_{n}$ and $b_{n}$, we define $a_{n}=\Theta(b_{n})$ if both $a_{n}=\mathcal{O}(b_{n})$ and $b_{n}=\mathcal{O}(a_{n})$. Then, using the definition of $\Theta(\cdot)$, the variance of $\sum_{a\in B_{t}}\mu_{a}Y_{a}$ has the order of $\Theta(p_{n})$ because $Y_{a}$ is $\ell$-dependent. Then the variance of $V_{n}$ has the order of $\Theta(t_{n}p_{n})=\Theta(n)$. Thus, $\|V_{n}-\mathbb{E}(V_{n})\|^{-3}$ has an order of $\Theta(n^{-3/2})$. Overall, we have that

[TABLE]

To complete the proof, we first replace $V_{n}=\sum_{t}\sum_{a\in B_{t}}\mu_{a}Y_{a}$ by $\tilde{T}_{n}=\sum_{k}\mu_{k}\tilde{X}_{k}$ then by $T_{n}=\sum_{k}X_{k}$. Since

[TABLE]

we have that

[TABLE]

Furthermore, one can get

[TABLE]

Using

[TABLE]

we get

[TABLE]

Also

[TABLE]

Letting $\delta=p_{n}^{-1/4}$ we have that

[TABLE]

Finally, use the above technique again with $\Delta_{1}=(T_{n}-\tilde{T}_{n})/\sqrt{n}$ and $\delta=\|\Delta_{1}\|^{1/2}$, we get

[TABLE]

Lemma A.3**.**

(Bernstein’s inequality) Let $X_{1},\dots,X_{n}$ be independent zero-mean random variables. Suppose $|X_{i}|\leq M$ a.s., for all $i$. Then for all positive $t$,

[TABLE]

Definition A.4**.**

Let $(U_{1},\dots,U_{k})$ be a random vector. Then the joint cumulant is defined as

[TABLE]

where $V_{1},\dots,V_{p}$ is a partition of the set $\{1,2,\dots,k\}$ and the sum is taken over all such partitions.

Lemma A.5**.**

Assume $\textrm{GMC}(\alpha)$ with $\alpha=k$ for some $k\in\mathbb{N}$, and $\operatorname*{\mathrm{sup}\vphantom{\mathrm{infsup}}}_{t}\mathbb{E}(|X_{t}|^{k})<\infty$ Then there exists a constant $C>0$ such that for all $u$ and $0\leq m_{1}\leq\dots\leq m_{k-1}$,

[TABLE]

where $X_{u,i}:=\tau\left(\frac{i-\lfloor n\rfloor/2\rfloor}{n}\right)X_{\lfloor uN\rfloor+i-\lfloor n/2\rfloor}$.

Proof.

Since $\tau(\cdot)$ is bounded, we have $\operatorname*{\mathrm{sup}\vphantom{\mathrm{infsup}}}_{u}\operatorname*{\mathrm{sup}\vphantom{\mathrm{infsup}}}_{i}\mathbb{E}(|X_{u,i}|^{k})<\infty$. We extend [WS04, Proposition 2] to the cases of locally stationary time series.

Given $1\leq l\leq k-1$, by multi-linearity of joint cumulants, we replace $X_{u,m_{i}}$ by independent $X_{u,m_{i}}^{\prime}$ for all $i\geq l$ as follows

[TABLE]

Note that $(X_{u,0},X_{u,m_{1}},\dots,X_{u,m_{l-1}})$ is independent with $(X_{u,m_{l}}^{\prime},\dots,X_{u,m_{k-1}}^{\prime})$. By [Ros85, pp.35], we have $B=0$. Suppose we have that

[TABLE]

for $l\leq i\leq k-1$ and some constant $C$ that does not depend on $l$. Then $|J|\leq C\rho^{(m_{l}-m_{l-1})/k}$ for any $1\leq l\leq k-1$. Then we get

[TABLE]

Next, we show Eq. 137. In particular, we show the case $i=l$ and the other cases can be proven similarly. Note that $\mathbb{E}(|X_{u,i}|^{k})$ is uniformly bounded, by the definition of joint cumulants in Definition A.4, we only need to show that for $V\subset\{0,\dots,k-1\}$ such that $l\notin V$, we have that

[TABLE]

Letting $|V|$ be the cardinality of the set $V$, then $|V|\leq k-1$, and we have

[TABLE]

By Hölder’s inequality and Jensen’s inequality

[TABLE]

∎

A.4. Proof of Lemma 8.4

Throughout this proof, we write $\tilde{X}_{k}^{[\ell]}$ as $\tilde{X}_{k}$ and $\tilde{X}_{u,i}^{[\ell]}$ as $\tilde{X}_{u,i}$ for short. First, letting $\alpha_{k}=a(k/B_{n})\cos(k\theta)$, we have that

[TABLE]

By the summability of cumulants of orders $2$ and $4$ [Ros85, page 185], one can get

[TABLE]

Therefore, we have $\operatorname*{\mathrm{sup}\vphantom{\mathrm{infsup}}}_{\theta}\operatorname*{\mathrm{sup}\vphantom{\mathrm{infsup}}}_{u}\|h_{n}(u,\theta)\|=(nB_{n})^{-1/2}\mathcal{O}(B_{n})$.

Next, note that by the assumption of GMC($4$) defined in Eq. 102, we have that

[TABLE]

Then we have

[TABLE]

Finally

[TABLE]

A.5. Proof of Lemma 8.5

We write $\tilde{X}_{k}^{[\ell]}$ as $\tilde{X}_{k}$ and $\tilde{X}_{u,i}^{[\ell]}$ as $\tilde{X}_{u,i}$ for short. To show Eq. 72, since $\alpha_{k}$ is bounded, letting $z_{n}=k_{n}(p_{n}+q_{n})+1-q_{n}$, we have that

[TABLE]

Since $\tilde{X}_{u,i}\tilde{X}_{u,i+j}$ is $2\ell_{n}$-dependent, if $|j|<\ell_{n}$, we have

[TABLE]

If $|j|\leq\ell_{n}$, since $\mathbb{E}(\tilde{X}_{u,i}\tilde{X}_{u,i+j}\tilde{X}_{u,i^{\prime}}\tilde{X}_{u,i^{\prime}+j})=0$ if $|i-i^{\prime}|>\ell_{n}$, we have that

[TABLE]

where we have used the assumption $\operatorname*{\mathrm{sup}\vphantom{\mathrm{infsup}}}_{i}\mathbb{E}(|X_{i}|^{4+\delta})<M$. Therefore, we get Eq. 72.

To show Eq. 73, we first define $\tilde{h}_{n}(u,\theta)$ by replacing $X_{i}$ by $\tilde{X}_{i}$. Then we can prove similarly to Eq. 145 that

[TABLE]

Therefore, it suffices to show $\operatorname*{\mathrm{sup}\vphantom{\mathrm{infsup}}}_{u}\mathbb{E}(\operatorname*{\mathrm{max}\vphantom{\mathrm{infsup}}}_{\theta}|\tilde{h}_{n}(u,\theta)|)=o(1)$. Using similar technique to Eq. 147 we can show that

[TABLE]

where we have used $\eta<\frac{1}{2}$ and $\sqrt{\ell_{n}}B_{n}/\sqrt{n}=\mathcal{O}((\log n)^{1/2}n^{\eta-1/2})=o(1)$.

To show Eq. 74, we note that GMC($4$) implies the absolute summability of cumulants up to the fourth order. Also, for zero-mean random variables $X,Y,Z,W$, the joint cumulants

[TABLE]

Therefore, letting $\mathcal{L}_{r}$ be the set of the indices $i$’s such that $Y_{u,i}$ belongs to the block corresponding to $U_{u,r}$, we have that

[TABLE]

where the first term is finite since the fourth cumulants are summable. For the second term (the last term can also be shown similarly), we use the condition Eq. 11, so that

[TABLE]

Then using $p_{n}=o(n)$, $B_{n}=o(p_{n})$ and $\operatorname*{\mathrm{sup}\vphantom{\mathrm{infsup}}}_{u}\sum_{k=-\infty}^{\infty}|r(u,k)|<\infty$, one can get

[TABLE]

To show Eq. 75, we note that

[TABLE]

From Lemma 8.8, we know that $\textrm{var}(U_{u,r}(\theta))\sim p_{n}B_{n}\sigma^{2}_{u}(\theta)$ and $\sigma^{2}_{u}(\theta)=[1+\eta(2\theta)]f^{2}(u,\theta)\int_{-1}^{1}a^{2}(t)\mathrm{d}t\geq f_{*}^{2}\int_{-1}^{1}a^{2}(t)\mathrm{d}t>0$. Thus, it suffices to show that

[TABLE]

By Lemma 8.6, applying similar inequalities as Eq. 82, we have that

[TABLE]

Finally, since $\mathbb{E}(U_{u,r}^{\prime})\leq\mathbb{E}(|U_{u,r}|)\leq[\mathbb{E}(|U_{u,r}|^{2+\delta/2})]^{\frac{1}{2+\delta/2}}$, using again similar inequalities as Eq. 82, we have that

[TABLE]

A.6. Proof of Lemma 8.6

In this proof, we write $\tilde{X}_{k}^{[\ell]}$ as $\tilde{X}_{k}$ and $\tilde{X}_{u,i}^{[\ell]}$ as $\tilde{X}_{u,i}$ for short. For simplicity, we first consider that $u$ and $i$ are fixed. Without loss of generality, we consider the first block sum ($i=1$) so

[TABLE]

We will first show that

[TABLE]

where $\alpha_{k}=a(k/B_{n})\cos(k\theta)$. Then we conclude that $\mathcal{O}(\ell_{n}\sqrt{p_{n}B_{n}})$ is also uniformly over $u$ and $i$ since the assumption $\operatorname*{\mathrm{sup}\vphantom{\mathrm{infsup}}}_{u}\operatorname*{\mathrm{sup}\vphantom{\mathrm{infsup}}}_{i}\mathbb{E}(|X_{u,i}|^{4+\delta})<M$. We first write by the triangle inequality

[TABLE]

Now consider two cases (i) $\ell_{n}=o(B_{n})$, then

[TABLE]

where the first term of the right hand side of Eq. 163 satisfies

[TABLE]

Continuing to divide the sum of $\sum_{k=-B_{n}}^{-\ell_{n}}\tilde{X}_{u,h+(j-1)\ell_{n}+k}\alpha_{k}$ into $\ell_{n}$ parts, then by $\operatorname*{\mathrm{sup}\vphantom{\mathrm{infsup}}}_{u,i}\mathbb{E}(|X_{u,i}|^{4+\delta})<M$, we have that

[TABLE]

which holds uniformly over $u$ and $i$. Similarly, for the second term of the right hand side of Eq. 163

[TABLE]

Note that the order $\mathcal{O}(\ell_{n}^{2}\sqrt{p_{n}/\ell_{n}})$ also holds uniformly over $u$ and $i$. This is because $\|\tilde{X}_{u,h+3j\ell_{n}}\tilde{X}_{u,h+3j\ell_{n}+k}\|_{2+\delta/2}$ is uniformly bounded, which can be shown using Cauchy–Schwarz’s inequality and $\operatorname*{\mathrm{sup}\vphantom{\mathrm{infsup}}}_{u}\operatorname*{\mathrm{sup}\vphantom{\mathrm{infsup}}}_{i}\mathbb{E}(|X_{u,i}|^{4+\delta})<M$. Therefore, we have proven that, for case (i), we have $\operatorname*{\mathrm{sup}\vphantom{\mathrm{infsup}}}_{u}\operatorname*{\mathrm{sup}\vphantom{\mathrm{infsup}}}_{i}\operatorname*{\mathrm{sup}\vphantom{\mathrm{infsup}}}_{\theta}\|U_{u,i}(\theta)\|_{2+\delta/2}=\mathcal{O}(\ell_{n}\sqrt{p_{n}B_{n}})$.

For the second case (ii) $B_{n}=\mathcal{O}(\ell_{n})$, we have that

[TABLE]

which is also uniform over $u$ and $i$.

A.7. Proof of Lemma 8.8

Using the property of cumulants in Eq. 152, similarly to Eqs. 153 and 154, one can get that

[TABLE]

By Lemma A.5, we have that

[TABLE]

which implies that the first term of the right hand side of Eq. 168 is finite.

Finally, according to [Ros84, Theorem 2, Eqs. (3.9)–(3.12)], one can show that

[TABLE]

Lemma A.6**.**

Let $\{X_{k}\}$ be $\ell$-dependent with $\mathbb{E}X_{k}=0$ and $X_{k}\in\mathcal{L}^{p}$ with $p\geq 2$. Let $W_{n}=\sum_{k=1}^{n}X_{k}$. Then for any $Q>0$, there exists $C_{1},C_{2}>0$ only depending on $Q$ such that

[TABLE]

Proof.

See [LW10, Lemma 2]. ∎

Lemma A.7**.**

Let $\{X_{t}\}$ be $\ell$-dependent with $\mathbb{E}X_{t}=0$, $|X_{t}|\leq M$ a.s., $\ell\leq n$, and $M\geq 1$. Let $S_{k,l}=\sum_{t=l+1}^{l+k}X_{t}\sum_{s=1}^{t-1}{\alpha}_{n,t-s}X_{s}$, where $l\geq 0$, $l+k\leq n$ and assume that $\operatorname*{\mathrm{max}\vphantom{\mathrm{infsup}}}_{1\leq t\leq n}|\alpha_{n,t}|\leq K_{0}$, $\operatorname*{\mathrm{max}\vphantom{\mathrm{infsup}}}_{1\leq t\leq n}\mathbb{E}X_{t}^{2}\leq K_{0}$, $\operatorname*{\mathrm{max}\vphantom{\mathrm{infsup}}}_{1\leq t\leq n}\mathbb{E}X_{t}^{4}\leq K_{0}$ for some $K_{0}>0$. Then for any $x\geq 1$, $y\geq 1$, and $Q>0$,

[TABLE]

where $C_{1},C_{2}>0$ are constants depending only on $Q$ and $K_{0}$.

Proof.

See [LW10, Proposition 3]. ∎

Lemma A.8**.**

Assume that $X_{k}\in\mathcal{L}^{p}$, with $p>1$, and $\mathbb{E}X_{k}=0$. Let $C_{p}=18p^{3/2}(p-1)^{-1/2}$ and $p^{\prime}=\operatorname*{\mathrm{min}\vphantom{\mathrm{infsup}}}(2,p)$. Let $\alpha_{1},\dots,\in\mathbb{C}$. Then under GMC, we have that

[TABLE]

and

[TABLE]

for some constant $C$.

Proof.

This lemma follows from [LW10, Lemma 1] with $\Theta_{\ell+1,p}=o(\sum_{j=\ell+1}^{\infty}\rho^{j})=o(\rho^{\ell})$. ∎

Lemma A.9**.**

Assume $\mathbb{E}X_{u,k}=0$, $\operatorname*{\mathrm{sup}\vphantom{\mathrm{infsup}}}_{u}\mathbb{E}|X_{u,k}|^{2p}<\infty$, $p\geq 2$. Let

[TABLE]

where $\alpha_{1},\dots,\in\mathbb{C}$. Then under GMC, we have that

[TABLE]

Proof.

For fixed $u$, if $\mathbb{E}|X_{u,k}|^{2p}<\infty$, the result follow from [LW10, Proposition 1] with $\Theta_{0,2p}=o(1)$ and $d_{\ell,2p}=\sum_{t=0}^{\infty}\operatorname*{\mathrm{min}\vphantom{\mathrm{infsup}}}\{o(\rho^{t}),o(\rho^{\ell})\}=o(\ell\rho^{\ell})$. Since we have $\operatorname*{\mathrm{sup}\vphantom{\mathrm{infsup}}}_{u}\mathbb{E}|X_{u,k}|^{2p}<\infty$ the proof of [LW10, Proposition 1] also holds uniformly over $u$. ∎

Lemma A.10**.**

Assume that $\mathbb{E}X_{u,k}=0$, $\operatorname*{\mathrm{sup}\vphantom{\mathrm{infsup}}}_{u}\mathbb{E}X_{u,k}^{4}<\infty$ and $\operatorname{GMC}(2)$. Let $\alpha_{j}=\beta_{j}\exp(ij\theta)$, where $i=\sqrt{-1}$, $\theta\in\mathbb{R}$, $\beta_{j}\in\mathbb{R}$, $1-n\leq j\leq-1$, $m\in\mathbb{N}$ and $\tilde{L}_{n,u}=\sum_{1\leq j<t\leq n}\alpha_{j-t}{\tilde{X}_{u,j}^{[\ell]}\tilde{X}_{u,t}^{[\ell]}}$. Define

[TABLE]

where $(\cdot)^{*}$ denotes the complex conjugate, $A_{u,k}=\sum_{t=0}^{\infty}\mathbb{E}({\tilde{X}_{u,t+k}^{[\ell]}}\,|\,\mathcal{F}_{u,k})\exp(ij\theta)$ where $\mathcal{F}_{u,k-1}:=\mathcal{F}_{\lfloor uN-n/2\rfloor+k-1}$. Then

[TABLE]

where

[TABLE]

Proof.

For fixed $u$, the result comes from [LW10, Proposition 2]. Since here we have assumed $\operatorname*{\mathrm{sup}\vphantom{\mathrm{infsup}}}_{u}\mathbb{E}X_{u,k}^{4}<\infty$, following the proof of [LW10, Proposition 2], the upper bound also holds uniformly over $u$. ∎

Lemma A.11**.**

Suppose that $\mathbb{E}X_{k}=0$, $\operatorname*{\mathrm{sup}\vphantom{\mathrm{infsup}}}_{u}\mathbb{E}X_{k}^{4}<\infty$, and $\operatorname{GMC}(2)$ holds, then

1. (1)

We have that

[TABLE]

uniformly on $(u_{1},u_{2},\theta_{1},\theta_{2})$ such that either $(u_{1},u_{2})\in\mathcal{U}^{2}$ or $(\theta_{1},\theta_{2})\in\Theta^{2}$ where $\mathcal{U}^{2}=\{(u_{1},u_{2}):\frac{n}{2N}\leq u_{1}\leq u_{2}\leq 1-\frac{n}{2N},|u_{1}-u_{2}|\geq\frac{n}{N}(1-1/(\log B_{n})^{2})\}$ and $\Theta^{2}=\{(\theta_{1},\theta_{2}):0\leq\theta_{1}<\theta_{2}\leq\pi-B_{n}^{-1}(\log B_{n})^{2},|\theta_{1}-\theta_{2}|\geq B_{n}^{-1}(\log B_{n})^{2}\}$. 2. (2)

For $\alpha_{n}>0$ with $\limsup\alpha_{n}<1$, we have that

[TABLE]

uniformly on $(u_{1},u_{2},\theta_{1},\theta_{2})$ such that either $(u_{1},u_{2})\in\mathcal{U}^{2}$ or $(\theta_{1},\theta_{2})\in\bar{\Theta}^{2}$ where $\mathcal{U}^{2}=\{(u_{1},u_{2}):\frac{n}{2N}\leq u_{1}\leq u_{2}\leq 1-\frac{n}{2N},|u_{1}-u_{2}|\geq\frac{n}{N}(1-1/(\log B_{n}))\}$ and $\bar{\Theta}^{2}=\{(\theta_{1},\theta_{2}):B_{n}^{-1}(\log B_{n})^{2}\leq\theta_{1}<\theta_{2}\leq\pi-B_{n}^{-1}(\log B_{n})^{2},|\theta_{1}-\theta_{2}|\geq B_{n}^{-1}\}$. 3. (3)

We have that

[TABLE]

uniformly on $\{(u,\theta):B_{n}^{-1}(\log B_{n})^{2}\leq\theta\leq\pi-B_{n}^{2}(\log B_{n})^{2},\frac{n}{2N}<u<1-\frac{n}{2N}\}$.

Proof.

Throughout this proof, we write $\tilde{X}_{k}^{[\ell]}$ as $\tilde{X}_{k}$ and $\tilde{X}_{u,i}^{[\ell]}$ as $\tilde{X}_{u,i}$ for simplicity.

1. (1)

By Lemma A.9 we approximate $g_{n}-\mathbb{E}g_{n}$ first by $\tilde{g}_{n}-\mathbb{E}\tilde{g}_{n}$. Then by Lemma A.10, we approximate $\tilde{g}_{n}-\mathbb{E}\tilde{g}_{n}$ by $M_{n}(u,\theta)$, where $M_{n}(u,\theta)=\sum_{t=1}^{n}D_{t}(u,\theta)^{*}\sum_{j=1}^{t-1}\alpha_{n,j-t}D_{j}(u,\theta)$. Then it is suffices to show that $|\mathbb{E}[M_{n}(u_{1},\theta_{1})-M_{n}^{*}(u_{1},\theta_{1})][M_{n}(u_{2},\theta_{2})-M_{n}^{*}(u_{2},\theta_{2})]|\leq C\frac{nB_{n}}{(\log B_{n})^{2}}$ and $|\mathbb{E}[M_{n}(u_{1},\theta_{1})+M_{n}^{*}(u_{1},\theta_{1})][M_{n}(u_{2},\theta_{2})+M_{n}^{*}(u_{2},\theta_{2})]|\leq C\frac{nB_{n}}{(\log B_{n})^{2}}$. We only prove the first inequality here, since the other inequality can be proved similarly. Define

[TABLE]

Since the martingale differences $\{D_{t}(u,\theta)\}$ are uncorrelated but not independent, we further define $N_{n}(u,\theta)=\sum_{t=1}^{n}D_{t}(u,\theta)^{*}\sum_{j=1}^{t-\ell-1}\alpha_{n,j-t}D_{j}(u,\theta)$, then $\|M_{n}(u,\theta)-N_{n}(u,\theta)\|=\mathcal{O}(\sqrt{n\ell})$ and $|r_{n}(u_{1},\theta_{1},u_{2},\theta_{2})|\leq|\tilde{r}_{n}(u_{1},\theta_{1},u_{2},\theta_{2})|+\mathcal{O}(\sqrt{n\ell(nB_{n})}+\sqrt{n\ell(B_{n}^{2})})$, where

[TABLE]

Since $\ell=\lfloor n^{\gamma}\rfloor$ where $\gamma$ is small enough, it suffices to show that $\tilde{r}_{n}(u_{1},\theta_{1},u_{2},\theta_{2})=\mathcal{O}(nB_{n}/(\log B_{n})^{2})$. Now we substitute $N_{n}(u,\theta)=\sum_{t=1}^{n}D_{t}(u,\theta)^{*}\sum_{j=1}^{t-\ell-1}\alpha_{n,j-t}D_{j}(u,\theta)$ to $\tilde{r}_{n}(u_{1},\theta_{1},u_{2},\theta_{2})$.

If $\theta_{1}\neq\theta_{2}$ and $u_{1}=u_{2}$, we have that

[TABLE]

Now it suffices to show that

[TABLE]

Since $|\theta_{1}-\theta_{2}|\geq B_{n}^{-1}(\log B_{n})^{2}$, using $1+2\sum_{k=1}^{n}\cos(k\theta)=\sin((n+1)\theta/2)/\sin(\theta/2)\leq 1/\sin(\theta/2)$, $\sin(x)=\Theta(x)$ when $x\to 0$, and denoting $j=t-s$, we have that

[TABLE]

If $\theta_{1}=\theta_{2}$ but $u_{1}\neq u_{2}$, using Eq. 189 and $n-N|u_{1}-u_{2}|\leq n/(\log B_{n})^{2}$, we have that

[TABLE] 2. (2)

When $\theta_{1}\neq\theta_{2}$, using [WN67, Lemma 3.2(ii)] with the assumption on the continuity of $a(\cdot)$ in Theorem 5.3, we have that

[TABLE]

If $\theta_{1}=\theta_{2}$ and $u_{1}\neq u_{2}$ then

[TABLE] 3. (3)

Since $\|D_{t}(u,\theta)\|^{2}=\sum_{j=\ell}^{\ell}\mathbb{E}(\tilde{X}_{u,t}\tilde{X}_{u,t+j})\exp(ij\theta)$, we have that

[TABLE]

∎

Lemma A.12**.**

Let $X_{1},\dots,X_{m}$ be independent mean zero $d$-dimensional random vectors such that $|X_{i}|\leq M$. If the underlying probability space is rich enough, one can define independent normally distributed mean zero random vectors $V_{1},\dots,V_{m}$ such that the covariance matrices of $V_{i}$ and $X_{i}$ are equal, for all $1\leq i\leq m$; furthermore

[TABLE]

Proof.

See [EM97, Fact 2.2]. ∎

Lemma A.13**.**

If $X$ and $Y$ have a bi-variate normally distributed distribution with expectations [math], unit variances, and correlation coefficient $r$, then

[TABLE]

uniformly for all $r$ such that $|r|\leq\delta$, for all $0<\delta<1$.

Proof.

See [Ber62, Lemma 2]. ∎

A.8. Proof of Lemma 8.9

By Markov’s inequality, we have that

[TABLE]

By Lemma A.9, $\mathbb{E}|g_{n}(u,\theta_{i})-\tilde{g}_{n}(u,\theta_{i})|=o(n^{1+\gamma}\rho^{\lfloor n^{\gamma}\rfloor})$ uniformly on $u$ and $\theta_{i}$. Since $D_{n}=B_{n}C_{n}$ is polynomial of $n$, the GMC assumption guarantees

[TABLE]

A.9. Proof of Lemma 8.10

Lemma A.14**.**

Let $X_{i},i=1,\dots,n$ be an arbitrary sequence of real-valued random variables with finite mean and variance. Then

[TABLE]

Proof.

See [Ave85, Theorem 2.1]. ∎

In this proof, we write $\tilde{X}_{k}^{[\ell]}$ and $\tilde{X}_{u,i}^{[\ell]}$ as $\tilde{X}_{k}$ and $\tilde{X}_{u,i}$ for simplicity. First of all, since $a(\cdot)$ has bounded support $[-1,1]$, we only need to consider the case that $|s-k|\leq B_{n}$. Furthermore, let $\alpha_{*}$ be an upper bound of $\alpha_{n,i}$ uniformly over $u$. Then, we have that

[TABLE]

Next, we show that the first term of the right hand side of Eq. 195 satisfies

[TABLE]

Similar arguments yield the same result for the second term of the right hand side of Eq. 195. Note that

[TABLE]

Applying Lemma A.14 and using $\ell$-independence and Hölder’s inequality, we have that uniformly on $u$

[TABLE]

Furthermore, uniformly on $u$, we also have that

[TABLE]

By the assumptions $p>4$ and $(p-1)\alpha>3/4$, we have that

[TABLE]

since we have assumed $C_{n}^{1/2}=o[(nB_{n})^{\alpha(p-1)-\frac{1}{2}}]$. Next, uniformly on $u$, the second term of the right hand side of Eq. 195 satisfies that

[TABLE]

Furthermore, uniformly on $u$, we have that

[TABLE]

Overall, we have that

[TABLE]

since we have assumed $C_{n}=o(B_{n}^{1+2\alpha(p-2)}n^{-2-2\gamma})$.

A.10. Proof of Lemma 8.11

We prove this lemma by first showing that

[TABLE]

and then showing

[TABLE]

To show Eq. 204, we note that $\{V_{j}\}$ are independent. Applying Lemma A.6, we have that

[TABLE]

Similar to the proof of Lemma 8.12, one can show $\sum_{j=1}^{k_{n}+1}\mathbb{E}V_{j}^{2}=\mathcal{O}(n^{1+\gamma}B_{n}^{1-\beta})$. Therefore, by choosing $\gamma$ close to zero and $Q$ large enough, we have that

[TABLE]

for any $c>0$. For the other term

[TABLE]

we apply Lemma A.7 with $M=(nB_{n})^{\alpha}$, $k=B_{n}+\ell$, $\ell=\lfloor n^{\gamma}\rfloor$ and $y=(\log B_{n})^{2}$, which yields

[TABLE]

where the second term of the right hand side is $\mathcal{O}(n^{-c})$ by choosing $Q$ large enough. Since $\alpha<1/4$ and $|\bar{X}_{i,\ell}|<(nB_{n})^{\alpha}$ almost surely, the last term of the right hand side converges to zero almost surely if

[TABLE]

which can be satisfied by choosing $\gamma$ close enough to zero. Therefore, by choosing $Q$ large enough so that $\mathcal{O}(C_{n}B_{n}n^{-c})=o(1)$ (Note that this only requires $C_{n}=o(n^{c}B_{n}^{-1})$ for some $c$, which is always satisfied when $C_{n}$ is polynomial of $n$), we have that

[TABLE]

which implies Eq. 204.

To prove Eq. 205, note that

[TABLE]

Therefore, other than using $p_{n}=B_{n}^{1+\beta}$ instead of $q_{n}=B_{n}+\ell$, the proof of Eq. 205 is essentially the same as the proof of Eq. 204.

A.11. Proof of Lemma 8.12

By Lemma A.3, we have that

[TABLE]

Note that $U_{j}=\sum_{k\in H_{j}}(\bar{Y}_{k,\ell}-\mathbb{E}\bar{Y}_{k,\ell})$, we first divide $\sum_{k\in H_{j}}(\bar{Y}_{k,\ell}-\mathbb{E}\bar{Y}_{k,\ell})$ into $\ell$ sums of sub-sequences. Note that $\bar{Y}_{k,\ell}=\bar{X}_{k,\ell}\sum_{s=1}^{k-1}\alpha_{n,k-s}\bar{X}_{s,\ell}=\bar{X}_{k,\ell}\sum_{s=\operatorname*{\mathrm{max}\vphantom{\mathrm{infsup}}}(1,k-B_{n})}^{k-1}\alpha_{n,k-s}\bar{X}_{s,\ell}$. Thus, one can get $\|\bar{U}_{j}\|^{2}=\mathcal{O}(\ell B_{n}^{2})$. Then using $\ell=\mathcal{O}(n^{\gamma})$ and $k_{n}=\lfloor n/(p_{n}+q_{n})\rfloor=\mathcal{O}(n/B_{n}^{1+\beta})$, one can get $\sum_{j=1}^{k_{n}}\mathbb{E}\bar{U}_{j}^{2}=\mathcal{O}(n^{1+\gamma}B_{n}^{1-\beta})=o(nB_{n})$ by choosing $\gamma$ and $\beta$ such that $n^{\gamma}B_{n}^{-\beta}=o(1)$.

Finally, we have that

[TABLE]

since $\log C_{n}+\log B_{n}=o(\log B_{n})^{4}$ when $C_{n}$ and $B_{n}$ are polynomials of $n$.

A.12. Proof of Lemma 8.13

1. (1)

We first show that for $|i_{1}-i_{2}|\geq(\log B_{n})^{2}/B_{n}$ or $|u_{1}-u_{2}|\geq\frac{n}{N}(1-1/(\log B_{n})^{2})$, we have that

[TABLE]

Note that $\sum_{j}\bar{U}_{j}$ can be approximated by $\bar{g}_{n}$. This is because according to the proof of Lemma 8.11, we have that

[TABLE]

Next, we can approximate $\bar{g}_{n}$ by $\tilde{g}_{n}$. This is because by Lemma 8.10 we have that

[TABLE]

Finally, we only need to show that

[TABLE]

which has been proved in Lemma A.11(i). 2. (2)

For convenience, we assume $\int a^{2}(t)\mathrm{d}t=1$. Select $d$ distinct tuples $(\theta_{i_{1}},u_{i}),i=1,\dots,d$ that $(\log B_{n})^{2}\leq i_{1}\leq\dots\leq i_{d}\leq B_{n}-(\log B_{n})^{2}$ and $u_{i}\in\mathcal{U},i=1,\dots,d$. Let $\bm{W}_{n}=\sum_{j=1}^{k_{n}}W_{j}$ where

[TABLE]

Note that by Lemma A.11(iii), we have that

[TABLE]

Together with Eq. 215, we have that

[TABLE]

Then we approximate $\bm{W}_{n}$ by $\bm{W}^{\prime}_{n}=\sum_{j=1}^{k_{n}}W^{\prime}_{j}$ using Lemma A.12, where $\{W^{\prime}_{j}\}$ are independent centered normally distributed random vectors. Then by Lemma A.12, we have $\textrm{Cov}(W_{j})=\textrm{Cov}(W^{\prime}_{j})$, for $1\leq j\leq k_{n}$, and

[TABLE]

Therefore, we have that

[TABLE] 3. (3)

Next, for $z=(z_{1},\dots,z_{d})$, we define the minimum of $\{z_{i}\}$ by $|z|_{d}:=\operatorname*{\mathrm{min}\vphantom{\mathrm{infsup}}}_{1\leq i\leq d}\{z_{i}\}$. Then we show that

[TABLE]

uniformly on distinct tuples of $\{(u_{j},\theta_{i_{j}}),j=1,\dots,d:(\log B_{n})^{2}\leq j_{1}\leq\cdots\leq j_{d}\leq B_{n}-(\log B_{n})^{2},\frac{n}{2N}<u_{j}<1-\frac{n}{2N}\}$ such that for any two tuples $(u_{j_{1}},\theta_{i_{j_{1}}})$ and $(u_{j_{2}},\theta_{i_{j_{2}}})$, if $u_{j_{1}}=u_{j_{2}}$ then $|\theta_{i_{j_{1}}}-\theta_{i_{j_{2}}}|\geq(\log B_{n})^{2}/B_{n}$; if $\theta_{i_{j_{1}}}=\theta_{i_{j_{2}}}$ then $|u_{j_{1}}-u_{j_{2}}|\geq\frac{n}{N}(1-1/(\log B_{n})^{2})$.

According to Eq. 222, we have that

[TABLE]

From Eq. 223, we have that

[TABLE]

so that for a standard normally distributed $R^{d}$-valued random vector, $\tilde{W}$, the tail probability of $\frac{\textrm{Cov}^{1/2}(\bm{W}^{\prime}_{n})}{\sqrt{nB_{n}}}\tilde{W}-2\pi\bm{I}_{d}\tilde{W}$ satisfies that

[TABLE]

Putting together the above results we can use $2\pi|\tilde{W}|_{d}$ (recall that we defined the minimum of $\{z_{i}\}$ by $|z|_{d}:=\operatorname*{\mathrm{min}\vphantom{\mathrm{infsup}}}_{1\leq i\leq d}\{z_{i}\}$) instead of $\frac{|\bm{W}^{\prime}_{n}|_{d}}{\sqrt{nB_{n}}}$ to bound the tail probability of $\frac{|\bm{W}_{n}|_{d}}{\sqrt{nB_{n}}}$:

[TABLE]

Using the following approximation of tail probability of a standard normally distributed random variable $Z$,

[TABLE]

we can get that

[TABLE]

Then we have shown that

[TABLE]

Similarly, using Lemma A.13 and Lemma A.11(ii), we can also have that

[TABLE]

for some $\delta>0$, uniformly on distinct tuples of $\{(u_{j},\theta_{i_{j}}),j=1,\dots,d:(\log B_{n})^{2}\leq j_{1}\leq\cdots\leq j_{d}\leq B_{n}-(\log B_{n})^{2},\frac{n}{2N}<u_{j}<1-\frac{n}{2N}\}$ such that for any two tuples $(u_{j_{1}},\theta_{i_{j_{1}}})$ and $(u_{j_{2}},\theta_{i_{j_{2}}})$with $j_{1}\leq j_{2}$, if $u_{j_{1}}=u_{j_{2}}$ then if $\theta_{i_{j_{1}}}=\operatorname*{\mathrm{min}\vphantom{\mathrm{infsup}}}_{j}\theta_{i_{j}}$ then $|\theta_{i_{j_{1}}}-\theta_{i_{j_{2}}}|\geq B_{n}^{-1}$; otherwise $|\theta_{i_{j_{1}}}-\theta_{i_{j_{2}}}|\geq(\log B_{n})^{2}/B_{n}$; if $\theta_{i_{j_{1}}}=\theta_{i_{j_{2}}}$ then $|u_{j_{1}}-u_{j_{2}}|\geq\frac{n}{N}(1-1/(\log B_{n})^{2})$. 4. (4)

Finally, we define

[TABLE]

and we show

[TABLE]

To this end, we define

[TABLE]

and

[TABLE]

Then by Bonferroni’s inequality, we have for every fixed $k$ and $u$ that

[TABLE]

Next following the proof of [Wat54, Theorem] and [WN67, Theorem 3.3] based on Eq. 224 and Eq. 232, we can show that

[TABLE]

as $n\to\infty$. As shown in [Wat54, pp.799], with Eq. 224 and Eq. 232, when $n\to\infty$, we have that

[TABLE]

Therefore, we have shown that

[TABLE]

Finally, we use the above techniques again to show

[TABLE]

which means we only need to show

[TABLE]

Letting $y_{n}^{2}/4\pi^{2}=2\log B_{n}+2\log C_{n}-\log(\pi\log B_{n}+\pi\log C_{n})+x$, as in Eq. 224, we have that

[TABLE]

A.13. Proof of Remark 4.2

First of all, by the assumption $\operatorname{GMC}(2)$

[TABLE]

Next, by the SLC condition, we know $r(u,k)$ is Lipschitz. Together with the Lipschitz condition of $\tau(\cdot)$, we have that

[TABLE]

Since $r(u,k)$ is twice continuously differentiable with respect to $u$, we have that

[TABLE]

Furthermore, since $\tau(\cdot)$ is an even function

[TABLE]

Therefore, we have that

[TABLE]

Therefore, by the locally quadratic property of $a(\cdot)$ at [math], we have that

[TABLE]

Then, using the fact that if $\theta\notin\{0,\pi\}$, we know that

[TABLE]

Then, for fixed $\theta\notin\{0,\pi\}$, we have that

[TABLE]

If $\operatorname*{\mathrm{sup}\vphantom{\mathrm{infsup}}}_{u}\sum_{k\in\mathbb{Z}}|r(u,k)|k^{2}<\infty$ and $B_{n}=o(n)$, then

[TABLE]

Finally, $B_{n}=o(N^{1/3})$ implies $\mathcal{O}(B_{n}/N)=o(1/B_{n}^{2})$. Also, $n=o(N^{2/3})$ and $B_{n}=o(N^{1/3})$ implies $\mathcal{O}(B_{n}^{2}n^{2}/N^{2})=o(1)$.

A.14. Proof of Theorem 5.5

For simplicity, we denote $\delta_{u,n}$ as $\delta_{u}$ and $\delta_{\theta,n}$ as $\delta_{\theta}$. First, we write

[TABLE]

Then by continuity we have that

[TABLE]

Letting $\hat{g}_{n}(u,u_{i},\theta,\theta_{j}):=\hat{f}_{n}(u,\theta)-\hat{f}_{n}(u_{i},\theta_{j})$, it suffices to show that

[TABLE]

Note that

[TABLE]

Then we can write

[TABLE]

Since by Theorem 4.1, we have that

[TABLE]

Therefore, the following result completes the proof.

Lemma A.15**.**

If $\delta_{u}=\mathcal{O}(\frac{n}{N(\log n)^{\alpha}})$ and $\delta_{\theta}=\mathcal{O}(\frac{1}{B_{n}(\log n)^{\alpha}})$ for some $\alpha>0$, then

[TABLE]

Proof.

See Section A.15. ∎

A.15. Proof of Lemma A.15

First, we pick any $(u_{0},\theta_{0})$ such that $|u_{0}-u|\leq\delta_{u}$ and $|\theta_{0}-\theta|\leq\delta_{\theta}$. Then

[TABLE]

Using $\tau\left(\frac{i-\lfloor u_{0}N\rfloor}{n}\right)=\tau\left(\frac{i-\lfloor uN\rfloor}{n}\right)+\mathcal{O}\left(\frac{\delta_{u}N}{n}\right)$, we have that

[TABLE]

Note that $\exp(\sqrt{-1}k\theta_{0})=\exp(\sqrt{-1}k\theta)[\exp(\sqrt{-1}k(\theta_{0}-\theta))]$ and $\cos(k\theta_{0})=\cos(k\theta)\cos(k(\theta_{0}-\theta))-\sin(k\theta)\sin(k(\theta_{0}-\theta))$. Therefore, we have that

[TABLE]

where we have used the fact that the GMC condition implies $\sum_{k=0}^{B_{n}}kr(u,k)=\mathcal{O}(\sum_{k=0}^{B_{n}}k\rho^{k})=\mathcal{O}(B_{n})$. Note that we have assumed that $f(u,\theta)>f_{*}>0$ uniformly over $u$ and $\theta$, so we can write $\mathcal{O}_{\mathbb{P}}(B_{n}\delta_{\theta})=(B_{n}\delta_{\theta})\mathcal{O}_{\mathbb{P}}(\hat{f}_{n}(u,\theta))$. Therefore, we have that

[TABLE]

which implies that

[TABLE]

In order to make it equal to $o_{\mathbb{P}}(1)$, we only need $\delta_{u}=o(n/N)$ and $\delta_{\theta}=o(1/B_{n})$. Therefore, choosing $\alpha>0$, $\delta_{u}=\mathcal{O}\left(\frac{n}{N(\log n)^{\alpha}}\right)$ and $\delta_{\theta}=\mathcal{O}\left(\frac{1}{B_{n}(\log n)^{\alpha}}\right)$ is sufficient.

A.16. Proof of Remark 2.9

By the triangle inequality and Hölder’s inequality, we have that

[TABLE]