The exact contents of citations.db main_text.text for this paper — one flattened LaTeX string, title through conclusion, appendix excluded, unmodified except for removing email addresses. This is what our citation measures are computed over.
128,029 characters
Theory of Evolutionary Spectra for Heteroskedasticity and Autocorrelation Robust Inference in Possibly Misspecified and Nonstationary Models
\pagebreak{}
\setcounter{page}{0}
\raggedbottom
\title{\textbf{\Large{}Theory of Evolutionary Spectra for Heteroskedasticity
and Autocorrelation Robust Inference in Possibly Misspecified and
Nonstationary Models}\textbf{}\thanks{This paper is based on the first chapter of my doctoral dissertation
at Boston University. I am grateful to Pierre Perron for his continuous
support and advice. I thank Whitney Newey and Tim Vogelsang for helpful
discussions. I also thank Marine Carrasco, Giuseppe Cavaliere, Bin
Chen, Andrew Chesher, Jean-Jacques Forneron, Massimo Franchi,
Raffaella Giacomini, Liudas Giraitis, Jesus Gonzalo, Chris Hansen,
Jungbin Hwang, Hiroaki Kaido, Ivana Komunjer, Oliver Linton, Adam
McCloskey, Zhongjun Qu, Myung Seo, Julius Vainora, Rasmus Varneskov,
Daniel Whilem, Dacheng Xiu, Yohei Yamamoto, and seminar participants
at BC-BU Econometrics Workshop, Cambridge, Chicago Booth, Connecticut,
KU Leuven, Georgetown and UCL for comments. }}
\maketitle
\begin{abstract}
{\footnotesize{}The literature on heteroskedasticity and autocorrelation
robust (HAR) inference is extensive but its usefulness relies on
stationarity of the relevant process, say $V_{t}$, usually a function
of the data and estimated model residuals. Yet, a large body of work
shows widespread evidence of various forms of nonstationarity in the
latter. Also, many testing problems are such that $V_{t}$ is stationary
under the null hypothesis but nonstationary under the alternative.
In either case, the consequences are possible size distortions and,
especially, a reduction in power which can be substantial (e.g., non-monotonic
power), since all such estimates are based on weighted sums of the
sample autocovariances of $V_{t}$, which are inflated. We propose
HAR inference methods valid under a broad class of nonstationary processes,
labelled Segmented Local Stationary, which possess a spectrum that
varies both over frequencies and time. It is allowed to change either
slowly and continuously and/or abruptly at some time points, thereby
encompassing most nonstationary models used in applied work. We introduce
a double kernel estimator (DK-HAC) that applies a smoothing over both
lagged autocovariances and time. The optimal kernels and bandwidth
sequences are derived under a mean-squared error criterion. The data-dependent
bandwidths rely on the ``plug-in'' approach using approximating
parametric models having time-varying parameters estimated with standard
methods applied to local data. Our method yields tests with good size
and power under both stationary and nonstationary, thereby encompassing
previous methods. In particular, the power gains are achieved without
notable size distortions, the exact size being as good as those delivered
by the best fixed-$b$ approach, when the latter works well.}{\footnotesize\par}
\end{abstract}
\indent {\bf{JEL Classification}}: C12, C13, C18, C22, C32, C51\\
\noindent {\bf{Keywords}}: Fixed-$b$, HAC standard errors, HAR, Long-run variance, Nonstationarity, Misspecification, Outliers, Segmented locally stationary.
\onehalfspacing
\thispagestyle{empty}
\pagebreak{}
\section{Introduction}
The literature on heteroskedasticity and autocorrelation robust (HAR)
inference is extensive and quite mature by now. For concreteness,
consider the linear model where $x_{t}$ is a vector of regressors
and $e_{t}$ is an unobservable disturbance, which can be serially
correlated. It is now common practice to use OLS and correct the
standard errors. This entails the estimation of the covariance matrix
(referred to as the long-run variance, LRV) of $V_{t}=x_{t}e_{t}$
or ($2\pi$ times) the spectral density of $V_{t}$ at frequency zero
when the latter is stationary (of course, in general, the relevant
process $V_{t}$ can be generated from a more complex model; e.g.,
a moment condition in a GMM context). Early important contributions
in econometrics are \citeauthor{newey/west:87} (\citeyear{newey/west:87};
\citeyear{newey/west:94}) and \citet{andrews:91} who proposed heteroskedasticity
and autocorrelation consistent (HAC) estimators with some optimal
properties. This approach aims at devising good estimate of the LRV
of $V_{t}$. An alternative method foregoes that aim and concentrates
on having a test with a pivotal non-normal limit distribution that
is obtained through an inconsistent estimate of the LRV of $V_{t}$
that keeps the bandwidth at a fixed fraction of the sample size. This
is the so-called fixed-$b$ HAR inference initiated by \citet{Kiefer/vogelsang/bunzel:00}
and \citeauthor{Kiefer/vogelsang:02} (\citeyear{Kiefer/vogelsang:02};
\citeyear{kiefer/vogelsang:05}). The drawback of this approach is
that the limit distribution changes depending on the context and critical
values are to be obtained numerically on a case by case basis. The
literature since then has focused on providing various refinements,
mostly to have tests having exact sizes closer to the nominal level.\footnote{The fixed-$b$ or post-HAC literature is vast; see \citet{dou:18},
\citet{lazarus/lewis/stock:17}, Lazarus et al. \citeyearpar{lazarus/lewis/stock/watson:18},
\nocite{Goncalves/vogelsang:11} \citet{dejong/davidson:00}, \citet{ibragimov/muller:10},
\citet{jansson:04}, \citeauthor{muller:07} (\citeyear{muller:07};
\citeyear{mueller:14}), \citet{phillips:05}, \citet{politis:11},
\citet{preinerstorfer/potscher:16}, \citet{potscher/preinerstorfer:18},
\citet{rho/vogelsang:2020}, \citet{robinson:98}, \citeauthor{sun:13}
\citeyearpar{sun:13,sun:14,sun:14a}, \citet{sun/phillips/jin:08}
and \citet{zhang/shao:13}.}
Most of this literature relies on stationarity with exception of the
consistency results in \citet{newey/west:87} and of some results
in \citet{andrews:91} which, however, do not provide accurate approximations.
Yet, another strand of the literature has argued convincingly that
the processes governing economic data $\{x_{t}\}$ and the errors
in the relevant regressions $\left\{ e_{t}\right\} $ are nonstationary.\footnote{By nonstationary we mean non-constant moments. As in the literature,
we consider processes whose sum of absolute autocovariances is finite.
That is, we rule out processes with unbounded second moments (e.g.,
unit root).} This can occur for several reasons: changes in the moments of $x_{t}$
induced by changes in the model parameters that govern the data {[}cf.
\citet{perron:89}, \citet{stock/watson:96} and the surveys of \citet{ng/wright:13}
and \citet{giacomini/rossi:15}{]}; changes in the moments of $e_{t}$
(think about the Great Moderation with the decline in variance for
many macroeconomic variables or the effects of the COVID-19 pandemic);
smooth changes in the distributions governing either processes that
arise from transitory dynamics; and so on. All these induces nonstationarity
in $\{V_{t}\}$, which then makes $\mathbb{E}(V_{t}V'_{t-k})$ depend
on both $k$ and $t$. Furthermore, even if the data and primitive
shocks $\{e_{t}\}$ are stationary, many HAR testing problems are
such that the relevant process $\{V_{t}\}$ is stationary under the
null hypothesis but is affected by changes in means (or other forms
of nonstationarity) under the alternative. This occurs, for instance,
when using tests involving structural breaks based on estimating the
model under the null hypothesis; e.g., popular tests for forecast
evaluation {[}e.g., \citet{diebold/mariano:95}{]}, tests for forecast
instability {[}cf. \citet{casini_CR_Test_Inst_Forecast}, \citet{giacomini/rossi:09}
and \citet{perron/yamamoto:18}{]}, tests for structural change {[}cf.
\citet{casini/perron_Oxford_Survey} and \citet{perron:06}{]}. When
standardized by classical HAC estimators such tests may suffer from
issues such as non-monotonic power, i.e., power that goes to zero
as the alternative gets farther away from the null value. Various
forms of misspecication and/or nonstationarity generate low frequency
contamination and make the series or residuals appear much more persistent.
As a consequence, HAC standard errors are too large and when used
as normalizing factors of test statistics, the tests lose power {[}see
\citet{casini/perron_Low_Frequency_Contam_Nonstat:2020} for formal
details{]}.\footnote{A partial list of works that present evidence of power issues with
HAR inference is \citet{altissimo/corradi:2003}, \citeauthor{casini/perron_Oxford_Survey}
(\citeyear{casini/perron_Oxford_Survey}, \citeyear{casini/perron_Lap_CR_Single_Inf},
\citeyear{casini/perron_SC_BP_Lap}), \citet{chan:2020}, \citet{crainiceanu/vogelsang:07},
\citet{juhl/xiao:09}, \citet{kim/perron:09}, \citet{martins/perron:16},
\citet{perron/yamamoto:18}, \citet{shao/zhang:2010}, \citet{vogeslang:99},
\citet{xu:2013}, \citet{zhang/lavitas:2018}.} This applies even more forcefully to the fixed-$b$ type methods
and to the recent refinements by \citet{lazarus/lewis/stock:17} and
Lazarus et al. \citeyearpar{lazarus/lewis/stock/watson:18}, since
they involve more lagged autocovariances (or long bandwidths) and,
hence, larger contaminations.
This points to the importance of extending the methods for HAR inference
so that they have the correct size and good power even under nonstationarity.
This is the aim of the paper. We first develop a theoretical framework
under which to analyze the statistical properties of our suggested
estimate. We introduce a class of nonstationary processes which possess
a spectrum that varies both over frequencies and time, thereby encompassing
the nonstationary models used in applied work. We work in an infill
asymptotic setting akin to the one used in nonparametric regression
{[}cf. \citet{robinson:89}{]}. For a process $V_{t}$, its spectrum
at frequency $\omega$ and time $u=t/T$, denoted by $f\left(u,\,\omega\right)$,
is allowed to change slowly yet continuously as well as to change
abruptly in $u$ at a finite number of time points; the latter allows
for structural breaks in the spectrum of $V_{t}$. We label this class
as Segmented Locally Stationary (SLS). It is related to the locally
stationary processes introduced by \citet{dahlhaus:96} that have
the characterizing property of behaving as a stationary process in
a small neighborhood of $u$. This is achieved via smoothness of
$f\left(u,\,\omega\right)$ in $u$. By allowing discontinuities across
some segments, we can deal with relevant features such as structural
change, regime switching-type and threshold models {[}cf. \citet{bai/perron:98},
\citeauthor{casini/perron_CR_Single_Break} (\citeyear{casini/perron_Oxford_Survey},
\citeyear{casini/perron_Lap_CR_Single_Inf}, \citeyear{casini/perron_SC_BP_Lap}
and \citeyear{casini/perron_CR_Single_Break}), \citet{hamilton:89}
and \citet{hansen:00ecma}{]}. The SLS class extends some of the analysis
of \citet{dahlhaus:96} to processes having a more general time-varying
spectrum.\footnote{A few authors used a notion of local stationarity that allows for
breaks {[}see, e.g., \citet{dahlhaus:2009} and \citet{last/shumway:08}{]}.
However, none of these works was concerned with HAR inference. \citet{dahlhaus:2009}
presented some results for local spectral density estimation and required
smoothness (see Example 4.2 there). \citet{last/shumway:08} considered
testing for change-points in a locally stationary series which under
the alternative hypothesis results in a piecewise locally stationary
series. Furthermore, our notion of SLS processes and related framework
are more general from a theoretical standpoint since we provide precise
definitions and establish theoretical results about the identification
of the local spectral density when there are discontinuities. } Our framework is of independent interest and can be useful in many
contexts in econometrics if one is interested in deriving the properties
of estimators or inference under nonstationarity.
Under this framework, we introduce a double kernel HAC (DK-HAC) estimator
in order to flexibly account for nonstationarity and we show that
it is robust to low frequency contamination and other misspecifications.
This entails an extension of the classical HAC estimators since in
addition to the usual smoothing procedure over lagged autocovariances,
it applies a second smoothing over time for each lagged autocovariance,
involving a second kernel and bandwidth. If $\{V_{t}\}$ is Segmented
Locally Stationary, $\mathbb{E}(V_{t}V'_{t-k})$ changes smoothly
in $t$, as long as $t$ is away from the change-points in the spectrum
$f\left(t/T,\,\omega\right)$. Thus, the smoothing over time yields
good estimates for the time path of $\mathbb{E}(V_{t}V'_{t-k})$ for
all $k$. We determine the optimal kernels and optimal values for
both bandwidth sequences under a mean-squared error (MSE) criterion.
We establish new MSE bounds that show how nonstationarity affects
the bias-variance trade-off and are more informative than previously
established MSE bounds. We use them to construct data-dependent bandwidths
relying on the ``plug-in'' approach. Unlike \citet{andrews:91},
our candidate parametric models have time-varying parameters which
can be estimated by applying standard methods to local data, akin
to using rolling regressions. The procedure depends on three elements:
the bandwidths for the smoothing over autocovariances and over time,
and a block size to separate the regimes. In this paper, we consider
a sequential bandwidth selection procedure by first deriving the optimal
bandwidth for smoothing over time, then conditioning on this to obtain
the optimal bandwidth for smoothing over autocovariances.
The DK-HAC estimators can result in HAR tests that are oversized when
there is high temporal dependence in the data, a well-known problem
for all methods, though for ours these distortions are relatively
minor compared to, e.g., the methods of \citet{newey/west:87} and
\citet{andrews:91}. Still, in order to improve the size control of
HAR tests, \citet{casini/perron_PrewhitedHAC}, using the theory of
this paper, propose a nonparametric nonlinear VAR prewhitened DK-HAC
estimators. This form of prewhitening differs from those discussed
previously {[}e.g., \citet{andrews/monahan:92} and \citet{rho/shao:13}{]}
in that it accounts explicitly for nonstationarity. HAR tests based
on prewhitened DK-HAC estimators have size control competitive to
fixed-$b$ HAR tests when the latter work well (i.e., under stationarity).
Notably, non-prewhitened and prewhitened DK-HAC have excellent power
properties even when existing HAR tests have serious issues with power.
\subsection*{Comparison to Existing Literature}
There are two main approaches to HAR inference differing on whether
the LRV estimator is consistent or not. The classical approach relies
on consistency, which results in HAC estimators {[}cf. \citeauthor{newey/west:87}
\citeyearpar{newey/west:87,newey/west:94}, \citet{andrews:91} and
\citet{hansen:92ecma}{]}, and on bandwidths chosen via MSE criterion.
Inference is standard because HAR tests follow asymptotically standard
distributions. The researcher then uses corrected standard errors
and asymptotic critical values. It was shown early that classical
HAC standard errors can result in oversized tests when there is substantial
temporal dependence. This stimulated a second approach based on inconsistent
LRV estimators that keep the bandwidth at some fixed fraction of the
sample size {[}cf. \citet{Kiefer/vogelsang/bunzel:00}{]}. Because
of the inconsistency, inference is nonstandard and HAR tests do not
asymptotically follow standard distributions. Critical values are
to be obtained numerically. Long bandwidths/fixed-$b$ methods require
stationarity and reduce the oversize problem of HAR tests. The bandwidth
choice is often based on testing-oriented criteria {[}e.g., \citet{sun/phillips/jin:08}{]}.
Our approach falls in the first category; we propose HAC estimators
and standard HAR inference.
We now compare in detail our approach to the existing literature.
We believe that a fair comparison has to consider the following four
criteria: (1) applicability to general HAR inference tests; (2) size
of HAR tests; (3) power of HAR tests; (4) theoretical validity under
stationarity/nonstationarity. In terms of (1), it is clear that the
HAC and DK-HAC estimators are generally and immediately applicable
to any HAR inference test and that they are simple to use in practice.
This also explains why the classical HAC estimators have become the
standard practice in econometrics. Methods that rely on long bandwidths/fixed-$b$
do not share the same property. They are not generally applicable
to HAR inference tests because a researcher would first need to derive
a new asymptotic non-standard fixed-$b$ distribution. This can be
unfeasible in non-standard testing problems {[}e.g., tests for parameter
instability, etc.{]}. Turning to (2), all existing HAR inference
tests are known to be oversized when there is strong serial dependence.
However, fixed-$b$ HAR tests (or versions thereof) are less oversized
than other tests based on the classical HAC estimators. The stronger
is the temporal dependence the larger is the difference in size between
the two approaches. Our prewhitened DK-HAC estimators are competitive
with fixed-$b$ HAR tests in controlling the size. Moving to (3),
prewhitened and non-prewhitened DK-HAC estimators have excellent power
under either stationarity or nonstationarity whereas existing methods
have serious problems with power under nonstationarity or under nonstationary
alternative hypotheses. These problems result in non-monotonic power
and little or no power in relevant circumstances especially in HAR
tests outside the stable linear regression model. Fixed-$b$ or long
bandwidths methods suffer most from these problems. Finally, turning
to (4), our method like the classical HAC approach is valid under
nonstationarity whereas methods using long bandwidths/fixed-$b$ are
only valid under stationarity {[}cf. \citet{casini_fixed_b_erp}{]}.
It should be mentioned that the fixed-$b$ approach is shown to achieve
(pointwise) higher-order refinements under stationarity while the
MSE-based optimality of the HAC or DK-HAC estimators pertains only
to the first-order but it holds under nonstationarity.\footnote{\citet{preinerstorfer/potscher:16} pointed out some limitations of
these approaches because optimality does not hold uniformly over DGPs.
They showed negative non-asymptotic results about size and power of
HAR tests when one allows for all correlation structures corresponding
to stationary Gaussian AR(1) processes. }
Recently, \citet{lazarus/lewis/stock:17} and Lazarus et al. \citeyearpar{lazarus/lewis/stock/watson:18}
made some progress on the applicability of fixed-$b$ methods. They
showed that the $t$-test using a LRV estimator based on equally-weighted
cosine (EWC) under fixed-$b$ asymptotics can achieve a $t$-distribution
with the degrees of freedom depending on the bandwidth choice.\footnote{\citet{hwang/sun:2017} proposed a modification to the trinity of
test statistics in the two-step GMM setting and showed that the modified
test statistics are asymptotically $F$ distributed under fixed-$b$
asymptotics and stationarity. However, following the logic of \citet{lazarus/lewis/stock:17},
procedures such as e.g., bootstrap-based autocorrelation robust tests,
modification of the test statistic, etc., are not as simple as HAC-based
inference and therefore may find less traction than HAC-based inference
in empirical work. } However, \citet{casini/perron_PrewhitedHAC} showed that EWC is oversized
relative to the original fixed-$b$ of \citet{Kiefer/vogelsang/bunzel:00}
and to the prewhitened DK-HAC when there is strong dependence. \citet{lazarus/lewis/stock:17}
relied on the Neyman-Pearson Lemma or simply ``apple-to-apple comparison''
to compare HAR tests. This is certainly a reasonable criterion. The
lemma suggests to compare the power of tests that have an empirical
size no greater than the significance level. However, the Neyman-Pearson
Lemma alone does not suffice to find the ``best'' test in this context
because all tests are oversized when there is strong dependence. Indeed,
it does not even apply in this context. We face a trade-off between
size and power. Our proposed method is competitive with the existing
methods which is least oversized {[}i.e., original fixed-$b$ of \citet{Kiefer/vogelsang/bunzel:00}{]}
and has excellent power even when existing HAR tests do not have any.
We believe that our method strikes a good balance with respect to
criteria (1)-(4).
Our approach is different from methods based on subsampling of $t$-statistics
{[}see, e.g., \citet{ibragimov/muller:10}{]}. The latter rely on
splitting the sample in subsamples and estimating the model within
each subsample. Under the assumption that the estimates from the subsamples
are asymptotically independent, the test statistic based on an average
of estimates across subsamples follows asymptotically a $t$-distribution.
In terms of point (1) above, this approach is not general enough compared
to the HAC/DK-HAC approach because this changes the test statistic
and its asymptotic distribution. Also, subsampling test statistics
which are not $t$-tests can be challenging/unfeasible and would require
at best extra work in general HAR inference contexts to derive the
new distribution. Finally, our simulation experience (not reported)
suggests that this method suffers from the same finite-sample issues
about size and power as the fixed-$b$ methods.
\subsection*{Related Work}
This paper is part of a set of papers on HAR inference by the author
and collegues. The current paper provides the core theoretical and
empirical elements that are used in all other papers, which can be
viewed as providing extensions or refinements. \citet{casini/perron_PrewhitedHAC}
used our theoretical framework to derive minimax MSE bounds for LRV
estimation that are sharper than previously established and extended
some of our theoretical results to general nonstationarity. As a
finite-sample refinement, they also developed a new prewhitening procedure
robust to nonstationarity for DK-HAC estimators. Even though the latter
procedure is included in our simulations, we established the corresponding
theoretical results in a separate paper because of the extent of the
work needed in the analysis. \citet{casini/perron_Low_Frequency_Contam_Nonstat:2020}
showed analytically that the poor finite-sample performance of existing
LRV-based HAR tests under nonstationarity and misspecification is
induced by low frequency contamination. Belotti et al. \citeyearpar{belotti/casini/catania/grassi/perron_HAC_Sim_Bandws}
used our theoretical framework to propose alternative data-dependent
bandwidths for DK-HAC estimators that are optimal under a global MSE
criterion. \citet{casini:change-point-spectra} considered change-point
detection in time series with evolutionary spectra. Initially, it
was intended to be used in the method suggested in this paper to improve
the finite-sample size and power properties, given that we work with
segmented locally stationary processes. However, it turns out that
the current method to select the blocks (see Section \ref{subsec Choice-of nT})
is able to handle even abrupt structural change.
\bigskip{}
The remainder of the paper is organized as follows. Section \ref{Section: Statistical Enviromnent}
introduces the statistical setting and the new HAC estimator. Section
\ref{Section HAC-Estimation-with Predetermined} presents consistency,
rates of convergence and the asymptotic MSE results for the DK-HAC
estimators. Asymptotically optimal kernels and bandwidths are derived
in Section \ref{Section Optimal-Kernels-and}. A data-dependent method
for choosing the bandwidths and its asymptotic properties are discussed
in Section \ref{Section Data-Dependent-Bandwidths}. Section \ref{Section Monte Carlo}
presents a Monte Carlo study. Section \ref{Section Conclusions}
concludes the paper. The supplemental materials {[}cf. \citet{casini_hac_supp}
and an additional supplement not for publication{]} contain some implementation
details and all mathematical proofs. The code to implement our methods
is provided in $\mathrm{\mathsf{Matlab}}$, $\mathrm{\mathrm{\mathsf{R}}}$
and $\mathrm{\mathrm{\mathsf{Stata}}}$ languages through a $\mathrm{\mathsf{Github}}$
repository.
\section{The Statistical Environment\label{Section: Statistical Enviromnent}}
To motivate our approach, consider the linear regression model estimated
by least-squares (LS): $y_{t}=x'_{t}\beta_{0}+e_{t}$ $(t=1,\ldots,\,T)$,
where $\beta_{0}\in\Theta\subset\mathbb{R}^{p}$, $y_{t}$ is an observation
on the dependent variable, $x_{t}$ is a $p$-vector of regressors
and $e_{t}$ is an unobserved disturbance. The LS estimator is given
by $\widehat{\beta}=(X'X)^{-1}X'Y$, where $Y=(y_{1},\ldots,\,y_{T})'$
and $X=(x_{1},\ldots,\,x_{T})'$. Classical inference about $\beta_{0}$
requires estimation of $\mathrm{Var}(\sqrt{T}(\widehat{\beta}-\beta_{0}))$
where
\begin{align*}
\mathrm{Var}(\sqrt{T}(\widehat{\beta}-\beta_{0})) & \triangleq\mathbb{E}\left[\left(T^{-1}\sum_{t=1}^{T}x_{t}x'_{t}\right)^{-1}T^{-1}\sum_{s=1}^{T}\sum_{t=1}^{T}e_{s}x_{s}(e_{t}x{}_{t})'\left(T^{-1}\sum_{t=1}^{T}x_{t}x'_{t}\right)^{-1}\right],
\end{align*}
where ``$\triangleq$'' is used for definitional equivalence. Consistent
estimation of $\mathrm{Var}(\sqrt{T}(\widehat{\beta}-\beta_{0}))$
relies on consistent estimation of $\mathrm{lim}_{T\rightarrow\infty}T^{-1}\sum_{s=1}^{T}\sum_{t=1}^{T}$
$\mathbb{E}(e_{s}x_{s}(e_{t}x{}_{t})')$. More generally, one needs
a consistent estimate of $J\triangleq\mathrm{lim}_{T\rightarrow\infty}J_{T}$
where $J_{T}=T^{-1}\sum_{s=1}^{T}\sum_{t=1}^{T}\mathbb{E}(V_{s}(\beta_{0})$
$V_{t}(\beta_{0})')$ with $V_{t}(\beta)$ being a random $p$-vector
for each $\beta\in\Theta$. For the linear regression model, $V_{t}(\beta)=(y_{t}-x'_{t}\beta)x_{t}$.
For HAR tests outside the regression model $V_{t}(\beta)$ takes different
forms.\footnote{If one suspects that $\beta_{0}$ may not be constant, one can use
appropriate tests for parameter instability. However, our discussions
and methods still apply because these tests are HAR inference tests
and one needs a LRV estimate of $J$ based on the appropriate $V_{t}.$ } Hence, our problem is to estimate $J$ when $\{V_{t}\}$ is a Segmented
Locally Stationary process, as defined in Section \ref{Subsec Segmented-Locally-Stationary}.\footnote{\citet{casini/perron_PrewhitedHAC} extended the results to the case
where $\left\{ V_{t}\right\} $ is generally nonstationary (i.e.,
$\left\{ V_{t}\right\} $ is a sequence of unconditionally heteroskedastic
random variables). } Such estimate can then be used to conduct HAR inference in the usual
way using the theory developed below. By a change of variables,
$J_{T}$ can be rewritten as
\begin{align*}
J_{T} & =\sum_{k=-T+1}^{T-1}\Gamma_{T,k},\qquad\mathrm{where}\qquad\Gamma_{T,k}=\begin{cases}
T^{-1}\sum_{t=k+1}^{T}\mathbb{E}(V_{t}V'_{t-k}) & \mathrm{for\,}k\geq0\\
T^{-1}\sum_{t=-k+1}^{T}\mathbb{E}(V_{t+k}V'_{t}) & \mathrm{for\,}k<0
\end{cases},
\end{align*}
and $V_{t}=V_{t}\left(\beta_{0}\right)$. The rest of this section
is structured as follows. In Section \ref{Subsec Segmented-Locally-Stationary}
we introduce a new class of nonstationary time series that we use
as the underlying framework for our theoretical analysis. Section
\ref{Subsection HAC-Estimation} presents the DK-HAC estimator. We
adopt the following notational conventions. The $j$th element of
a vector $x$ is indicated by $x^{\left(j\right)}$ while the $\left(j,\,l\right)$th
element of a matrix $X$ is indicated by $X^{\left(j,\,l\right)}$.
$\mathrm{tr}(\cdot)$ denotes the trace and $\otimes$ denotes the
tensor product. The $p^{2}\times p^{2}$ matrix $C_{pp}$ is a commutation
matrix that transforms $\mathrm{vec}\left(A\right)$ into $\mathrm{vec}\left(A'\right)$,
i.e., $C_{pp}=\sum_{j=1}^{p}\sum_{l=1}^{p}\iota_{j}\iota_{l}'\otimes\iota_{l}\iota_{j}'$,
where $\iota_{j}$ is the $j$th elementary $p$-vector. $\lambda_{\max}\left(A\right)$
denotes the largest eigenvalue of $A$. $W$ and $\widetilde{W}$
are used for $p^{2}\times p^{2}$ weight matrices. $\mathbb{C}$
is used for the set of complex numbers. $\overline{A}$ is used for
the complex conjugate of $A\in\mathbb{C}$. Let $0=\lambda_{0}<\ldots<\lambda_{m+1}=1$.
A function $G\left(\cdot,\,\cdot\right):\,\left[0,\,1\right]\times\mathbb{R}\rightarrow\mathbb{C}$
is said to be piecewise (Lipschitz) continuous with $m+1$ segments
if it is (Lipschitz) continuous within each segment (e.g., it is piecewise
Lipschitz continuous if for each $j=1,\ldots,\,m+1$ it satisfies
$\sup_{u\neq v}\left|G\left(u,\,\omega\right)-G\left(v,\,\omega\right)\right|\leq K\left|u-v\right|$
for any $\omega\in\mathbb{R}$ with $\lambda_{j-1}<u,\,v\leq\lambda_{j}$
for some $K<\infty$). We define $G_{j}\left(u,\,\omega\right)=G\left(u,\,\omega\right)$
for $\lambda_{j-1}<u\leq\lambda_{j}$. If we say piecewise Lipschitz
continuous with index $\vartheta>0$, then the above inequality is
replaced by $\sup_{u\neq v}\left|G\left(u,\,\omega\right)-G\left(v,\,\omega\right)\right|\leq K\left|u-v\right|^{\vartheta}$.
A function $G\left(\cdot,\,\cdot\right):\,\left[0,\,1\right]\times\mathbb{R}\rightarrow\mathbb{C}$
is said to be left-differentiable at $u_{0}$ if $\partial G\left(u_{0},\omega\right)/\partial_{-}u\triangleq\lim_{u\rightarrow u_{0}^{-}}\left(G\left(u_{0},\,\omega\right)-G\left(u,\,\omega\right)\right)/\left(u_{0}-u\right)$
exists $\forall\omega\in\mathbb{R}$.
\subsection{\label{Subsec Segmented-Locally-Stationary}Segmented Locally Stationary
Processes}
Suppose $\left\{ V_{t}\right\} _{t=1}^{T}$ is defined on an abstract
probability space $\left(\Omega,\,\mathscr{F},\,\mathbb{P}\right)$,
where $\Omega$ is the sample space, $\mathscr{F}$ is the $\sigma$-algebra
and $\mathbb{P}$ is a probability measure. In order to introduce
a framework to analyze time series models with a time-varying spectrum
it is necessary to introduce an infill asymptotic setting whereby
we rescale the original discrete time horizon $\left[1,\,T\right]$
by dividing each $t$ by $T.$ Letting $u=t/T$ and $T\rightarrow\infty,$
this defines a new time scale $u\in\left[0,\,1\right]$ which we interpret
as saying that as $T\rightarrow\infty$ we observe more and more realizations
of $V_{t}$ close to time $t$, i.e., we observe the rescaled process
$V_{Tu}$ on the interval $\left[u-\varepsilon,\,u+\varepsilon\right]$,
where $\varepsilon>0$ is a small number.
In order to define a general class of nonstationary processes, we
shall start from processes that have a time-varying spectral representation
specified by:
\begin{align}
V_{t,T} & =\mu\left(t/T\right)+\int_{-\pi}^{\pi}\exp\left(i\omega t\right)A\left(t/T,\,\omega\right)d\xi\left(\omega\right),\label{Eq. 2.2 Definition V_t}
\end{align}
where $i\triangleq\sqrt{-1}$, $\mu\left(t/T\right)$ is the trend
function, $A\left(t/T,\,\omega\right)$ is the transfer function and
$\xi\left(\omega\right)$ is some stochastic process whose properties
are specified below. Observe that this representation is similar to
the spectral representation of stationary processes {[}see \citet{anderson:71},
\citet{brillinger:75}, \citet{hannan:70} and \citet{priestley:85}{]}.
We shall see that the main difference is that $A\left(t/T,\,\omega\right)$
and $\mu\left(t/T\right)$ are not constant in $t$.\footnote{In HAR inference, a minimal assumption on $V_{t}$ under the null
hypothesis is that it has zero mean (i.e., $\mu\left(t/T\right)=0$
for all $t$). However, in this subsection we allow for arbitrary
$\mu\left(t/T\right)$ so as to introduce a general framework, also
applicable under various alternative hypotheses in both within and
outside the regression model. } \citet{dahlhaus:96} used the time-varying spectral representation
to define the so-called locally stationary processes which are characterized,
broadly speaking, by smoothness conditions on $\mu\left(\cdot\right)$
and $A\left(\cdot,\,\cdot\right)$. Locally stationary processes have
been used widely in both statistics and economics, though in the latter
field they are best known as time-varying parameter processes {[}see,
e.g., \citet{cai:07} and \citet{chen/hong:12}{]}. The smoothness
restrictions exclude many prominent models that account for time variation
in the parameters. For example, structural change and regime switching-type
models do not belong to this class because parameter changes occur
suddenly at a particular point in time. We propose a class of nonstationarity
processes which allow both continuous and discontinuous changes in
the parameters. Stationarity and local stationarity are recovered
as special cases.
\begin{defn}
\label{Definition Segmented-Locally-Stationary}A sequence of stochastic
processes $\{V_{t,T}\}_{t=1}^{T}$ is called Segmented Locally Stationary
(SLS)\textbf{ }with $m_{0}+1$ regimes, transfer function $A^{0}$
and trend $\mu$, if there exists a representation
\begin{align}
V_{t,T} & =\mu_{j}\left(t/T\right)+\int_{-\pi}^{\pi}\exp\left(i\omega t\right)A_{j,t,T}^{0}\left(\omega\right)d\xi\left(\omega\right),\qquad\qquad\left(t=T_{j-1}^{0}+1,\ldots,\,T_{j}^{0}\right),\label{Eq. Spectral Rep of SLS}
\end{align}
for $j=1,\ldots,\,m_{0}+1$, where by convention $T_{0}^{0}=0$ and
$T_{m_{0}+1}^{0}=T$ and the following holds:
(i) $\xi\left(\omega\right)$ is a stochastic process on $\left[-\pi,\,\pi\right]$
with $\overline{\xi\left(\omega\right)}=\xi\left(-\omega\right)$
and
\begin{align*}
\mathrm{cum}\left\{ d\xi\left(\omega_{1}\right),\ldots,\,d\xi\left(\omega_{r}\right)\right\} & =\varphi\left(\sum_{j=1}^{r}\omega_{j}\right)g_{r}\left(\omega_{1},\ldots,\,\omega_{r-1}\right)d\omega_{1}\ldots d\omega_{r},
\end{align*}
where $\mathrm{cum}\left\{ \cdot\right\} $ is the cumulant of
$r$th order, $g_{1}=0,\,g_{2}\left(\omega\right)=1$, $\left|g_{r}\left(\omega_{1},\ldots,\,\omega_{r-1}\right)\right|\leq M_{r}<\infty$
and $\varphi\left(\omega\right)=\sum_{j=-\infty}^{\infty}\delta\left(\omega+2\pi j\right)$
is the period $2\pi$ extension of the Dirac delta function $\delta\left(\cdot\right)$.
(ii) There exists a constant $K>0$ and a piecewise continuous function
$A:\,\left[0,\,1\right]\times\mathbb{R}\rightarrow\mathbb{C}$ such
that, for each $j=1,\ldots,\,m_{0}+1$, there exists a $2\pi$-periodic
function $A_{j}:\,(\lambda_{j-1}^{0},\,\lambda_{j}^{0}]\times\mathbb{R}\rightarrow\mathbb{C}$
with $A_{j}\left(u,\,-\omega\right)=\overline{A_{j}\left(u,\,\omega\right)}$,
$\lambda_{j}^{0}\triangleq T_{j}^{0}/T$ and for all $T,$
\begin{align}
A\left(u,\,\omega\right)=A_{j}\left(u,\,\omega\right) & \,\mathrm{\,for\,}\,\lambda_{j-1}^{0}<u\leq\lambda_{j}^{0},\label{Eq A(u) =00003D Ai}\\
\sup_{1\leq j\leq m_{0}+1}\sup_{T_{j-1}^{0}<t\leq T_{j}^{0},\,\omega}\left|A_{j,t,T}^{0}\left(\omega\right)-A_{j}\left(t/T,\,\omega\right)\right| & \leq KT^{-1}.\label{Eq. 2.4 Smothenss Assumption on A}
\end{align}
(iii) $\mu_{j}\left(t/T\right)$ is piecewise continuous.
\end{defn}
The smoothness properties of $A$ in $u$ guarantees that $V_{t,T}$
has a piecewise locally stationary behavior. Later we will require
additional smoothness properties for $A$.
\begin{example}
(i) Suppose $X_{t}$ is a stationary process with spectral representation
$X_{t}=\int_{-\pi}^{\pi}\exp\left(i\omega t\right)$ $A\left(\omega\right)d\xi\left(\omega\right),$
and $\mu,\,\sigma:\,\left[0,\,1\right]\rightarrow\mathbb{R}$ are
piecewise continuous. Then, $V_{t,T}=\mu_{j}\left(t/T\right)+\sigma_{j}\left(t/T\right)X_{t}$,
with $T_{j-1}^{0}<t\leq T_{j}^{0}\,(1\leq j\leq m_{0}+1)$ is a SLS
process with $m_{0}+1$ regimes where $A_{j,t,T}^{0}\left(\omega\right)=A_{j}\left(t/T,\,\omega\right)=\sigma_{j}\left(t/T\right)A\left(\omega\right)$.
Within each segment, $V_{t,T}$ is locally stationary. When $t=Tu$
is away from the change-points, as $T\rightarrow\infty$ more and
more realizations of $V_{Tu,T}$ with $u\in\left[u-\varepsilon,\,u+\varepsilon\right]$
are observed, that is, realizations with amplitude close to $\sigma_{j}\left(u\right)$
for the appropriate $j$.
(ii) Suppose $e_{t}$ is an $i.i.d.$ sequence and $V_{t,T}=\sum_{k=0}^{\infty}a_{j,k}\left(t/T\right)e_{t-k},\,T_{j-1}^{0}<t\leq T_{j}^{0}$
$(1\leq j\leq m_{0}+1)$. Then, $V_{t,T}$ is SLS with $A_{j,t,T}^{0}\left(\omega\right)=A_{j}\left(t/T,\,\omega\right)=\sum_{k=0}^{\infty}a_{j,k}\left(t/T\right)$
$\exp\left(-i\omega k\right).$
(iii) Autoregressive processes with time-varying coefficients, known
as TVAR, augmented with structural breaks are SLS. In this case,
we do not have the exact relationship $A_{j,t,T}^{0}\left(\omega\right)=A_{j}\left(t/T,\,\omega\right)$
but only the approximate relationship \eqref{Eq. 2.4 Smothenss Assumption on A}.
\end{example}
If there is only a single regime (i.e., $m_{0}=0$) then $V_{t,T}$
is locally stationary {[}cf. \citet{dahlhaus:96}{]}. If $\mu$ and
$A^{0}$ do not depend on $t$, then $V_{t,T}$ is stationary and
the spectral representation of stationary processes applies. However,
$m_{0}=0$ rules out structural change and regime switching models.
With $m_{0}>0$, we propose a framework where parameter variation
can occur either smoothly or abruptly, both being relevant for economic
data.\footnote{Some authors have used alternative notions of local stationarity that
allow for discontinuities (i.e., piecewise locally stationary) and
have established some results in other contexts which are not related
to HAR inference {[}see, e.g., \citet{dahlhaus:2009}, \citet{last/shumway:08}
and \citet{zhou:2013}{]}. In particular, our framework is more general
because we also define (and work with) the covariance between observations
belonging to different regimes whereas previous works considered only
the covariance between observations belonging to the same regime thereby
using smoothness which restricts the framework substantially.}
Let $\left\lfloor \cdot\right\rfloor $ denote the largest smaller
integer function and let $\mathcal{T}\triangleq\{T_{1}^{0},\,\ldots,\,T_{m_{0}}^{0}\}$.
We define the spectrum of $V_{t,T}$ in \eqref{Eq. 2.2 Definition V_t}
(for fixed $T$) as
\begin{align*}
f_{j,T}\left(u,\,\omega\right) & \triangleq\begin{cases}
\left(2\pi\right)^{-1}\sum_{s=-\infty}^{\infty}\mathrm{Cov}\left(V_{\left\lfloor Tu-3\left|s\right|/2\right\rfloor ,T},\,V_{\left\lfloor Tu-\left|s\right|/2\right\rfloor ,T}\right)\exp\left(-i\omega s\right), & Tu\in\mathcal{T},\,u=\lambda_{j}^{0}\\
\left(2\pi\right)^{-1}\sum_{s=-\infty}^{\infty}\mathrm{Cov}\left(V_{\left\lfloor Tu-s/2\right\rfloor ,T},\,V_{\left\lfloor Tu+s/2\right\rfloor ,T}\right)\exp\left(-i\omega s\right), & Tu\notin\mathcal{T},\,u\in(\lambda_{j-1}^{0},\,\lambda_{j}^{0})
\end{cases}
\end{align*}
with $A_{1,t,T}^{0}\left(\omega\right)=A_{1}\left(0,\,\omega\right)$
for $t<1$ and $A_{m_{0}+1,t,T}^{0}\left(\omega\right)=A_{m_{0}+1}\left(1,\,\omega\right)$
for $t>T$. Our definition coincides with the Wigner-Ville spectrum
{[}cf. \citet{martin/flandrin:85}{]} when there are no change-points
(i.e., $m_{0}=0$). Below we show that $f_{j,T}\left(u,\,\omega\right)$
tends in mean-squared to $f_{j}\left(u,\,\omega\right)\triangleq\left|A_{j}\left(u,\,\omega\right)\right|^{2}$
for $T_{j-1}^{0}/T<u=t/T\leq T_{j}^{0}/T$ which is the spectrum that
corresponds to the spectral representation. Therefore, we call $f_{j}\left(u,\,\omega\right)$
the time-varying spectral density matrix of the process.
\begin{assumption}
\label{Assumption Smothness of A}$A\left(u,\,\omega\right)$ is \textcolor{red}{
}piecewise Lipschitz continuous in the first component and uniformly
Lipschitz continuous in the second component, with index $\vartheta>1/2$
for both.
\end{assumption}
\begin{thm}
\label{Theorem 2.2 in Dal}Assume $V_{t,T}$ is Segmented Locally
Stationary with $m_{0}+1$ regimes and Assumption \ref{Eq. 2.2 Definition V_t}
holds. Then, for all $u\in\left(0,\,1\right)$, $\int_{-\pi}^{\pi}\sum_{j=1}^{m_{0}+1}\left|f_{j,T}\left(u,\,\omega\right)-f_{j}\left(u,\,\omega\right)\right|^{2}d\omega=o\left(1\right)$.
\end{thm}
Let $f\left(u,\,\omega\right)=f_{j}\left(u,\,\omega\right)$ if $Tu\in(T_{j-1}^{0},\,T_{j}^{0}]$
so as to suppress the subscript $j$ from $f$. It is well-known that,
even when $m_{0}=0,$ the spectral representation \eqref{Eq. Spectral Rep of SLS}
is not unique {[}cf. \citet{priestley:85}, Chapter 11.1{]}. A consequence
of Theorem \ref{Theorem 2.2 in Dal} is that $\{f_{j}\left(u,\,\omega\right)=\left|A_{j}\left(u,\,\omega\right)\right|^{2},\,j=1,\ldots,\,m_{0}+1\}$
is uniquely determined from the whole triangular array $\{V_{t,T}\}$.
For $Tu\notin\mathcal{T}$ with $T_{j-1}^{0}/T<u=t/T<T_{j}^{0}/T$,
only the realizations of $V_{t,T}$ in the time interval $u\in\left[u-n/T,\,u+n/T\right]$
with $n\rightarrow\infty$ contribute to $f_{j}\left(u,\,\omega\right)$.
Since this interval is fully contained in a segment $j$ where $A_{j}\left(u,\,\omega\right)$
is smooth, and given that the length of this interval tends to zero,
$V_{t,T}$ becomes ``asymptotically stationary'' on this interval.
The length of the interval in which $V_{t,T}$ can be considered
stationary is given by $n\ln n/T^{\vartheta}\rightarrow0$ . For $Tu\in\mathcal{T}$,
the arguments are different. Suppose $Tu=T_{j}^{0}$. The spectrum
$f_{j,T}\left(u,\,\omega\right)$ is defined in such a way that only
observations prior to $T_{j}^{0}$ are used in order to construct
an approximation to $f_{j}\left(u,\,\omega\right)$. Since the length
of this interval tends to zero and $A_{j}\left(u,\,\omega\right)$
is left-Lipschitz continuous, then those observations become ``asymptotically
stationary'' and thus provide the same information about $f_{j}\left(u,\,\omega\right)$.
Given $f\left(u,\,\omega\right),$ we can define the local covariance
of $V_{t,T}$ at rescaled time $u$ with $Tu\notin\mathcal{T}$ and
lag $k\in\mathbb{Z}$ as
\[
c\left(u,\,k\right)\triangleq\int_{-\pi}^{\pi}e^{i\omega k}f\left(u,\,\omega\right)d\omega.
\]
The same definition is also used when $Tu\in\mathcal{T}$ and $k\geq0$.
For $Tu\in\mathcal{T}$ and $k<0$ it is defined as $c\left(u,\,k\right)\triangleq\int_{-\pi}^{\pi}e^{i\omega k}A\left(u,\,\omega\right)A\left(u-k/T,\,-\omega\right)d\omega$.
\subsection{\label{Subsection HAC-Estimation}DK-HAC Estimation}
In model \eqref{Eq. Spectral Rep of SLS}, if $m_{0}=0$ and $A^{0}$
is constant in its first argument, then $\left\{ V_{t,T}\right\} $
is second-order stationary. Its spectral density matrix is then equal
to $f\left(\omega\right)\triangleq\left(2\pi\right)^{-1}\sum_{k=-\infty}^{\infty}\Gamma\left(k\right)e^{-i\omega k}$
where $\Gamma\left(k\right)\triangleq\mathbb{E}(V_{t,T}V'_{t-k,T})$.
When evaluated at frequency $\omega=0$ it plays a prominent role
because $\lim_{T\rightarrow\infty}J_{T}=2\pi f\left(0\right).$ Nonstationarity
implies that the spectral density is time-varying since $\mathbb{E}(V_{t}V'_{t-k})$
now depends on $k$ as well as on $t$. The SLS processes introduced
above accommodate this property because they have a time-varying spectrum
$f\left(u,\,\omega\right)$. Accordingly, we introduce the notation
$\Gamma_{u}\left(k\right)\triangleq\mathbb{E}(V_{Tu,T}V'_{Tu-k,T})$
where $u=t/T$. We show below that $\Gamma_{u}\left(k\right)=c\left(u,\,k\right)+O\left(T^{-1}\right)$
uniformly in $1\leq j\leq m+1$, $Tu\leq T_{j}^{0}$ and $k\in\mathbb{Z}.$
Under the rescaling $u=t/T,\,u\in\left[0,\,1\right]$, the limit of
$J_{T}$ for SLS processes is given by,
\begin{align*}
J\triangleq\lim_{T\rightarrow\infty}J_{T} & =\int_{0}^{1}c\left(u,\,0\right)du+\sum_{k=1}^{\infty}\int_{0}^{1}\left(c\left(u,\,k\right)+c\left(u,\,k\right)'\right)du.
\end{align*}
Using the definition of $f\left(u,\,\omega\right)$ it can be shown
that $J=2\pi\int_{0}^{1}f\left(u,\,0\right)du$. \citet{dahlhaus:2009}
discussed how to estimate $f\left(u,\,\omega\right)$ for the scalar
case under smoothness in both arguments using the smoothed local periodogram.
Our goals are to estimate $J$ using a time-domain method and to relax
the smoothness assumption in $u$. This is different from Dahlhaus'
work that considered local problems (i.e., estimation of $f\left(u,\,\omega\right)$
under smoothness) and not full-sample problems (i.e., estimation of
$J$). The class of estimators of $J$ relies on double kernel smoothing
over lags and time,
\begin{align*}
\widehat{J}_{T} & =\widehat{J}_{T}\left(b_{1,T},\,b_{2,T}\right)\triangleq\frac{T}{T-p}\sum_{k=-T+1}^{T-1}K_{1}\left(b_{1,T}k\right)\widehat{\Gamma}\left(k\right),\,\,\mathrm{with}\,\,\\
\widehat{\Gamma}\left(k\right) & \triangleq\frac{n_{T}}{T-n_{T}}\sum_{r=0}^{\left\lfloor \left(T-n_{T}\right)/n_{T}\right\rfloor }\widehat{c}_{T}\left(rn_{T}/T,\,k\right),
\end{align*}
where $K_{1}\left(\cdot\right)$ is a real-valued kernel in the class
$\boldsymbol{K}_{1}$ defined below, $b_{1,T}$ is a bandwidth sequence
discussed below, $n_{T}\rightarrow\infty$ satisfying the conditions
given below, and
\begin{align}
\widehat{c}_{T}\left(rn_{T}/T,\,k\right) & \triangleq\begin{cases}
\left(Tb_{2,T}\right)^{-1}\sum_{s=k+1}^{T}K_{2}^{*}\left(\frac{\left(\left(r+1\right)n_{T}-\left(s-k/2\right)\right)/T}{b_{2,T}}\right)\widehat{V}_{s}\widehat{V}'_{s-k}, & k\geq0\\
\left(Tb_{2,T}\right)^{-1}\sum_{s=-k+1}^{T}K_{2}^{*}\left(\frac{\left(\left(r+1\right)n_{T}-\left(s+k/2\right)\right)/T}{b_{2,T}}\right)\widehat{V}_{s+k}\widehat{V}'_{s}, & k<0
\end{cases},\label{eq: Def. chat}
\end{align}
with $K_{2}^{*}$ being a real-valued kernel and $b_{2,T}$ is a
bandwidth sequence discussed below. $\widehat{c}_{T}\left(u,\,k\right)$
is an estimate of the local autocovariance $c\left(u,\,k\right)$
of lag $k$ at time $u=rn_{T}/T$. Estimation of $c\left(u,\,k\right)$
for locally stationary processes was considered by \citet{dahlhaus:12}.
For positive semi-definiteness, it is necessary that $K_{2}^{*}$
takes the following form:
\[
K_{2}^{*}\left(\frac{\left(r+1\right)n_{T}-\left(s-k/2\right)}{Tb_{2,T}}\right)=\left(K_{2}\left(\frac{\left(r+1\right)n_{T}-s}{Tb_{2,T}}\right)K_{2}\left(\frac{\left(r+1\right)n_{T}-\left(s-k\right)}{Tb_{2,T}}\right)\right)^{1/2}\quad\mathrm{for}\,k\geq0,
\]
\[
K_{2}^{*}\left(\frac{\left(r+1\right)n_{T}-\left(s+k/2\right)}{Tb_{2,T}}\right)=\left(K_{2}\left(\frac{\left(r+1\right)n_{T}-s}{Tb_{2,T}}\right)K_{2}\left(\frac{\left(r+1\right)n_{T}-\left(s+k\right)}{Tb_{2,T}}\right)\right)^{1/2}\quad\mathrm{for}\,k<0.
\]
Setting $K_{2}\left(x\right)=(\int_{0}^{1}h\left(x\right)^{2}dx)^{-1}h\left(x+1/2\right)^{2}$
and $N_{T}=Tb_{2,T}$, we see that positive semi-definiteness requires
the use of a data taper $h\left(\cdot\right)$ with length $N_{T}$.
This follows because we need each $\widehat{V}_{t}$ $\left(t=1,\ldots,\,T\right)$
to be assigned the same weight across different $k$ for any given
$r$. Then, letting $\widehat{V}_{t}^{\ifmmode\begingroup\defbold{bold}
\text{\ifx\math@versionbold\bfseries\fi\textdegree}\endgroup\else\textdegree\fi}=\ensuremath{\left(K_{2}\left(\left(\left(r+1\right)n_{T}-t\right)/Tb_{2,T}\right)\right)^{1/2}\widehat{V}_{t}}$
we can use the same arguments as in \citet{andrews:91} applied now
to $\widehat{V}_{t}^{\ifmmode\begingroup\defbold{bold}
\text{\ifx\math@versionbold\bfseries\fi\textdegree}\endgroup\else\textdegree\fi}$ to show that $J_{T}$
is positive semi-definite for the appropriate choice of $K_{1}$.
The estimator $\widehat{J}_{T}$ involves two kernels: $K_{1}$ smooths
the lagged sample autocovariances, akin to the classical HAC estimators,
while $K_{2}$ applies smoothing over time. The factor $T/\left(T-p\right)$
is an optional small-sample degrees of freedom adjustment. In Section
\ref{Section HAC-Estimation-with Predetermined}-\ref{Section Optimal-Kernels-and},
we consider estimators $\widehat{J}_{T}$ for which $b_{1,T}$ and
$b_{2,T}$ are given sequences. In Section \ref{Section Data-Dependent-Bandwidths},
we consider adaptive estimators $\widehat{J}_{T}$ for which $b_{1,T}$
and $b_{2,T}$ are data-dependent. Observe that the optimal $b_{2,T}$
actually depends on the properties of $\left\{ V_{t,T}\right\} $
in any given block. Since the order of $b_{2,T}\left(\cdot\right)$
is the same across blocks, we omit this notation for the developments
of the asymptotic results. However, when we determine the data-dependent
estimate of $b_{2,T}\left(\cdot\right)$, we will estimate $b_{2,T}\left(rn_{T}/T\right)$
for each $r$. We consider the following class of kernels {[}cf.
\citet{andrews:91}{]},
\begin{align}
\boldsymbol{K}_{1} & =\{K_{1}\left(\cdot\right):\,\mathbb{R}\rightarrow\left[-1,\,1\right]:\,K_{1}\left(0\right)=1,\,K_{1}\left(x\right)=K_{1}\left(-x\right),\,\forall x\in\mathbb{R}\label{Eq. (2.6) K1 Kernel class}\\
& \quad{\textstyle \int\nolimits _{-\infty}^{\infty}}K_{1}^{2}\left(x\right)dx<\infty,\,K_{1}\left(\cdot\right)\,\mathrm{is\,continuous\,at\,0\,and\,at\,all\,but\,finite\,numbers\,of\,points}\}.\nonumber
\end{align}
Examples of kernels in $\boldsymbol{K}_{1}$ include the Truncated,
Bartlett, Parzen, Quadratic Spectral (QS) and Tukey-Hanning kernel.
We shall show below that the QS kernel has certain optimality properties:
\begin{align*}
K_{1}^{\mathrm{QS}}\left(x\right) & =\frac{25}{12\pi^{2}x^{2}}\left(\frac{\sin\left(6\pi x/5\right)}{6\pi x/5}-\cos\left(6\pi x/5\right)\right).
\end{align*}
\section{\label{Section HAC-Estimation-with Predetermined}HAC Estimation
with Predetermined Bandwidths}
In Section \ref{Subsection: Estimation-of-Local} we present some
asymptotic properties of $\widehat{c}\left(\cdot,\,\cdot\right)$.
We use them in Section \ref{Subsec: Large-Sample-Results-on} in order
to establish consistency, rate of convergence and MSE properties of
predetermined bandwidths HAC estimators. Let $\widetilde{J}_{T}$
denote the pseudo-estimator identical to $\widehat{J}_{T}$ but based
on $\{V_{t,T}\}=\{V_{t,T}(\beta_{0})\}$ rather than on $\{\widehat{V}_{t,T}\}=\{V_{t,T}(\widehat{\beta})\}$.
We first require some smoothness of $A\left(u,\,\cdot\right)$ in
$u$.
\begin{assumption}
\label{Assumption Smothness of A (for HAC)}(i) $\left\{ V_{t,T}\right\} $
is a mean-zero SLS process with $m_{0}+1$ regimes; (ii) $A\left(u,\,\omega\right)$
is \textcolor{red}{ }twice continuously differentiable in $u$ at
all $u\neq\lambda_{j}^{0}$ $(j=1,\ldots,\,m_{0}+1)$ with uniformly
bounded derivatives $\left(\partial/\partial u\right)A\left(u,\,\cdot\right)$
and $\left(\partial^{2}/\partial u^{2}\right)A\left(u,\,\cdot\right)$,
and Lipschitz continuous in the second component with index $\vartheta=1$;
(iii) $\left(\partial^{2}/\partial u^{2}\right)A\left(u,\,\cdot\right)$
is Lipschitz continuous at all $u\neq\lambda_{j}^{0}$ $(j=1,\ldots,\,m_{0}+1)$;
(iv) $A\left(u,\,\omega\right)$ is twice left-differentiable in $u$
at $u=\lambda_{j}^{0}$, $(j=1,\ldots,\,m_{0}+1)$ with uniformly
bounded derivatives $\left(\partial/\partial_{-}u\right)A\left(u,\,\cdot\right)$
and $\left(\partial^{2}/\partial_{-}u^{2}\right)A\left(u,\,\cdot\right),$
and has piecewise Lipschitz continuous derivative $\left(\partial^{2}/\partial_{-}u^{2}\right)A\left(u,\,\cdot\right)$.
\end{assumption}
We also need to impose conditions on the temporal dependence of $V_{t}=V_{t,T}$.
Let
\begin{align*}
\kappa_{V,t}^{\left(a,b,c,d\right)}\left(u,\,v,\,w\right) & \triangleq\kappa^{\left(a,b,c,d\right)}\left(t,\,t+u,\,t+v,\,t+w\right)-\kappa_{\mathscr{N}}^{\left(a,b,c,d\right)}\left(t,\,t+u,\,t+v,\,t+w\right)\\
& \triangleq\mathbb{E}(V_{t}^{\left(a\right)}V_{t+u}^{\left(b\right)}V_{t+v}^{\left(c\right)}V_{t+w}^{\left(d\right)})-\mathbb{E}(V_{\mathscr{N},t}^{\left(a\right)}V_{\mathscr{N},t+u}^{\left(b\right)}V_{\mathscr{N},t+v}^{\left(c\right)}V_{\mathscr{N},t+w}^{\left(d\right)}),
\end{align*}
where $\left\{ V_{\mathscr{N},t}\right\} $ is a Gaussian sequence
with the same mean and covariance structure as $\left\{ V_{t}\right\} $.
$\kappa_{V,t}^{\left(a,b,c,d\right)}\left(u,\,v,\,w\right)$ is the
time-$t$ fourth-order cumulant of $(V_{t}^{\left(a\right)},\,V_{t+u}^{\left(b\right)},\,V_{t+v}^{\left(c\right)},$
$\,V_{t+w}^{\left(d\right)})$ while $\kappa_{\mathscr{N}}^{\left(a,b,c,d\right)}(t,\,t+u,$
$\,t+v,\,t+w)$ is the time-$t$ centered fourth moment of $V_{t}$
if $V_{t}$ were Gaussian.
\begin{assumption}
\label{Assumption A - Dependence}(i) $\sum_{k=-\infty}^{\infty}\sup_{u\in\left[0,\,1\right]}$
$\left\Vert c\left(u,\,k\right)\right\Vert <\infty$, $\sum_{k=-\infty}^{\infty}\sup_{u\in\left[0,\,1\right]}\left\Vert \left(\partial^{2}/\partial u^{2}\right)c\left(u,\,k\right)\right\Vert <\infty$
and $\sum_{k=-\infty}^{\infty}\sum_{j=-\infty}^{\infty}\sum_{l=-\infty}^{\infty}\sup_{u\in\left[0,\,1\right]}|\kappa_{V,\left\lfloor Tu\right\rfloor }^{\left(a,b,c,d\right)}$
$\left(k,\,j,\,l\right)|<\infty$ for all $a,\,b,\,c,\,d\leq p$.
(ii) For all $a,\,b,\,c,\,d\leq p$ there exists a function $\widetilde{\kappa}_{a,b,c,d}:\,\left[0,\,1\right]\times\mathbb{Z}\times\mathbb{Z}\times\mathbb{Z}\rightarrow\mathbb{R}$
such that $\sup_{u\in\left(0,\,1\right)}|\kappa_{V,\left\lfloor Tu\right\rfloor }^{\left(a,b,c,d\right)}\left(k,\,s,\,l\right)$
$-\widetilde{\kappa}_{a,b,c,d}\left(u,\,k,\,s,\,l\right)|\leq KT^{-1}$
for some constant $K$; the function $\widetilde{\kappa}_{a,b,c,d}\left(u,\,k,\,s,\,l\right)$
is twice differentiable in $u$ at all $u\neq\lambda_{j}^{0}$, $(j=1,\ldots,\,m_{0}+1)$
with uniformly bounded derivatives $\left(\partial/\partial u\right)\widetilde{\kappa}_{a,b,c,d}\left(u,\cdot,\cdot,\cdot\right)$
and $\left(\partial^{2}/\partial u^{2}\right)\widetilde{\kappa}_{a,b,c,d}\left(u,\cdot,\cdot,\cdot\right)$,
and twice left-differentiable in $u$ at $u=\lambda_{j}^{0}$ $(j=1,\ldots,\,m_{0}+1)$
with uniformly bounded derivatives $\left(\partial/\partial_{-}u\right)\widetilde{\kappa}_{a,b,c,d}\left(u,\cdot,\cdot,\cdot\right)$
and $\left(\partial^{2}/\partial_{-}u^{2}\right)\widetilde{\kappa}_{a,b,c,d}$
$\left(u,\cdot,\cdot,\cdot\right),$ and piecewise Lipschitz continuous
derivative $\left(\partial^{2}/\partial_{-}u^{2}\right)\widetilde{\kappa}_{a,b,c,d}\left(u,\cdot,\cdot,\cdot\right)$.
\end{assumption}
If $\left\{ V_{t,T}\right\} $ is stationary then the cumulant condition
of Assumption \ref{Assumption A - Dependence}-(i) reduces to the
standard one used in the time series literature {[}see also Assumption
A in \citet{andrews:91}{]}. We do not require fourth-order stationarity
but only that the time-$t=Tu$ fourth order cumulant is locally constant
in a neighborhood of $u$. One can show that $\alpha$-mixing and
moment conditions imply that the cumulant condition of Assumption
\ref{Assumption A - Dependence} holds.
\subsection{\label{Subsection: Estimation-of-Local}Estimation of the Local Covariance}
Let $\widetilde{c}_{T}\left(u,\,k\right)$ denote the estimator that
uses $\left\{ V_{t,T}\right\} $. We consider the following
class of kernels:
\begin{align}
\boldsymbol{K}_{2} & =\{K_{2}\left(\cdot\right):\,\mathbb{R}\rightarrow\left[0,\,\infty\right],\,K_{2}\left(x\right)=K_{2}\left(1-x\right),\,{\textstyle \int}K_{2}\left(x\right)dx=1,\label{Eq. K2 Kernel class}\\
& \qquad\qquad K_{2}\left(x\right)=0,\,\mathrm{for\,}\,x\notin\left[0,\,1\right],\,K_{2}\left(\cdot\right)\,\mathrm{is\,continuous}\}.\nonumber
\end{align}
\begin{lem}
\label{Lemma Rate of Convergence of ctilde - SLS}Suppose that Assumption
\ref{Assumption Smothness of A (for HAC)}-\ref{Assumption A - Dependence}
hold. If $b_{2,T}\rightarrow0$ and\textbf{ }$Tb_{2,T}^{5}\rightarrow\eta\in(0,\,\infty)$,
then $\widetilde{c}_{T}\left(u_{0},\,k\right)-c\left(u_{0},\,k\right)=O_{\mathbb{P}}(\sqrt{Tb_{2,T}})$
for all $u_{0}\in\left(0,\,1\right)$.
\end{lem}
\subsection{\label{Subsec: Large-Sample-Results-on}Results on DK-HAC Estimation
with Predetermined Bandwidths}
Following \citet{parzen:57}, we define $K_{1,q}\triangleq\lim_{x\downarrow0}\left(1-K_{1}\left(x\right)\right)/\left|x\right|^{q}$
for $q\in[0,\,\infty);$ $q$ increases with the smoothness of $K_{1}\left(\cdot\right)$
with the largest value being such that $K_{1,q}<\infty$. When $q$
is an even integer, $K_{1,q}=-\left(d^{q}K_{1}\left(x\right)/dx^{q}\right)|_{x=0}/q!$
and $K_{1,q}<\infty$ if and only if $K_{1}\left(x\right)$ is $q$
times differentiable at zero. We define the index of smoothness
of $f\left(u,\,\omega\right)$ at $\omega=0$ by $f^{\left(q\right)}\left(u,\,0\right)\triangleq\left(2\pi\right)^{-1}\sum_{k=-\infty}^{\infty}\left|k\right|^{q}c\left(u,\,k\right)$,
for $q\in[0,\,\infty)$. If $q$ is even, then $f^{\left(q\right)}\left(u,\,0\right)=\left(-1\right)^{q/2}\left(d^{q}f\left(u,\,\omega\right)/d\omega^{q}\right)|_{\omega=0}$.
Further, $||f^{\left(q\right)}\left(u,\,0\right)||<\infty$ if and
only if $f\left(u,\,\omega\right)$ is $q$ times differentiable at
$\omega=0$. We define
\begin{align}
\mathrm{MSE}\left(Tb_{1,T}b_{2,T},\,\widetilde{J}_{T},\,W\right) & =Tb_{1,T}b_{2,T}\mathbb{E}\left[\mathrm{vec}\left(\widetilde{J}_{T}-J_{T}\right)'W\mathrm{vec}\left(\widetilde{J}_{T}-J_{T}\right)\right].\label{Eq: 3.5 MSE}
\end{align}
\begin{thm}
\label{Theorem MSE J}Suppose $K_{1}\left(\cdot\right)\in\boldsymbol{K}_{1}$,
$K_{2}\left(\cdot\right)\in\boldsymbol{K}_{2}$, Assumption \ref{Assumption Smothness of A (for HAC)}-\ref{Assumption A - Dependence}
hold, $b_{1,T},\,b_{2,T}\rightarrow0$, $n_{T}\rightarrow\infty,\,n_{T}/T\rightarrow0$
and $1/Tb_{1,T}b_{2,T}\rightarrow0$. We have: (i)~
\begin{align*}
\lim_{T\rightarrow\infty} & Tb_{1,T}b_{2,T}\mathrm{Var}\left[\mathrm{vec}\left(\widetilde{J}_{T}\right)\right]\\
& =4\pi^{2}\int K_{1}^{2}\left(y\right)dy\int_{0}^{1}K_{2}^{2}\left(x\right)dx\left(I+C_{pp}\right)\left(\int_{0}^{1}f\left(u,\,0\right)du\right)\otimes\left(\int_{0}^{1}f\left(v,\,0\right)dv\right).
\end{align*}
(ii) If $1/Tb_{1,T}^{q}b_{2,T}\rightarrow0$, $n_{T}/Tb_{1,T}^{q}\rightarrow0$
and $b_{2,T}^{2}/b_{1,T}^{q}\rightarrow0$ for some $q\in[0,\,\infty)$
for which $K_{1,q},$ $||\int_{0}^{1}f^{\left(q\right)}\left(u,\,0\right)du||\in[0,\,\infty)$,
then $\lim_{T\rightarrow\infty}b_{1,T}^{-q}\mathbb{E}(\widetilde{J}_{T}-J_{T})=-2\pi K_{1,q}\int_{0}^{1}f^{\left(q\right)}\left(u,\,0\right)du.$
(iii) If $n_{T}/Tb_{1,T}^{q}\rightarrow0$, $b_{2,T}^{2}/b_{1,T}^{q}\rightarrow0$
and $Tb_{1,T}^{2q+1}b_{2,T}\rightarrow\gamma\in\left(0,\,\infty\right)$
for some $q\in[0,\,\infty)$ for which $K_{1,q},\,||\int_{0}^{1}f^{\left(q\right)}\left(u,\,0\right)du||\in[0,\,\infty)$,
then
\begin{align*}
\lim_{T\rightarrow\infty} & \mathrm{MSE}\left(Tb_{1,T}b_{2,T},\,\widetilde{J}_{T},\,W\right)=4\pi^{2}\left[\gamma K_{1,q}^{2}\mathrm{vec}\left(\int_{0}^{1}f^{\left(q\right)}\left(u,\,0\right)du\right)'W\mathrm{vec}\left(\int_{0}^{1}f^{\left(q\right)}\left(u,\,0\right)du\right)\right.\\
& \quad\left.+\int K_{1}^{2}\left(y\right)dy\int K_{2}^{2}\left(x\right)dx\,\mathrm{tr}\left(W\left(I_{p^{2}}+C_{pp}\right)\left(\int_{0}^{1}f\left(u,\,0\right)du\right)\otimes\left(\int_{0}^{1}f\left(v,\,0\right)dv\right)\right)\right].
\end{align*}
\end{thm}
If $b_{2,T}^{2}/b_{1,T}^{q}\rightarrow\nu<\infty$ replaces $b_{2,T}^{2}/b_{1,T}^{q}\rightarrow0$
in part (ii), then the asymptotic bias for the case of locally stationary
processes becomes
\begin{align}
\lim_{T\rightarrow\infty}b_{1,T}^{-q}\mathbb{E}(\widetilde{J}_{T}-J_{T}) & =-2\pi K_{1,q}\int_{0}^{1}f^{\left(q\right)}\left(u,\,0\right)du+\frac{\nu}{2}\int_{0}^{1}x^{2}K_{2}\left(x\right)\sum_{k=-\infty}^{\infty}\int_{0}^{1}\frac{\partial^{2}}{\partial u^{2}}c\left(u,\,k\right)du.\label{eq Asymptotic MSE full}
\end{align}
For the general case of SLS processes the term involving $\left(\partial^{2}/\partial^{2}u\right)c\left(u,\,k\right)$
is different. The second summand on the right-hand side of \eqref{eq Asymptotic MSE full}
cancels when $\int_{0}^{1}\left(\partial^{2}/\partial^{2}u\right)c\left(u,\,k\right)du=0$.
The latter occurs when the process is stationary. \citet{dahlhaus:12}
presented MSE results for a pointwise estimate of $f\left(u,\,\omega\right)$
under continuity in both components by applying smoothing over $u$
and $\omega.$ His results depends on the local behavior of $f\left(u,\,\omega\right)$
at time $u$ and frequency $\omega$ whereas in our problem the MSE
results depend on properties of the full time path of $f\left(u,\,0\right)$.
The theorem suggests that the optimal choice of $b_{1,T}$ hinges
on the degree of nonstationary in the data, a feature that does not
appear from the corresponding results in the literature. The results
are derived as $n_{T}\rightarrow\infty$. It is possible and indeed
easier to keep $n_{T}$ fixed, in which case the results are unchanged.
However, the case with $n_{T}$ fixed can have some disadvantages
when the spectrum is discontinuous because then the estimator would
be often dealing with observations from different regimes, which as
explained above might lead to low frequency contamination. We now
move to the results concerning $\widehat{J}_{T}$.
\begin{assumption}
\label{Assumption B}(i) $\sqrt{T}(\widehat{\beta}-\beta_{0})=O_{\mathbb{P}}\left(1\right)$;
(ii) $\sup_{u\in\left[0,\,1\right]}\mathbb{E}||V_{\left\lfloor Tu\right\rfloor }||^{2}<\infty$;
(iii) $\sup_{u\in\left[0,\,1\right]}\mathbb{E}\sup_{\beta\in\Theta}$
$||\left(\partial/\partial\beta'\right)V_{\left\lfloor Tu\right\rfloor }\left(\beta\right)||^{2}<\infty$;
(iv) $\int_{-\infty}^{\infty}\left|K_{1}\left(y\right)\right|dy,$
$\int_{0}^{1}\left|K_{2}\left(x\right)\right|dx<\infty.$
\end{assumption}
Assumption \ref{Assumption B}-(i,iii) is the same as Assumption B
in \citet{andrews:91}. As remarked above, we interpret $\beta_{0}$
as the pseudo-true parameter $\beta^{*}$ when the model is misspecified.
Part (iv) of the assumption is satisfied by most commonly used kernels.
In order to obtain rate of convergence results we replace Assumption
\ref{Assumption A - Dependence} with the following assumptions.
\begin{assumption}
\label{Assumption C Andrews 91}(i) Assumption \ref{Assumption A - Dependence}
holds with $V_{t,T}$ replaced by
\begin{align*}
\left(V'_{t},\,\mathrm{vec}\left(\left(\frac{\partial}{\partial\beta'}V_{t}\left(\beta_{0}\right)\right)-\mathbb{E}\left(\frac{\partial}{\partial\beta'}V_{t}\left(\beta_{0}\right)\right)\right)'\right)' & .
\end{align*}
(ii) $\sup_{u\in\left[0,\,1\right]}\mathbb{E}(\sup_{\beta\in\Theta}||\left(\partial^{2}/\partial\beta\partial\beta'\right)V_{\left\lfloor Tu\right\rfloor }^{\left(a\right)}\left(\beta\right)||^{2})<\infty$
for all $a=1,\ldots,\,p$.
\end{assumption}
\begin{assumption}
\label{Assumption W_T and unbounded kernel and Cumulant 8}Let $W_{T}$
denote a $p^{2}\times p^{2}$ weight matrix such that $W_{T}\overset{\mathbb{P}}{\rightarrow}W$.
\end{assumption}
\begin{thm}
\label{Theorem 1 -Consistency and Rate}Suppose $K_{1}\left(\cdot\right)\in\boldsymbol{K}_{1}$,
$K_{2}\left(\cdot\right)\in\boldsymbol{K}_{2}$, $b_{1,T},\,b_{2,T}\rightarrow0$,\textbf{
}$n_{T}\rightarrow\infty,\,n_{T}/Tb_{1,T}\rightarrow0,$ and $1/Tb_{1,T}b_{2,T}\rightarrow0$.
We have:
(i) If Assumption \ref{Assumption Smothness of A (for HAC)}-\ref{Assumption B}
hold, $\sqrt{T}b_{1,T}\rightarrow\infty$, $b_{2,T}/b_{1,T}\rightarrow0$
then $\widehat{J}_{T}-J_{T}\overset{\mathbb{P}}{\rightarrow}0$ and
$\widehat{J}_{T}-\widetilde{J}_{T}\overset{\mathbb{P}}{\rightarrow}0$.
(ii) If Assumption \ref{Assumption Smothness of A (for HAC)}, \ref{Assumption B}-\ref{Assumption C Andrews 91}
hold, $n_{T}/Tb_{1,T}^{q}\rightarrow0$, $1/Tb_{1,T}^{q}b_{2,T}\rightarrow0$,
$b_{2,T}^{2}/b_{1,T}^{q}\rightarrow0$ and $Tb_{1,T}^{2q+1}b_{2,T}\rightarrow\gamma\in\left(0,\,\infty\right)$
for some $q\in[0,\,\infty)$ for which $K_{1,q},\,||\int_{0}^{1}f^{\left(q\right)}\left(u,\,0\right)du||\in[0,\,\infty)$,
then $\sqrt{Tb_{1,T}b_{2,T}}(\widehat{J}_{T}-J_{T})=O_{\mathbb{P}}\left(1\right)$
and $\sqrt{Tb_{1,T}}(\widehat{J}_{T}-\widetilde{J}_{T})=o_{\mathbb{P}}\left(1\right).$
(iii) Under the conditions of part (ii) and Assumption \ref{Assumption W_T and unbounded kernel and Cumulant 8},
\begin{align*}
\lim_{T\rightarrow\infty}\mathrm{MSE}\left(Tb_{1,T}b_{2,T},\,\widehat{J}_{T},\,W_{T}\right)=\lim_{T\rightarrow\infty}\mathrm{MSE}\left(Tb_{1,T}b_{2,T},\,\widetilde{J}_{T},\,W\right) & .
\end{align*}
\end{thm}
The consistency result of $\widehat{J}_{T}$ in part (i) applies to
kernels $K_{1}\left(\cdot\right)$ with unbounded support and to bandwidths
$b_{1,T}$ and $b_{2,T}$ such that $1/b_{1,T}b_{2,T}$ grows at rate
$o(\sqrt{T/b_{2,T}})$. Part (ii) yields the consistency of $\widehat{J}_{T}$
with $b_{1,T}$ only required to be $o\left(Tb_{2,T}\right)$. This
rate is slower than the corresponding rate $o\left(T\right)$ of the
classical kernel HAC estimators as shown by \citet{andrews:91} in
his Theorem 1-(b). However, this property is of little practical import
because optimal growth rates typically are less than $T^{1/2}$; for
the QS kernel HAC estimator the optimal growth rate is $T^{1/5}$
while it is $T^{1/3}$ for the Newey-West HAC estimator. Part (ii)
of the theorem presents the rate of convergence of $\widehat{J}_{T}$
which is $\sqrt{Tb_{2,T}b_{1,T}}$. In Section \ref{Section Optimal-Kernels-and},
we compare the rate of convergence of $\widehat{J}_{T}$ with that
of the classical HAC estimators when the respective optimal bandwidths
are used.
\section{Optimal Kernels, \label{Section Optimal-Kernels-and} Bandwidths
and Choice of $n_{T}$}
In this section, we show the optimality of quadratic-type kernels
under MSE criterion.\footnote{Besides \citet{andrews:91} and \citet{newey/west:87} in the context
of LRV estimation, the MSE-optimality criterion was also used more
recently by \citet{whilelm:2015} in a GMM context to determine the
optimal bandwidth of the nonparametric estimator of the optimal weighting
matrix. } For $K_{1},$ the result states that the QS kernel minimizes the
asymptotic MSE for any $K_{2}\left(\cdot\right)$. Let
\begin{align*}
\mathrm{MSE} & (b_{2,T}^{-4},\,\widehat{c}_{T}\left(u_{0},\,k,\,\right),\,\widetilde{W}_{T})\\
& \triangleq b_{2,T}^{-4}\mathbb{E}\left[\mathrm{vec}\left(\widehat{c}_{T}\left(u_{0},\,k\right)-c\left(u_{0},\,k\right)\right)\right]'\widetilde{W}_{T}\left[\mathrm{vec}\left(\widehat{c}_{T}\left(u_{0},\,k\right)-c\left(u_{0},\,k\right)\right)\right],
\end{align*}
where $\widetilde{W}_{T}$ is some $p\times p$ positive semidefinite
matrix. The optimal bandwidths $b_{1,T}^{\mathrm{opt}}$ and $b_{2,T}^{\mathrm{opt}}$
satisfy the following sequential MSE criterion:
\begin{align}
\mathrm{MSE} & \left(Tb_{1,T}^{\mathrm{opt}}\overline{b}_{2,T}^{\mathrm{opt}},\,\widehat{J}_{T}\left(b_{1,T}^{\mathrm{opt}},\,\overline{b}_{2,T}^{\mathrm{opt}}\right),\,W_{T}\right)\leq\mathrm{MSE}\left(Tb_{1,T}^{\mathrm{opt}}\overline{b}_{2,T}^{\mathrm{opt}},\,\widehat{J}_{T}\left(b_{1,T},\,\overline{b}_{2,T}^{\mathrm{opt}}\right),\,W_{T}\right)\label{eq (MSE criterio)}\\
& \quad\mathrm{where}\quad\overline{b}_{2,T}^{\mathrm{opt}}=\int_{0}^{1}b_{2,T}^{\mathrm{opt}}\left(u\right)du\quad\nonumber \\
& \quad\mathrm{and}\quad\quad b_{2,T}^{\mathrm{opt}}\left(u\right)=\underset{b_{2,T}}{\mathrm{argmin}}\,\,\,\mathrm{MSE}\left(b_{2,T}^{-4},\,\widehat{c}_{T}\left(u_{0},\,k\right),\,\widetilde{W}_{T}\right).\nonumber
\end{align}
The first inequality above has to hold as $T\rightarrow\infty$.
The above criterion determines the globally optimal $b_{1,T}^{\mathrm{opt}}$
given the integrated locally optimal $b_{2,T}^{\mathrm{opt}}\left(u\right)$.
Thus, $b_{1,T}^{\mathrm{opt}}$ and $\overline{b}_{2,T}^{\mathrm{opt}}$
need not be the same as the bandwidths $(\widetilde{b}_{1,T}^{\mathrm{opt}},\,\widetilde{b}_{2,T}^{\mathrm{opt}})$
that jointly minimize the global asymptotic MSE,
\begin{equation}
\lim_{T\rightarrow\infty}\mathrm{MSE}\left(Tb_{1,T}b_{2,T},\,\widehat{J}_{T}\left(b_{1,T},\,b_{2,T}\right),\,W_{T}\right).\label{Eq. (Global MSE)}
\end{equation}
Theorem \ref{Theorem MSE J}-(ii) states that, under the condition
$b_{2,T}^{2}/b_{1,T}^{q}\rightarrow0$, the bias only depends on the
smoothing over lagged autocovariances but not on $b_{2,T}$. Then,
the global solution $\widetilde{b}_{2,T}^{\mathrm{opt}}$ would be
trivial: $b_{2,T}$ affects the MSE only through the variance term
and optimality requires to set the bandwidth as large as possible.
In contrast, the MSE criterion \eqref{eq (MSE criterio)} based on
the MSE given in Theorem \ref{Theorem MSE J}-(iii) leads to a unique
solution which can be obtained analytically. Under the condition
$b_{2,T}^{2}/b_{1,T}^{q}\rightarrow\nu<\infty$, Belotti et al. \citeyearpar{belotti/casini/catania/grassi/perron_HAC_Sim_Bandws}\nocite{belotti/casini/catania/grassi/perron:2020}
determined the bandwidths $(\widetilde{b}_{1,T}^{\mathrm{opt}},\,\widetilde{b}_{2,T}^{\mathrm{opt}})$
that jointly minimize \eqref{Eq. (Global MSE)}. They showed that
$\widetilde{b}_{1,T}^{\mathrm{opt}},\,\widetilde{b}_{2,T}^{\mathrm{opt}}=O(T^{-1/6})$
while the optimal bandwidths $(b_{1,T}^{\mathrm{opt}},\,\overline{b}_{2,T}^{\mathrm{opt}})$
from \eqref{eq (MSE criterio)} satisfy $b_{1,T}^{\mathrm{opt}}=O(T^{-4/25})$
and $\overline{b}_{2,T}^{\mathrm{opt}}=O(T^{-1/5})$. Thus, the criterion
\eqref{eq (MSE criterio)} leads to a slightly shorter block length
relative to the global criterion \eqref{Eq. (Global MSE)} (i.e.,
$T\overline{b}_{2,T}^{\mathrm{opt}}<T\widetilde{b}_{2,T}^{\mathrm{opt}}$).
A shorter bock length is beneficial if there is substantial nonstationarity
and implies less sensitivity to low frequency contamination from not
properly accounting for nonstationarity {[}cf. Casini et al. \citeyearpar{casini/perron_Low_Frequency_Contam_Nonstat:2020}{]}.
For a throughout comparison between the two criteria see Belotti
et al. \citeyearpar{belotti/casini/catania/grassi/perron_HAC_Sim_Bandws}.\nocite{belotti/casini/catania/grassi/perron:2020}
\subsection{Optimal $K_{2}\left(\cdot\right)$ and $b_{2,T}$}
Let $F\left(K_{2}\right)\triangleq\int_{0}^{1}K_{2}^{2}\left(x\right)dx$,
$H\left(K_{2}\right)=(\int_{0}^{1}x^{2}K_{2}\left(x\right)dx)^{2}$,
and for any $k\in\mathbb{Z}$,
\begin{align*}
D_{1}\left(u_{0}\right) & \triangleq\mathrm{\,vec}\left(\partial^{2}c\left(u_{0},\,k\right)/\partial u^{2}\right)'\widetilde{W}\mathrm{\,vec}\left(\partial^{2}c\left(u_{0},\,k\right)/\partial u^{2}\right),\quad\\
D_{2}\left(u_{0}\right) & \triangleq\mathrm{tr}\widetilde{W}\left(I_{p^{2}}+C_{pp}\right)\sum_{l=-\infty}^{\infty}c\left(u_{0},\,l\right)\otimes\left[c\left(u_{0},\,l\right)+c\left(u_{0},\,l+2k\right)\right].
\end{align*}
\begin{prop}
\label{Proposition: Optimal Local Covariance}Suppose Assumption \ref{Assumption Smothness of A (for HAC)},
\ref{Assumption B}-\ref{Assumption C Andrews 91} hold and $\widetilde{W}_{T}\overset{\mathbb{P}}{\rightarrow}\widetilde{W}$.
We have for all $a,\,b\leq p,$
\begin{align*}
\mathrm{MSE} & \left(1,\,\widehat{c}_{T}^{\left(a,b\right)}\left(u_{0},\,k\right),\,1\right)\\
& =\frac{1}{4}b_{2,T}^{4}\left(\int_{0}^{1}xK_{2}\left(x\right)dx\right)^{2}\left(\frac{\partial^{2}}{\partial^{2}u}c^{\left(a,b\right)}\left(u_{0},\,k\right)\right)^{2}\\
& \quad+\frac{1}{Tb_{2,T}}\int_{0}^{1}K_{2}^{2}\left(x\right)dx\sum_{l=-\infty}^{\infty}c^{\left(a,b\right)}\left(u_{0},\,l\right)\left[c^{\left(a,b\right)}\left(u_{0},\,l\right)+c^{\left(a,b\right)}\left(u_{0},\,l+2k\right)\right]\\
& \quad+\frac{1}{Tb_{2,T}}\int_{0}^{1}K_{2}^{2}\left(x\right)dx\sum_{h_{1}=-\infty}^{\infty}\kappa_{V,\left\lfloor Tu_{0}\right\rfloor }^{\left(a,b,a,b\right)}\left(-k,\,h_{1},\,h_{1}-k\right)+o\left(b_{2,T}^{4}\right)+O\left(1/\left(b_{2,T}T\right)^{2}\right).
\end{align*}
$\mathrm{MSE}(b_{2,T}^{-4},\,\widehat{c}_{T}\left(u_{0},\,k\right)-c\left(u_{0},\,k\right),\,\widetilde{W}_{T})$
is minimized with
\[
b_{2,T}^{\mathrm{opt}}\left(u_{0}\right)=[H\left(K_{2}^{\mathrm{opt}}\right)D_{1}\left(u_{0}\right)]^{-1/5}\left(F\left(K_{2}^{\mathrm{opt}}\right)\left(D_{2}\left(u_{0}\right)+D_{3}\left(u_{0}\right)\right)\right)^{1/5}T^{-1/5},
\]
where $D_{3}\left(u_{0}\right)$ depends on $\widetilde{\kappa}$
(for $p=1$, $D_{3}\left(u_{0}\right)=\sum_{h_{1}=-\infty}^{\infty}\kappa_{V,\left\lfloor Tu_{0}\right\rfloor }\left(-k,\,h_{1},\,h_{1}-k\right)$),
and $K_{2}^{\mathrm{opt}}\left(x\right)=6x\left(1-x\right),\,0\leq x\leq1$.
In addition if $V_{t}$ is Gaussian, then $D_{3}\left(u_{0}\right)=0$,
for $u_{0}\in\left(0,\,1\right)$.
\end{prop}
The optimal kernel $K_{2}^{\mathrm{opt}}\left(x\right)$ is a transformation
of the Epanechnikov kernel. Optimality of quadratic kernels under
a MSE criterion has been shown in many contexts {[}cf. \citet{epanechnikov:69}
and \citet{priestley:85}{]}. The optimal bandwidth sequence decreases
at rate $T^{-1/5}$ which is the same optimal rate derived in the
context of parameter estimation of locally stationary processes {[}see
e.g., \citet{Dahlhaus/Giraitis:98}{]}. The term $D_{1}\left(u_{0}\right)$
is due to nonstationary, while the term $D_{2}\left(u_{0}\right)$
measures the variability of $\widehat{c}_{T}\left(u_{0},\,k\right)$.
The bandwidth $b_{2,T}^{\mathrm{opt}}$ converges to zero at a slower
rate as the process becomes closer to stationary (i.e., as the square
root of $D_{1}\left(u_{0}\right)$ decreases).
\subsection{Optimal $K_{1}\left(\cdot\right)$ }
We next determine the optimal kernel $K_{1}$ and the optimal bandwidth
sequence $b_{1,T}$ given any $K_{2}$ and any $b_{2,T}$ of order
$O(T^{-1/5})$, i.e., the same order of $b_{2,T}^{\mathrm{opt}}\left(u\right)$
for any $u\in\left[0,\,1\right]$. Let $\widehat{J}_{T}^{\mathrm{QS}}$
denote $\widehat{J}_{T}$ when the latter is based on the QS kernel.
For some results below, we consider a subset of $\boldsymbol{K}_{1}$.
Let $\boldsymbol{\widetilde{K}}_{1}=\{K_{1}\left(\cdot\right)\in\boldsymbol{K}_{1}|\,\widetilde{K}\left(\omega\right)\geq0\,\forall\,\omega\in\mathbb{R}\}\bigr\}$
where $\widetilde{K}\left(\omega\right)=\left(2\pi\right)^{-1}\int_{-\infty}^{\infty}K_{1}\left(x\right)e^{-ix\omega}dx.$
The function $\widetilde{K}\left(\omega\right)$ is referred to as
the spectral window generator. The set $\widetilde{\boldsymbol{K}}_{1}$
contains all kernels $K_{1}$ that necessarily generate positive semidefinite
estimators in finite samples.
We adopt the notation $\widehat{J}_{T}\left(b_{1,T}\right)=\widehat{J}_{T}\left(b_{1,T},\,b_{2,T},\,K_{2}\right)$
to denote the estimator $\widehat{J}_{T}$ that uses $b_{1,T},\,b_{2,T}=\overline{b}_{2,T}^{\mathrm{opt}}+o\left(T^{-1/5}\right)$
and $K_{2}^{\mathrm{}}\left(\cdot\right)$. We then compare two kernels
$K_{1}$ using comparable\textcolor{red}{{} }bandwidths $b_{1,T}$ which
are defined as follows. Given $K_{1}\left(\cdot\right)\in\widetilde{\boldsymbol{K}}_{1}$,
the QS kernel $K_{1}^{\mathrm{QS}}\left(\cdot\right)$, and a bandwidth
sequence $\left\{ b_{1,T}\right\} $ to be used with the QS kernel,
define a comparable bandwidth sequence $\left\{ b_{1,T,K_{1}}\right\} $
for use with $K_{1}\left(\cdot\right)$ such that both kernel/bandwidth
combinations have the same asymptotic variance when scaled by the
same factor $Tb_{1,T}b_{2,T}$. This means that
\begin{align*}
\lim_{T\rightarrow\infty}\mathrm{MSE} & (Tb_{1,T}b_{2,T},\,\widehat{J}_{T}^{\mathrm{QS}}(b_{1,T})-\mathbb{E}(\widetilde{J}_{T}^{\mathrm{QS}}(b_{1,T}))+J_{T},\,W_{T})\\
=\lim_{T\rightarrow\infty}\mathrm{MSE} & (Tb_{1,T}b_{2,T},\,\widehat{J}_{T}^{\mathrm{}}(b_{1,T,K_{1}})-\mathbb{E}(\widetilde{J}_{T}(b_{1,T,K_{1}}))+J_{T},\,W_{T}).
\end{align*}
This definition yields $b_{1,T,K_{1}}=b_{1,T}/(\int K_{1}^{2}\left(x\right)dx)$
and $b_{1,T,\mathrm{QS}}=b_{1,T}$ since $\int\left(K_{1}^{\mathrm{QS}}\right)^{2}\left(x\right)dx=1$.
\begin{thm}
\label{Theorem Optimal Kernels}Suppose Assumption \ref{Assumption Smothness of A (for HAC)},
\ref{Assumption B}-\ref{Assumption W_T and unbounded kernel and Cumulant 8}
hold, $\int_{0}^{1}||f^{\left(2\right)}\left(u,\,0\right)||du<\infty$,
$b_{2,T}\rightarrow0$, $b_{2,T}^{5}T\rightarrow\eta\in\left(0,\,\infty\right)$,
$(\mathrm{vec}(\int_{0}^{1}f^{\left(q\right)}\left(u,\,0\right)du))'W\mathrm{vec}(\int_{0}^{1}f^{\left(q\right)}\left(u,\,0\right)du)>0$
and $W$ is positive semidefinite. For any bandwidth sequence $\left\{ b_{1,T}\right\} $
such that $b_{2,T}/b_{1,T}\rightarrow0$, $n_{T}/Tb_{1,T}^{2}\rightarrow0$
and $Tb_{1,T}^{5}b_{2,T}\rightarrow\gamma\in\left(0,\,\infty\right)$,
and for any kernel $K_{1}\left(\cdot\right)\in\boldsymbol{\widetilde{K}}_{1}$
used to construct $\widehat{J}_{T}^{\mathrm{}}$, the QS kernel is
preferred to $K_{1}\left(\cdot\right)$ in the sense that
\begin{align*}
\lim_{T\rightarrow\infty} & \left(\mathrm{MSE}\left(Tb_{1,T}b_{2,T},\,\widehat{J}_{T}^{\mathrm{}}\left(b_{1,T,K_{1}}\right),\,W_{T}\right)-\mathrm{MSE}\left(Tb_{1,T}b_{2,T},\,\widehat{J}_{T}^{\mathrm{QS}}\left(b_{1,T}\right),\,W_{T}\right)\right)\\
& =4\gamma\pi^{2}\left(\mathrm{vec}\left(\int_{0}^{1}f^{\left(2\right)}\left(u,\,0\right)du\right)\right)'W\mathrm{vec}\left(\int_{0}^{1}f^{\left(2\right)}\left(u,\,0\right)du\right)\int_{0}^{1}\left(K_{2}^{\mathrm{opt}}\left(x\right)\right)^{2}dx\\
& \quad\times\left[K_{1,2}^{2}\left(\int K_{1}^{2}\left(y\right)dy\right)^{4}-\left(K_{1,2}^{\mathrm{QS}}\right)^{2}\right]\geq0.
\end{align*}
The inequality is strict if $K_{1}\left(x\right)\neq K_{1}^{\mathrm{QS}}\left(x\right)$
with positive Lebesgue measure.
\end{thm}
The requirement $\int_{0}^{1}||f^{\left(2\right)}\left(u,\,0\right)||du<\infty$
is not stringent and it reduces to the one used by \citet{andrews:91}
when $\left\{ V_{t,T}\right\} $ is stationary. If $\int_{0}^{1}||f^{\left(q\right)}\left(u,\,0\right)||du<\infty$
only for some $1\leq q<2$, one can show that any kernel with $K_{1,q}=0$
has smaller asymptotic MSE than a kernel with $K_{1,q}>0$. The QS,
Parzen, and Tukey-Hanning kernels have $K_{1,q}=0$ for $1\leq q<2$,
whereas the Bartlett has $K_{1,q}>0$ for $1\leq q<2$. Thus, the
asymptotic superiority of the former kernels over the Bartlett kernel
holds even if $\int_{0}^{1}||f^{\left(q\right)}\left(u,\,0\right)||du<\infty$
only for $1\leq q<2$.
\subsection{Optimal Predetermined Bandwidth Sequence $b_{1,T}$ }
We now present the predetermined bandwidth sequence that minimizes
the asymptotic MSE given $b_{2,T}=O(b_{2,T}^{\mathrm{opt}})$ and
$K_{2}=K_{2}^{\mathrm{opt}}.$ This optimality result applies to
each kernel $K_{1}\left(\cdot\right)\in\boldsymbol{K}_{1}$ for which
$K_{1,q}\in\left(0,\,\infty\right)$ for some $q\in\left(0,\,\infty\right)$.
Thus, most commonly used kernels are allowed with the exception of
the truncated kernel. Let
\begin{align*}
\phi\left(q\right) & =\frac{\mathrm{vec}\left(\int_{0}^{1}f^{\left(q\right)}\left(u,\,0\right)du\right)'W\mathrm{vec}\left(\int_{0}^{1}f^{\left(q\right)}\left(u,\,0\right)du\right)}{\mathrm{tr}W\left(I_{p^{2}}+C_{pp}\right)\left(\int_{0}^{1}f\left(u,\,0\right)du\right)\otimes\left(\int_{0}^{1}f\left(v,\,0\right)dv\right)}.
\end{align*}
The optimal bandwidth is $b_{1,T}^{\mathrm{opt}}=(2qK_{1,q}^{2}\phi\left(q\right)Tb_{2,T}^{\mathrm{\mathrm{opt}}}/(\int K_{1}^{2}\left(y\right)dy\int_{0}^{1}K_{2}^{2}\left(x\right)dx))^{-1/\left(2q+1\right)}$,
where $\phi\left(q\right)$ is a function of the\textcolor{red}{{} }unknown
spectral density $f\left(\cdot,\,\cdot\right)$. Hence, the optimal
bandwidth $b_{1,T}^{\mathrm{opt}}$ is unknown in practice, and we
consider data-dependent estimates of $\phi\left(q\right)$ in Section
\ref{Section Data-Dependent-Bandwidths}.
\begin{condition}
\label{Condition b_T}$b_{1,T},\,b_{2,T}\rightarrow0$ with $b_{2,T}/b_{1,T}\rightarrow0$,
and $Tb_{1,T}^{2q+1}b_{2,T}\rightarrow\gamma\in\left(0,\,\infty\right)$
for some $q\in[0,\,\infty)$ for which $K_{1,q},\,||\int_{0}^{1}f^{\left(q\right)}\left(u,\,0\right)du||\in[0,\,\infty)$,
where $b_{2,T}^{\mathrm{}}=O(T^{-1/5})$.
\end{condition}
\begin{cor}
\label{Corollary 1 -Optimal b1 }Suppose Assumption \ref{Assumption Smothness of A (for HAC)},
\ref{Assumption B}-\ref{Assumption W_T and unbounded kernel and Cumulant 8}
hold, $||\int_{0}^{1}f^{\left(q\right)}\left(u,\,\omega\right)du||<\infty$,
$\phi\left(q\right)\in\left(0,\,\infty\right)$, and $W$ is positive
definite. Consider $K_{1}\left(\cdot\right)\in\boldsymbol{K}_{1}$
for which $K_{1,q}\in\left(0,\,\infty\right)$ for some $q\in\left(0,\,\infty\right)$.
Then, $\{b_{1,T}^{\mathrm{opt}}\}$ is optimal among the sequences
$\{b_{1,T}\}$ that satisfy Condition \ref{Condition b_T} in the
sense that,
\begin{align*}
\lim_{T\rightarrow\infty} & \Biggl(\mathrm{MSE}\left(\left(Tb_{2,T}^{\mathrm{}}\right)^{2q/\left(2q+1\right)},\,\widehat{J}_{T}^{\mathrm{}}\left(b_{1,T},\,b_{2,T}\right),\,W_{T}\right)\\
& -\mathrm{MSE}\left(\left(Tb_{2,T}^{\mathrm{}}\right)^{2q/\left(2q+1\right)},\,\widehat{J}_{T}^{\mathrm{}}\left(b_{1,T}^{\mathrm{opt}},\,b_{2,T}\right),\,W_{T}\right)\Biggr)\geq0.
\end{align*}
The inequality is strict unless $b_{1,T}=b_{1,T}^{\mathrm{opt}}+o((Tb_{2,T})^{-1/\left(2q+1\right)})$.
\end{cor}
In Corollary \ref{Corollary 1 -Optimal b1 }, $q=2$ for the QS kernel
and so $b_{1,T}^{\mathrm{opt}}=0.6584(\phi\left(2\right)Tb_{2,T}^{\mathrm{\mathrm{opt}}})^{-1/5}(\int_{0}^{1}K_{2}^{2}\left(y\right)dy)^{1/5}.$
For $K_{2}\left(y\right)=K_{2}^{\mathrm{opt}}\left(y\right),$ the
latter reduces to,
\begin{align}
b_{1,T}^{\mathrm{opt}} & =0.6828(\phi\left(2\right)Tb_{2,T}^{\mathrm{\mathrm{opt}}})^{-1/5}.\label{eq: b1 (opt)}
\end{align}
The optimal bandwidth is of order $T^{-4/25}.$ Thus, the optimal
bandwidth sequence decreases to zero at a slower rate than the optimal
bandwidth sequence for the QS kernel-based HAC estimator of \citet{andrews:91},
for which the rate is of order $T^{-1/5}$. The slower rate is due
to the fact that our estimator smooths the spectrum over time through
$K_{2}\left(\cdot\right)$ and this restricts the smoothing of $K_{1}\left(\cdot\right)$.
In particular, the optimal choice of $b_{1,T}$ hinges on the degree
of nonstationary through $b_{2,T}^{\mathrm{\mathrm{opt}}}$. The
more nonstationary are the data, the smaller is $b_{2,T}^{\mathrm{\mathrm{opt}}}$
and the large is $b_{1,T}^{\mathrm{\mathrm{opt}}}$ which means that
less weight is given to $\widehat{\Gamma}\left(k\right)$ for $k\neq0$.
In contrast, the optimal choice of $b_{1,T}$ for the methods proposed
in the literature is independent of the degree of nonstationarity.
When $b_{1,T}$ and $b_{2,T}$ are chosen optimally, the convergence
rate from Theorem \ref{Theorem 1 -Consistency and Rate} reduces to
$T^{8/25}$. Thus, the rate is slower than the corresponding one
for the QS kernel HAC estimator considered in \citet{andrews:91}.
However, it is misleading to compare our DK-HAC estimator with the
classical HAC estimators only on the basis of the rate of convergence.
In fact, the DK-HAC estimators account flexibly for nonstationarity
and are robust to low frequency contamination induced by nonstationarity/misspecification
whereas the classical HAC estimators are not in general {[}cf. Casini
et al. \citeyearpar{casini/perron_Low_Frequency_Contam_Nonstat:2020}{]}.
\subsection{Choice of $n_{T}$\label{subsec Choice-of nT}}
Our MSE analysis does not indicate an optimal value for $n_{T}$.
It only suggests growth rate bounds. When $K_{1}^{\mathrm{QS}}$
is used, $n_{T}$ cannot grow faster than $T^{2/3}$. We set $n_{T}=T^{0.66}$
for the $\mathrm{QS}$ kernel. That is, we choose $n_{T}$ to be
the largest possible value allowed by the condition. Our sensitivity
analysis (not reported) suggests that choosing a smaller $n_{T}$
might result in excessive overlapping of regimes when the process
is SLS (i.e., $m_{0}>0$). See Belotti et al. \citeyearpar{belotti/casini/catania/grassi/perron_HAC_Sim_Bandws}
for more details.
\section{Data-Dependent Bandwidths\label{Section Data-Dependent-Bandwidths}}
In this section we consider estimators $\widehat{J}_{T}$ that use
bandwidths $b_{1,T}$ and $b_{2,T}$ whose values are determined via
data-dependent methods. We use the ``plug-in'' method which
is characterized by plugging-in estimates of unknown quantities into
an asymptotic formula for an optimal bandwidth parameter (i.e., the
expressions for $b_{1,T}^{\mathrm{opt}}$ and $b_{2,T}^{\mathrm{opt}}$
from Section \ref{Section Optimal-Kernels-and}). Section \ref{subsec:Implementation}
explains how to construct the automatic bandwidths while Section \ref{subsec:Theoretical-Results}
presents the corresponding theoretical results.
\subsection{\label{subsec:Implementation}Implementation}
Let us begin with $b_{1,T}^{\mathrm{opt}}$ and then move to $b_{2,T}^{\mathrm{opt}}$.
The first step for the construction of data-dependent bandwidth parameters
is to specify $p$ univariate parametric models for the elements
of $V_{t}=(V_{t}^{\left(1\right)},\ldots,\,V_{t}^{\left(p\right)})'$.
The second step involves the estimation of the parameters. In our
context, the logical estimation methods to use are local (weighted)
least-squares (LS) (i.e., LS method applied to rolling windows) and
nonparametric kernel methods. In a third step, we replace the unknown
parameters in $\phi\left(q\right)$ with corresponding estimates.
Such estimate $\widehat{\phi}\left(q\right)$ is then substituted
into the expression for $b_{1,T}^{\mathrm{opt}}$ to yield the data-dependent
bandwidth $\widehat{b}_{1,T}$:
\begin{align}
\widehat{b}_{1,T} & =\left(2qK_{1,q}^{2}\widehat{\phi}\left(q\right)T\widehat{\overline{b}}_{2,T}/\left(\int K_{1}^{2}\left(y\right)dy\int_{0}^{1}K_{2}^{2}\left(x\right)dx\right)\right)^{-1/\left(2q+1\right)},\label{Eq. (6.1) Andrews 91}
\end{align}
where $\widehat{\overline{b}}_{2,T}=\left(n_{T}/T\right)\sum_{r=1}^{\left\lfloor T/n_{T}\right\rfloor -1}\widehat{b}_{2,T}\left(rn_{T}/T\right)$.
$\widehat{\overline{b}}_{2,T}$ is an average of the estimates $\widehat{b}_{2,T}\left(\cdot\right)$.
Since $b_{2,T}$ depends on $u$, it is more efficient to estimate
it for each block as its optimal value can change over time. In practice,
a reasonable candidate for an approximating parametric model is the
class of first order autoregressive {[}AR(1){]} models for $\{V_{t}^{\left(r\right)}\},\,r=1,\ldots,\,p$
(with different parameters for each $r$) or a first order vector
autoregressive {[}VAR(l){]} model for $\{V_{t}\}$. These classes
were also used by \citet{andrews:91}. However, in our context it
is reasonable to allow the parameters of the AR(1) model to be time-varying.
For parsimony, we consider a time-varying AR(1) with no breaks in
$f\left(u,\,\omega\right)$, i.e., $V_{t}^{\left(r\right)}=a_{1}\left(t/T\right)V_{t-1}^{\left(r\right)}+u_{t}^{\left(r\right)},$
where the $u_{t}^{\left(r\right)}$ need not be independent across
$r$.
The use of $p$ univariate parametric models requires a simple form
for the weight matrix $W$. In particular, $W$ has to be a diagonal
matrix which in turn implies that $\phi\left(q\right)$ reduces to
\begin{align*}
\phi\left(q\right) & =2^{-1}\sum_{r=1}^{p}W^{\left(r,r\right)}\left(\int_{0}^{1}f^{\left(q\right)\left(r,r\right)}\left(u,\,0\right)du\right)^{2}/\sum_{r=1}^{p}W^{\left(r,r\right)}\left(\int_{0}^{1}f^{\left(r,r\right)}\left(u,\,0\right)du\right)^{2}.
\end{align*}
The usual choice is $W^{\left(r,r\right)}=1$ for all $r$ except
that which corresponds to an intercept for which it is set to zero.
An estimate of $f^{\left(r,r\right)}\left(u,\,0\right)$ is $\widehat{f}^{\left(r,r\right)}\left(u,\,0\right)=\left(2\pi\right)^{-1}(\widehat{\sigma}^{\left(r\right)}\left(u\right))^{2}(1-\widehat{a}_{1}^{\left(r\right)}\left(u\right))^{-2}$
while $f^{\left(2\right)\left(r,r\right)}\left(u,\,0\right)$ can
be estimated by $\widehat{f}^{\left(2\right)\left(r,r\right)}\left(u,\,0\right)=3\pi^{-1}$
$((\widehat{\sigma}^{\left(r\right)}\left(u\right))^{2}\widehat{a}_{1}^{\left(r\right)}\left(u\right))(1-\widehat{a}_{1}^{\left(r\right)}\left(u\right))^{-4}$
where $\widehat{a}_{1}^{\left(r\right)}\left(u\right)$ and $\widehat{\sigma}^{\left(r\right)}\left(u\right)$
are the LS estimates computed using local data to the left of $u=t/T$:
\begin{align}
\widehat{a}_{1}^{\left(r\right)}\left(u\right) & =\frac{\sum_{j=\left\lfloor Tu\right\rfloor -n_{2,T}+1}^{\left\lfloor Tu\right\rfloor }\widehat{V}_{j}^{\left(r\right)}\widehat{V}_{j-1}^{\left(r\right)}}{\sum_{j=\left\lfloor Tu\right\rfloor -n_{2,T}+1}^{\left\lfloor Tu\right\rfloor }\left(\widehat{V}_{j-1}^{\left(r\right)}\right)^{2}},\qquad\label{eq: alpha and sigma}\\
\widehat{\sigma}^{\left(r\right)}\left(u\right) & =\left(\sum_{j=\left\lfloor Tu\right\rfloor -n_{2,T}+1}^{\left\lfloor Tu\right\rfloor }\left(\widehat{V}_{j}^{\left(r\right)}-\widehat{a}_{1}^{\left(r\right)}\left(u\right)\widehat{V}_{j-1}^{\left(r\right)}\right)^{2}\right)^{1/2},\nonumber
\end{align}
where $n_{2,T}\rightarrow\infty$. Then, for the QS kernel $K_{1}$,
\begin{align*}
\widehat{\phi}\left(2\right) & =\sum_{r=1}^{p}W^{\left(r,r\right)}\left(18\left(\frac{n_{3,T}}{T}\sum_{j=0}^{\left\lfloor T/n_{3,T}\right\rfloor -1}\frac{\left(\widehat{\sigma}^{\left(r\right)}\left(\left(jn_{3,T}+1\right)/T\right)\widehat{a}_{1}^{\left(r\right)}\left(\left(jn_{3,T}+1\right)/T\right)\right)^{2}}{\left(1-\widehat{a}_{1}^{\left(r\right)}\left(\left(jn_{3,T}+1\right)/T\right)\right)^{4}}\right)^{2}\right)/\\
& \quad\sum_{r=1}^{p}W^{\left(r,r\right)}\left(\frac{n_{3,T}}{T}\sum_{j=0}^{\left\lfloor T/n_{3,T}\right\rfloor -1}\frac{\left(\widehat{\sigma}^{\left(r\right)}\left(\left(jn_{3,T}+1\right)/T\right)\right)^{2}}{\left(1-\widehat{a}_{1}^{\left(r\right)}\left(\left(jn_{3,T}+1\right)/T\right)\right)^{2}}\right)^{2}.
\end{align*}
For most of the results below we can take $n_{3,T}=n_{2,T}=n_{T}.$
After plugging-in $\widehat{\phi}\left(2\right)$ into the formula
\eqref{eq: b1 (opt)}, we have $\widehat{b}_{1,T}=0.6828(\widehat{\phi}\left(2\right)T\widehat{\overline{b}}_{2,T})^{-1/5}.$
We now propose a data-dependent procedure for $b_{2,T}\left(u_{r}\right),$
where $u_{r}=rn_{T}/T$ for $r=1,\ldots,\,\left\lfloor \left(T-n_{T}\right)/n_{T}\right\rfloor $.
We assume that the parameters of the approximating time-varying
AR(1) models change slowly such that the smoothness of $f\left(\cdot,\,\omega\right)$
and thus of $c\left(\cdot,\,k\right)$ is the same as the one that
would arise if $a_{1}\left(u\right)=0.8\left(\cos1.5+\cos4\pi u\right)$
and $\sigma\left(u\right)=\sigma=1$ for all $u\in\left[0,\,1\right]$
{[}cf. \citet{dahlhaus:12}{]}. The reason for imposing this condition
is that it is otherwise difficult to estimate $\left(\partial^{2}/\partial u^{2}\right)c\left(u,\,k\right)$,
which enters $D_{1}\left(u\right)$, from the data. Under the above
specification, the exact expression of $D_{1}\left(u\right)$ can
be computed analytically:
\begin{align*}
D_{1}\left(u\right) & \triangleq(\int_{-\pi}^{\pi}\left[\frac{3}{\pi}\left(1+0.8\left(\cos1.5+\cos4\pi u\right)\exp\left(-i\omega\right)\right)^{-4}\left(0.8\left(-4\pi\sin\left(4\pi u\right)\right)\right)\exp\left(-i\omega\right)\right.\\
& \quad\left.-\frac{1}{\pi}\left|1+0.8\left(\cos1.5+\cos4\pi u\right)\exp\left(-i\omega\right)\right|^{-3}\left(0.8\left(-16\pi^{2}\cos\left(4\pi u\right)\right)\right)\exp\left(-i\omega\right)\right]d\omega)^{2}.
\end{align*}
An estimate of $D_{1}\left(u\right)$ is given by
\begin{align*}
\widehat{D}_{1}\left(u\right) & \triangleq(\left[S_{\omega}\right]^{-1}\sum_{s\in S_{\omega}}\left[\frac{3}{\pi}\left(1+0.8\left(\cos1.5+\cos4\pi u\right)\exp\left(-i\omega_{s}\right)\right)^{-4}\left(0.8\left(-4\pi\sin\left(4\pi u\right)\right)\right)\exp\left(-i\omega_{s}\right)\right.\\
& \quad\left.-\frac{1}{\pi}\left|1+0.8\left(\cos1.5+\cos4\pi u\right)\exp\left(-i\omega_{s}\right)\right|^{-3}\left(0.8\left(-16\pi^{2}\cos\left(4\pi u\right)\right)\right)\exp\left(-i\omega_{s}\right)\right])^{2},
\end{align*}
where $\left[S_{\omega}\right]$ is the cardinality of $S_{\omega}$
and $\omega_{s+1}>\omega_{s}$ with $\omega_{1}=-\pi,\,\omega_{\left[S_{\omega}\right]}=\pi.$
In our simulations we use $S_{\omega}=\left\{ -\pi,\,-3,\,-2,\,-1,\,0,\,1,\,2,\,3,\,\pi\right\} $.
Note that we have computed $\widehat{D}_{1}\left(u\right)$ for
$k=0$ because it makes the computation simpler. Further, this is
consistent with our sequential MSE criterion because $k=0$ is the
only lag for which $K_{1}(0)=1$ for all $K_{1}$ so that the choice
of $K_{1}$ does not influence $b_{2}^{\mathrm{opt}}\left(\cdot\right)$.
It remains to derive an estimate of $D_{2}\left(u\right)$ since $F\left(K_{2}\right)$
and $H\left(K_{2}\right)$ can be computed for a given $K_{2}\left(\cdot\right)$.
We assume that the innovations of the approximating time-varying
AR(1) model satisfy $\mathbb{E}(u_{t}^{\left(r\right)})=0$, $\mathbb{E}((u_{t}^{\left(r\right)})^{2})=\sigma^{2}$
and $\mathbb{E}((u_{t}^{\left(r\right)})^{4})=3\sigma^{4}$ so that
$D_{3}\left(u\right)=0$ for all $u\in\left(0,\,1\right)$. That is,
the term involving the cumulant drops from $b_{2}^{\mathrm{opt}}\left(u\right)$.
In practice this is convenient because it is complex to deal with
consistent estimation of cumulant terms. Note also that $D_{3}\left(u\right)=0$
if $u_{t}$ is Gaussian. Since $c\left(u,\,k\right)$ can be consistently
estimated by $\widehat{c}_{T}\left(u,\,k\right)$, an estimate of
$D_{2}\left(u\right)$ is given by
\[
\widehat{D}_{2}\left(u_{0}\right)\triangleq p^{-1}\sum_{r=1}^{p}\sum_{l=-\left\lfloor T^{4/25}\right\rfloor }^{\left\lfloor T^{4/25}\right\rfloor }\widehat{c}_{T}^{\left(r,r\right)}\left(u_{0},\,l\right)\left[2\widehat{c}_{T}^{\left(r,r\right)}\left(u_{0},\,l\right)\right],
\]
where the number of summands grows at the same rate as $(b_{1,T}^{\mathrm{opt}})^{-1}$;
a different choice is allowed as long as it grows at a slower rate
than $T^{2/5}$. Hence, the estimate of the optimal $b_{2,T}$ is
given by
\begin{align*}
\widehat{b}_{2,T}\left(u_{r}\right) & =1.6786\left(\widehat{D}_{1}\left(u_{r}\right)\right){}^{-1/5}(\widehat{D}_{2}\left(u_{r}\right))^{1/5}T^{-1/5},\qquad\mathrm{where}\qquad u_{r}=rn_{T}/T.
\end{align*}
\subsection{\label{subsec:Theoretical-Results}Theoretical Results}
Next, we establish consistency, rate of convergence and asymptotic
MSE results for the estimator $\widehat{J}_{T}(\widehat{b}_{1,T},\,\widehat{\overline{b}}_{2,T})$
that uses the data-dependent bandwidths $\widehat{b}_{1,T}$ and $\widehat{\overline{b}}_{2,T}$.
As in \citet{andrews:91}, we need to restrict the class of admissible
kernels to the following class:
\begin{align}
\boldsymbol{K}_{3} & =\{K_{3}\left(\cdot\right)\in\boldsymbol{K}_{1}:\,\left(i\right)\,\left|K_{1}\left(x\right)\right|\leq C_{1}\left|x\right|^{-b}\,\mathrm{with\,}b>\max\left(1+1/q,\,3\right)\label{Eq. (2.6) K3 Kernel class}\\
& \quad\,\,\,\mathrm{for}\,\left|x\right|\in\left[\overline{x}_{L},\,D_{T}h_{T}\overline{x}_{U}\right],\,b_{1,T}^{2}h_{T}\rightarrow\infty,\,D_{T}>0,\,\overline{x}_{L},\,\overline{x}_{U}\in\mathbb{R},\,1\leq\overline{x}_{L}<\overline{x}_{U},\,\mathrm{and}\,\nonumber \\
& \quad\,\,\,\mathrm{with\,}b>1+1/q\,\mathrm{\,for}\,\left|x\right|\notin\left[\overline{x}_{L},\,D_{T}h_{T}\overline{x}_{U}\right],\,\mathrm{and\,some\,}C_{1}<\infty,\,\mathrm{where}\,q\in\left(0,\,\infty\right)\nonumber \\
& \quad\,\,\,\mathrm{is\,such\,that\,}K_{1,q}\in\left(0,\,\infty\right),\,\left(ii\right)\,\left|K_{1}\left(x\right)-K_{1}\left(y\right)\right|\leq C_{2}\left|x-y\right|\,\forall x,\,y\in\mathbb{R}\,\mathrm{for\,some}\nonumber \\
& \quad\,\,\,\mathrm{costant}\,C_{2}<\infty,\,\mathrm{and}\,(iii)\,q<34/4\}.\nonumber
\end{align}
Let $\widehat{\theta}$ denote the estimator of the parameter of the
approximate (time-varying) parametric model(s) introduced above. For
example, with univariate AR(1) approximating parametric models, $\widehat{\theta}=(\int_{0}^{1}\widehat{a}_{1}^{\left(1\right)}\left(u\right)du,\,\int_{0}^{1}(\widehat{\sigma}^{\left(1\right)}\left(u\right))^{2}du,\ldots,\,\int_{0}^{1}\widehat{a}_{1}^{\left(p\right)}\left(u\right)du,\,\int_{0}^{1}(\widehat{\sigma}^{\left(p\right)}\left(u\right))^{2}du)'$.
Let $\theta^{*}$ denote the probability limit of $\widehat{\theta}$.
$\widehat{\phi}\left(q\right)$ is the value of $\phi\left(q\right)$
with $\widehat{\theta}$ instead of $\theta$. Its probability limit
is denoted by $\phi_{\theta^{*}}$.
\begin{assumption}
\label{Assumption E-F-G}(i) $\widehat{\phi}\left(q\right)=O\mathbb{_{P}}\left(1\right)$
and $1/\widehat{\phi}\left(q\right)=O\mathbb{_{P}}\left(1\right)$;
(ii) $\inf\{T/n_{3,T},\,\sqrt{n_{2,T}}\}(\widehat{\phi}\left(q\right)-\phi_{\theta^{*}})=O_{\mathbb{P}}\left(1\right)$
for some $\phi_{\theta^{*}}\in\left(0,\,\infty\right)$ where $n_{2,T}/T+n_{3,T}/T\rightarrow0,$
$n_{2,T}^{10/6}/T\rightarrow[c_{2},\,\infty),$ $n_{3,T}^{10/6}/T\rightarrow[c_{3},\,\infty)$
with $0<c_{2},\,c_{3}<\infty$; (iii) $\sup_{u\in\left[0,\,1\right]}\lambda_{\max}(\Gamma_{u}\left(k\right))\leq C_{3}k^{-l}$
for all $k\geq0$ for some $C_{3}<\infty$ and some $l>\max\left\{ 2,\,1+48q/\left(46+20q\right),\,1+q/\left(3/4+q/2\right)\right\} $,
where $q$ is as in $\boldsymbol{K}_{3}$; (iv) uniformly in $u\in\left[0,\,1\right]$,
$\widehat{D}_{1}\left(u\right),\,\widehat{D}_{2}\left(u\right)=O\mathbb{_{P}}\left(1\right)$
and $1/\widehat{D}_{1}\left(u\right),\,1/\widehat{D}_{2}\left(u\right)=O\mathbb{_{P}}\left(1\right)$;
(v) $|\omega_{s+1}-\omega_{s}|=O\left(T^{-1}\right)$ and $\left[S_{\omega}\right]=O\left(T\right)$;
(vi) $\sqrt{Tb_{2,T}\left(u\right)}(\widehat{D}_{2}\left(u\right)-D_{2}\left(u\right))=O_{\mathbb{P}}\left(1\right)$
for all $u\in\left[0,\,1\right]$; (vii) $\boldsymbol{K}_{2}$ includes
kernels that satisfy $|K_{2}\left(x\right)-K_{2}\left(y\right)|\leq C_{4}\left|x-y\right|$
for all $x,\,y\in\mathbb{R}$ and some constant $C_{4}<\infty$.
\end{assumption}
Parts (i)-(ii) of Assumption \ref{Assumption E-F-G} are the nonparametric
analogue to Assumption E and F, respectively, in \citet{andrews:91}.
Part (iii) is satisfied if $\left\{ V_{t}\right\} $ is strong mixing
with mixing numbers that are less stringent than those sufficient
for the cumulant condition in Assumption \ref{Assumption A - Dependence}-(i).
Part (iv) and (vi) extend (i)-(ii) to $\widehat{D}_{1}$ and $\widehat{D}_{2}$.
Part (v) is needed to apply the convergence of Riemann sums. Part
(vi) follows from the asymptotic results about $\widehat{c}_{T}\left(u,\,k\right)$.
Part (vii) requires $K_{2}$ to satisfy Lipschitz continuity. Note
that $\phi_{\theta^{*}}$ coincides with the optimal value $\phi\left(q\right)$
only when the approximate parametric model indexed by $\theta^{*}$
corresponds to the true data-generating mechanism.
Let $b_{\theta_{1},T}=(2qK_{1,q}^{2}\phi_{\theta^{*}}T\overline{b}_{\theta_{2},T}/\int K_{1}^{2}\left(y\right)dy\int_{0}^{1}K_{2}^{2}\left(x\right)dx)^{-1/\left(2q+1\right)},$
where $\overline{b}_{\theta_{2},T}\triangleq\int_{0}^{1}b_{2,T}^{\mathrm{opt}}\left(u\right)du$.
The asymptotic properties of $\widehat{J}_{T}(\widehat{b}_{1,T},\,\widehat{\overline{b}}_{2,T})$
are shown to be equivalent to those of $\widehat{J}_{T}(b_{\theta_{1},T},\,b_{\theta_{2},T})$.
\begin{thm}
\label{Theorem 3 Andrews 91}Suppose $K_{1}\left(\cdot\right)\in\boldsymbol{K}_{3}$,
$q$ is as in $\boldsymbol{K}_{3}$, $K_{2}\left(\cdot\right)\in\boldsymbol{K}_{2}$,\textbf{
}$n_{T}\rightarrow\infty,$ $n_{T}/Tb_{\theta_{1},T}\rightarrow0,$
and $||\int_{0}^{1}f^{\left(q\right)}\left(u,\,0\right)du||<\infty$.
Then,
(i) If Assumption \ref{Assumption Smothness of A (for HAC)}-\ref{Assumption B}
and \ref{Assumption E-F-G}-(i,iv,vii) hold, $n_{3,T}=n_{2,T}=n_{T},$
and $q>1/2$, then $\widehat{J}_{T}(\widehat{b}_{1,T},\,\widehat{\overline{b}}_{2,T})-J_{T}\overset{\mathbb{P}}{\rightarrow}0$.
(ii) If Assumption \ref{Assumption Smothness of A (for HAC)}, \ref{Assumption B}-\ref{Assumption C Andrews 91}
and \ref{Assumption E-F-G}-(ii,iii,v,vi,vii) hold and $n_{T}/Tb_{\theta_{1},T}^{2}\rightarrow0,$
then $\sqrt{Tb_{\theta_{1},T}b_{\theta_{2},T}}$ $(\widehat{J}_{T}(\widehat{b}_{1,T},\,\widehat{\overline{b}}_{2,T})-J_{T})=O_{\mathbb{P}}\left(1\right)$
and $\sqrt{Tb_{\theta_{1},T}b_{\theta_{2},T}}(\widehat{J}_{T}(\widehat{b}_{1,T},\,\widehat{\overline{b}}_{2,T})-\widehat{J}_{T}(b_{\theta_{1},T},\,b_{\theta_{2},T}))=o_{\mathbb{P}}\left(1\right)$.
(iii) Let $\gamma_{\theta}=2qK_{1,q}^{2}\phi_{\theta}/(\int K_{1}^{2}\left(y\right)dy\int_{0}^{1}K_{2}^{2}\left(x\right)dx)$.
If Assumption \ref{Assumption Smothness of A (for HAC)}, \ref{Assumption B}-\ref{Assumption W_T and unbounded kernel and Cumulant 8}
and \ref{Assumption E-F-G}-(ii,iii,v,vi,vii) hold, then
\begin{align*}
\lim_{T\rightarrow\infty} & \mathrm{MSE}\left(T^{4q/10\left(2q+1\right)},\,\widehat{J}_{T}\left(\widehat{b}_{1,T},\,\widehat{\overline{b}}_{2,T}\right),\,W_{T}\right)\\
& =\lim_{T\rightarrow\infty}\mathrm{MSE}\left(Tb_{\theta_{1},T}b_{\theta_{2},T},\,\widehat{J}_{T}\left(b_{\theta_{1},T},\,b_{\theta_{2},T}\right),\,W_{T}\right)\\
& =4\pi^{2}\left[\gamma_{\theta}K_{1,q}^{2}\mathrm{vec}\left(\int_{0}^{1}f^{\left(q\right)}\left(u,\,0\right)du\right)'W\mathrm{vec}\left(\int_{0}^{1}f^{\left(q\right)}\left(u,\,0\right)du\right)\right]\\
& \quad+\int K_{1}^{2}\left(y\right)dy\int K_{2}^{2}\left(x\right)dx\,\mathrm{tr}\left(W\left(I_{p^{2}}-C_{pp}\right)\left(\int_{0}^{1}f\left(u,\,0\right)du\right)\otimes\left(\int_{0}^{1}f\left(v,\,0\right)dv\right)\right).
\end{align*}
\end{thm}
When the chosen parametric model indexed by $\theta$ is correct,
it follows that $\phi_{\theta^{*}}=\phi\left(q\right)$ and $\widehat{\phi}\left(q\right)\overset{\mathbb{P}}{\rightarrow}\phi\left(q\right)$.
The theorem then implies that $\widehat{J}_{T}(\widehat{b}_{1,T},\,\widehat{\overline{b}}_{2,T})$
exhibits the same optimality properties presented in Theorem \ref{Theorem Optimal Kernels}
and Corollary \ref{Corollary 1 -Optimal b1 }. We omit the details.
\section{\label{Section Monte Carlo}Small-Sample Evaluations}
We conduct a Monte Carlo analysis to evaluate the properties of HAR
inference based on the HAC estimator $\widehat{J}_{T}$. We consider
HAR tests in the linear regression model as well as HAR tests used
in the forecast evaluation literature, namely the \citeauthor{diebold/mariano:95}'s
\citeyearpar{diebold/mariano:95} test\nocite{diebold/mariano:95}
and the forecast breakdown test of \citet{giacomini/rossi:09}. The
linear regression models have an intercept and a stochastic regressor.
We focus on the $t$-statistics $t_{r}=\sqrt{T}(\widehat{\beta}^{\left(r\right)}-\beta_{0}^{\left(r\right)})/\sqrt{\widehat{J}_{X,T}^{\left(r,r\right)}}$
where
\begin{align*}
\widehat{J}_{X,T}= & \left(T^{-1}\sum_{t=1}^{T}x_{t}x'_{t}\right)^{-1}\widehat{J}_{T}\left(T^{-1}\sum_{t=1}^{T}x_{t}x'_{t}\right)^{-1},
\end{align*}
is a consistent estimate of the limit of $\mathrm{Var}(\sqrt{T}(\widehat{\beta}-\beta_{0}))$
and $r=1,\,2$. $t_{1}$ is the $t$-statistic for the parameter associated
with the intercept while $t_{2}$ is associated with the stochastic
regressor $x_{t}$. Results for the $F$-test are qualitatively similar
{[}see \citet{casini_diss}{]}. Six basic regression models are considered.
We run a $t$-test on the intercept in model M1 and M5 whereas a $t$-test
on the coefficient of $x_{t}$ is run in model M2-M4 and M6. The
models are based on,
\begin{align}
y_{t} & =\beta_{0}^{\left(1\right)}+\delta+\beta_{0}^{\left(2\right)}x_{t}+e_{t},\qquad\qquad t=1,\ldots,\,T,\label{eq: Model P1}
\end{align}
for the $t$-test on the intercept (i.e., $t_{1}$) and
\begin{align}
y_{t} & =\beta_{0}^{\left(1\right)}+\left(\beta_{0}^{\left(2\right)}+\delta\right)x_{t}+e_{t},\qquad\qquad t=1,\ldots,\,T,\label{eq Model P1 beta2}
\end{align}
for the $t$-test on $\beta_{0}^{\left(2\right)}$ (i.e., $t_{2}$)
where $\delta=0$ under the null hypotheses. In Model M1 $e_{t}=0.5e_{t-1}+u_{t},\,u_{t}\sim\mathrm{i.i.d.\,}\mathscr{N}\left(0,\,0.5\right),$
$x_{t}\sim\mathrm{i.i.d.\,}\mathscr{N}\left(1,\,1\right)$, $\beta_{0}^{\left(1\right)}=0$
and $\beta_{0}^{\left(2\right)}=1.$\footnote{For the results with AR coefficient 0.9 see Table 1 in \citet{casini/perron_PrewhitedHAC}
and footnote \ref{fn:10} below.} Model M2 involves $e_{t}=0.8e_{t-1}+u_{t},\,u_{t}\sim\mathrm{i.i.d.\,}\mathscr{N}\left(0,\,1\right),$
$x_{t}\sim\mathrm{i.i.d.\,}\mathscr{N}\left(1,\,1\right)$, and $\beta_{0}^{\left(1\right)}=\beta_{0}^{\left(2\right)}=0.$
In Model M3 we have segmented locally stationary errors $e_{t}=\rho_{t}e_{t-1}+u_{t},\,u_{t}\sim\mathrm{i.i.d.\,}\mathscr{N}\left(0,\,1\right),\,\rho_{t}=\max\left\{ 0,\,-1\left(\cos\left(1.5-\cos\left(5t/T\right)\right)\right)\right\} $
for $t<4T/5$ and $e_{t}=0.9e_{t-1}+u_{t},\,u_{t}\sim\mathrm{\mathrm{i.i.d.\,}}\mathscr{N}\left(0,\,1\right)$
for $t\geq4T/5,$ and $x_{t}=0.4x_{t-1}+u_{X,t},\,u_{X,t}\sim\mathrm{i.i.d.\,}\mathscr{N}\left(0,\,1\right)$.
Note that $\rho_{t}$ varies smoothly between 0 and 0.8071. Model
M4 involves some misspecification that induces nonstationarity in
the errors,
\begin{align*}
y_{t} & =\beta_{0}^{\left(1\right)}+\left(\beta_{0}^{\left(2\right)}+\delta\right)x_{t}+w_{t}\mathbf{1}\left\{ t\geq4T/5\right\} +e_{t},\qquad\qquad t=1,\ldots,\,T,
\end{align*}
where $e_{t}=\rho_{t}e_{t-1}+u_{t},\,u_{t}\sim\mathrm{i.i.d.\,}\mathscr{N}\left(0,\,1\right),\,\rho_{t}$
as in M3, $x_{t}\sim\mathrm{i.i.d.\,}\mathscr{N}\left(1,\,1\right)$,
and $w_{t}\sim\mathrm{i.i.d.\,}\mathscr{N}\left(2,\,1\right)$ independent
from $x_{t}$. Model M5 involves misspecification under $H_{1}$
via a smooth change in the coefficient $\beta_{0}^{\left(2\right)}$
toward the end of the sample. This situation is very common in practice
and it is motivated by the model for the variable ``cay'' from \citet{bianchi/lettau/ludvigson:18}
(cf. Figure 3 in their paper). The model is given by
\begin{align*}
y_{t} & =\beta_{0}^{\left(1\right)}+\delta+\left(\beta_{0}^{\left(2\right)}+d_{t}\mathbf{1}\left\{ t\geq4.5T/5\right\} \right)x_{t}+e_{t},\qquad\qquad t=1,\ldots,\,T,
\end{align*}
where $d_{t}=1.5\delta\left(t-4.5T/5\right)/T$, $e_{t}=\rho_{t}e_{t-1}+u_{t},\,u_{t}\sim\mathrm{i.i.d.\,}\mathscr{N}\left(0,\,1\right),\,\rho_{t}=0.8(\cos(1.5-\cos(t/$
$2T)))$ for $t\in\left\{ 1,\ldots,\,T/2-1\right\} \cup\left\{ T/2+T/4+1,\ldots,\,T\right\} $
and $e_{t}=0.2e_{t-1}+2u_{t},\,u_{t}\sim\mathrm{\mathrm{i.i.d.\,}}\mathscr{N}\left(0,\,1\right)$
for $T/2\leq t\leq T/2+T/4$, and $x_{t}=2+0.5x_{t-1}+u_{X,t},\,u_{X,t}\sim\mathrm{i.i.d.}\,\mathscr{N}\left(0,\,1\right)$.
That is, $\rho_{t}$ varies smoothly between 0 and 0.7021. Model M6
is given by \eqref{eq Model P1 beta2} where $e_{t}=\rho_{t}e_{t-1}+u_{t},\,u_{t}\sim\mathrm{i.i.d.\,}\mathscr{N}\left(0,\,1\right),\,\rho_{t}=\max\left\{ 0,\,0.3\left(\cos\left(1.5-\cos\left(t/5T\right)\right)\right)\right\} $
for $t\in\left\{ 1,\ldots,\,T/2-1\right\} \cup\{T/2+4,\ldots,$ $T-16\}$
and $e_{t}=0.99e_{t-1}+2u_{t},\,u_{t}\sim\mathrm{\mathrm{i.i.d.\,}}\mathscr{N}\left(0,\,1\right)$
for $T/2\leq t\leq T/2+3$ and $e_{t}=0.9e_{t-1}+2u_{t},\,u_{t}\sim\mathrm{\mathrm{i.i.d.\,}}\mathscr{N}\left(0,\,1\right)$
for $T-15\leq t\leq T$, and $x_{t}\sim\mathrm{i.i.d.}\,\mathscr{N}\left(1,\,1\right)$.
Note that $\rho_{t}\in\left[0,\,0.2633\right]$.
Next, we move to the forecast evaluation tests. The Diebold-Mariano
test statistic is defined as $t_{\mathrm{DM}}\triangleq\sqrt{T_{n}}\,\overline{d}_{L}/\sqrt{\widehat{J}_{d_{L},T}}$,
where $\overline{d}_{L}$ is the average of the loss differentials
between two competing forecast models, $\widehat{J}_{d_{L},T}$ is
an estimate of the asymptotic variance of the the loss differential
series and $T_{n}$ is the number of observations in the out-of-sample.
Throughout we use the quadratic loss. In model M7, we consider an
out-of-sample forecasting exercise with a fixed scheme where, given
a sample of $T$ observations, $0.5T$ observations are used for the
in-sample and the remaining half is used for prediction. The true
model for the target variable is given by $y_{t}=\beta_{0}^{\left(1\right)}+\beta_{0}^{\left(2\right)}x_{t-1}^{(0)}+e_{t}$
where $x_{t-1}^{(0)}\sim\mathrm{i.i.d.}\,\mathscr{N}\left(1,\,1\right)$,
$e_{t}=0.3e_{t-1}+u_{t}$ with $u_{t}\sim\mathrm{i.i.d.\,}\mathscr{N}\left(0,\,1\right)$
and we set $\beta_{0}^{\left(1\right)}=\beta_{0}^{\left(2\right)}=1.$
The two competing models both involve an intercept but differ on the
predictor used in place of $x_{t}^{(0)}$. The first forecast model
uses $x_{t}^{(1)}$ while the second uses $x_{t}^{(2)}$ where $x_{t}^{(1)}$
and $x_{t}^{(2)}$ are independent $\mathrm{i.i.d.}\,\mathscr{N}\left(1,\,1\right)$
sequences, both independent from $x_{t}^{(0)}$. Each forecast model
generates a sequence of $\tau\left(=1\right)$-step ahead out-of-sample
losses $L_{t}^{(i)}$ $\left(i=1,\,2\right)$ for $t=T/2+1,\ldots,\,T-\tau.$
Then $d_{t}\triangleq L_{t}^{(2)}-L_{t}^{(1)}$ denotes the loss differential
at time $t$. The Diebold-Mariano test rejects the null of equal predictive
ability when (after normalization) $\overline{d}$ is sufficiently
far from zero.
Finally, we consider model M8 which we use to investigate the performance
of a $t$-test for forecast breakdown {[}cf. \citet{giacomini/rossi:09}{]}.
Suppose we want to forecast a variable $y_{t}$ following the equation
$y_{t}=\beta_{0}^{\left(1\right)}+\beta_{0}^{\left(2\right)}x_{t-1}+e_{t}$
where $x_{t}\sim\mathrm{i.i.d.\,}\mathscr{N}\left(1,\,1.5\right)$
and $e_{t}=0.3e_{t-1}+u_{t}$ with $u_{t}\sim\mathrm{i.i.d.\,}\mathscr{N}\left(0,\,1\right)$.
For a given forecast model and forecasting scheme, the test of \citet{giacomini/rossi:09}
(GR) detects a forecast breakdown when the average of the out-of-sample
losses differs significantly from the average of the in-sample losses.
The in-sample is used to obtain estimates of $\beta_{0}^{\left(1\right)}$
and $\beta_{0}^{\left(2\right)}$ which are in turn used to construct
out-of-sample forecasts $\widehat{y}_{t}=\widehat{\beta}_{0}^{\left(1\right)}+\widehat{\beta}_{0}^{\left(2\right)}x_{t-1}$.
We set $\beta_{0}^{\left(1\right)}=\beta_{0}^{\left(2\right)}=1.$
We consider a fixed forecasting scheme and one-step ahead forecasts.
The GR's (2009) test statistic is defined as $t^{\mathrm{GR}}\triangleq\sqrt{T_{n}}\overline{SL}/\sqrt{\widehat{J}_{SL}}$
where $\overline{SL}\triangleq T_{n}^{-1}\sum_{t=T_{m}}^{T-\tau}SL_{t+\tau}$,
$SL_{t+\tau}$ is the surprise loss at time $t+\tau$, i.e., the difference
between the time $t+\tau$ out-of-sample loss and in-sample loss,
$SL_{t+\tau}=L_{t+\tau}-\overline{L}_{t+\tau}$, $T_{n}$ is the sample
size in the out-of-sample, $T_{m}$ is the sample size in the in-sample
and $\widehat{J}_{SL}$ is an HAC estimator. We restrict attention
to $\tau=1.$
Throughout our study we consider the following LRV estimators: $\widehat{J}_{T}$
with automatic bandwidths; $\widehat{J}_{T}$ with automatic bandwidths
and the prewhitening of \citet{casini/perron_PrewhitedHAC}; Andrews's
(1991) HAC estimator with automatic bandwidth; Andrews's (1991) HAC
estimator with automatic bandwidth and the prewhitening procedure
of \citet{andrews/monahan:92}; Newey and West's (1987) HAC estimator
with the automatic bandwidth as proposed in \citet{newey/west:94};
Newey and West's (1987) HAC estimator with the automatic bandwidth
as proposed in \citet{newey/west:94} and the prewhitening procedure;
Newey-West with the fixed-$b$ method of \citet{Kiefer/vogelsang/bunzel:00}.\footnote{\label{fn:10}To save space, we do not report results for the Empirical
Weighted Periodogram (EWP) or Empirical Weighted Cosine (EWC) of Lazarus
et al. \citeyearpar{lazarus/lewis/stock:17} and Lazarus et al. \citeyearpar{lazarus/lewis/stock/watson:18},
respectively. Their performance is similar to the method of \citet{Kiefer/vogelsang/bunzel:00}.
The LRV estimator of \citet{Kiefer/vogelsang/bunzel:00} leads to
HAR tests that have better size control. \citet{casini/perron_PrewhitedHAC}
showed that EWC leads to oversized tests when there is strong dependence
in the data relative to the fixed-$b$ method of \citet{Kiefer/vogelsang/bunzel:00}
and to the prewhitened DK-HAC. The power properties of tests normalized
by the EWP and EWC are similar to those using the method of \citet{Kiefer/vogelsang/bunzel:00}.} \citet{casini/perron_PrewhitedHAC} proposed three forms of prewhitening:
(1) $\widehat{J}_{T,\mathrm{pw},1}$ uses a stationary model to whiten
the data; (2) $\widehat{J}_{T,\mathrm{pw},\mathrm{SLS}}$ uses a nonstationary
model to whiten the data; (3) $\widehat{J}_{T,\mathrm{pw},\mathrm{SLS},\mu}$
is the same as $\widehat{J}_{T,\mathrm{pw},\mathrm{SLS}}$ but it
adds a time-varying intercept in the VAR to whiten the data. For model
M7 we also report results using $\widehat{J}_{T}$ and $\widehat{J}_{T,\mathrm{pw},\mathrm{SLS}}$
with the pre-test for breaks in the spectrum as developed in \citet{casini:change-point-spectra}.
We do not report the results for the pre-test for model M1-M6 and
M8 because they are equivalent to those without the pre-test.
For all versions of $\widehat{J}_{T}$ we use $K_{1}^{\mathrm{opt}}$
and $K_{2}^{\mathrm{opt}}$. We set $n_{T}=T^{0.66}$ as explained
in Section \ref{subsec Choice-of nT} and $n_{2,T}=n_{3,T}=n_{T}.$
We consider the following sample sizes for M1-M6: $T=200,\,400$.
Simulation results for additional data-generating processes involving
ARMA, ARCH and heteroskedastic errors are not discussed here because
the results are qualitatively equivalent {[}see, e.g., \citet{casini_diss}
and \citet{casini/perron_PrewhitedHAC}{]}. The significance level
is $\alpha=0.05$ throughout the study.
\subsection{Empirical Sizes of HAR Inference Tests}
Table \ref{Table S1-S2}-\ref{Table Size Forecasting DM-GR} report
the rejection rates for model M1-M8. We begin with the $t$-test in
the linear regression models. As a general pattern, we confirm previous
evidence that the Newey-West's (1987) and Andrews' (1991) HAC estimators
lead to $t$-tests that are oversized when the data are stationary
{[}cf. model M1-M2{]}. The same problem occurs for the Newey-West
(1987) HAC estimator using the usual ``rule'' to determine the number
of lags (not reported). For extreme temporal dependence, simulations
in \citet{casini/perron_PrewhitedHAC} showed that the size distortions
can be even larger especially for the $t$-test on the intercept.
Prewhitening is often effective in helping the HAC estimators to better
control the size under stationarity. However, the simulation results
in \citet{casini/perron_PrewhitedHAC} and in the literature show
that the prewhitened HAC estimators can lead to oversized tests when
there is high serial dependence. The rejection rates of tests normalized
by the Newey-West estimator with fixed-$b$ are the most accurate
in model M1-M2 for $T=200$. Overall, the results in the literature
along with those in \citet{casini_diss} and \citet{casini/perron_PrewhitedHAC}
showed that under stationarity the original fixed-$b$ method of KVB
is in general the least oversized across different degrees of dependence
among all existing methods. Table \ref{Table S1-S2} shows that for
the $t$-test on the intercept the non-prewhitened DK-HAC leads to
HAR tests that are oversized while they are accurate for the $t$-test
on the coefficient on the stochastic regressor. The table also shows
that the prewhitened DK-HAC estimators are competitive with the KVB's
fixed-$b$ in controlling the size. $\widehat{J}_{T,\mathrm{pw},1}$
is the most accurate among the DK-HAC estimators. Since $\widehat{J}_{T,\mathrm{pw},1}$
uses a stationarity VAR model to whiten the data, it works better
than $\widehat{J}_{T,\mathrm{pw},\mathrm{SLS}}$ and $\widehat{J}_{T,\mathrm{pw},\mathrm{SLS},\mu}$
when stationarity actually holds which is consistent with the results
of Table \ref{Table S1-S2}.
Turning to nonstationarity, Table \ref{Table S3-S4} casts concerns
about the performance of existing methods in this context. For both
model M3 and M4, existing LRV estimators lead to HAR tests that have
either size equal or close to zero. The methods that use long bandwidths
(i.e., many lags) such as KVB's fixed-$b$ suffer most from this problem
relative to the classical HAC estimators. This is consistent with
the argument in \citet{casini/perron_Low_Frequency_Contam_Nonstat:2020}
who showed analytically that nonstationarity induces a positive bias
for each sample autocovariance. That bias is constant across different
lags. Since existing LRV estimators are weighted sum of sample autocovariances,
the larger the bandwidth (i.e., the more lagged autocovariances are
included) the larger the positive bias. Thus, LRV estimators are inflated
and HAR tests have rejection rates lower than the significance level.
This mechanism has consequences for power as well, as we show below
that traditional HAR tests have low power. In model M3-M4 the non-prewhitened
DK-HAC and the prewhitened DK-HAC (except $\widehat{J}_{T,\mathrm{pw},1}$)
perform well. $\widehat{J}_{T,\mathrm{pw},1}$ suffers from the same
problem as the existing estimators because it uses stationarity and
when this is violated its performance is affected. In model M5, the
classical HAC estimators yield HAR tests that are oversized. Also
the non-prewhitened DK-HAC is oversized. In contrast, the KVB's fixed-$b$
and the prewhitened DK-HAC have rejection rates close to the significance
level. In model M6, the KVB's fixed-$b$ HAR tests tend to be undersized
whereas the HAC and DK-HAC estimators lead to tests that control the
size more accurately.
Turning to the HAR tests for forecast evaluations, Table \ref{Table Size Forecasting DM-GR}
shows that for model M7 the KVB's fixed-$b$ HAR test has size essentially
equal to zero while the classical HAC estimators yield HAR tests that
are somewhat oversized. In contrast, the tests normalized by the prewhitened
DK-HAC estimators have most accurate rejection rates. In model M8,
the KVB's fixed-$b$ HAR tests are well-sized whereas the classical
HAC estimators lead to tests that are severely undersized. The DK-HAC
estimators control the size reasonably well.
In summary, the prewhitened DK-HAC estimators yield $t$-tests in
regression models with rejection rates that are relatively close to
the nominal size. The non-prewhitened DK-HAC can lead to oversized
tests for the $t$-tests on the intercept if there is high dependence.
Our results confirm the oversize problem of the HAR tests normalized
by the classical HAC estimators documented in the literature under
stationarity. The Fixed-$b$ HAR tests control the size well when
the data are stationary but can show severe undersized issues under
nonstationarity, a problem also affecting the tests normalized by
the classical HAC estimators. Thus, with regards to size control,
the prewhitened DK-HAC estimators are competitive with fixed-$b$
methods under stationarity and they also perform well when the data
are nonstationary.
\subsection{Empirical Power of HAR Inference Tests}
For model M1-M6 we report the values of the power in Table \ref{Table M1 Power}-\ref{Table Power M6}.
The sample size is $T=200$. Power functions for the Diebold-Maraino
and for the forecast breakdown test are presented next. For model
M1, the non-prewhitened HAC and DK-HAC lead to tests that have the
highest power but they were more oversized than the other methods.
The KVB's fixed-$b$ LRV leads to $t$-tests that sacrify some power
relative to the prewhitened HAC and DK-HAC estimators. In model M2,
a similar picture arises. HAR tests normalized by either classical
HAC or DK-HAC estimators have similarly good power while HAR tests
based on KVB's fixed-$b$ have relatively less power. In model M3,
the prewhitening HAC estimators and $\widehat{J}_{T,\mathrm{pw},1}$
(which uses a stationary model to whiten the data) have low power.
The best power is achieved by tests normalized by Andrews' (1991)
HAC estimator and $\widehat{J}_{T}$, followed by $\widehat{J}_{T,\mathrm{pw},\mathrm{SLS}}$
and $\widehat{J}_{T,\mathrm{pw},\mathrm{SLS},\mu}$. The KVB's fixed-$b$
leads to relatively less power than the latter methods. The Newey-West's
(1987) estimator leads to tests that have good power but they were
shown to be oversized. Similar comments apply to model M4. Here Andrews'
(1991) HAC estimator leads to tests that have better power for small
to medium breaks while tests based on $\widehat{J}_{T}$ have better
power for large breaks. In model M5, prewhitening HAC estimators and
KVB's fixed-$b$ lead to HAR tests that have non-monotonic power and
reach zero as $\delta$ increases. This does not occur for the classical
HAC estimators which, however, were oversized. HAR tests based on
$\widehat{J}_{T,\mathrm{QS}}$, $\widehat{J}_{T,\mathrm{pw},\mathrm{SLS}}$
and $\widehat{J}_{T,\mathrm{pw},\mathrm{SLS},\mu}$ perform best for
this model. $\widehat{J}_{T,\mathrm{pw},1}$ results in HAR tests
that have lower power relative to the tests based on the other DK-HAC
because stationarity is violated. In model M6, all HAR tests enjoy
monotonic power with small differences.
Next, let us move to the evaluation of the power properties of the
$t$-tests used in the forecasting literature. We begin with the Diebold-Mariano
test. For this test, the separation between the null and alternative
hypotheses does not depend on the value of a single parameter. Thus,
the data-generating mechanism is different from the one under the
null. The two competing forecast models are as follows: the first
model uses the actual true data-generating process while the second
model differs in that in place of $x_{t-1}^{(0)}$ it uses $x_{t-1}^{(2)}=x_{t-1}^{(0)}+u_{X_{2},t}$
for $t\leq3T/4$ and $x_{t-1}^{(2)}=\delta+x_{t-1}^{(0)}+u_{X_{2},t}$
for $t>3T/4$ with $u_{X_{2},t}\sim\mathrm{i.i.d.\,}\mathscr{N}\left(0,\,1\right)$.
Evidently, the null hypothesis of equal predictive ability should
be rejected by the Diebold-Mariano test whenever $\delta>0$. Table
\ref{Table Power DM Test} reports the power for several values of
$\delta.$ The HAR tests based on existing estimators have lower power
relative to the $\widehat{J}_{T}$ DK-HAC estimators for small values
of $\delta$. When we raise $\delta$ the tests based on the HAC estimators
of \citet{andrews:91} and Newey and West (1987), and KVB's fixed-$b$
method display non-monotonic power gradually converging to zero. In
contrast, the DK-HAC estimators lead to tests that have monotonic
power that reach and maintain unit power. The only exception is the
test based on $\widehat{J}_{T,\mathrm{pw},1}$ that has lower power
because stationarity is violated. The table also reports $\widehat{J}_{T}$
and $\widehat{J}_{T,\mathrm{pw},\mathrm{SLS}}$ with the pre-test
for breaks in the spectrum {[}cf. \citet{casini:change-point-spectra}{]}
that is used for choosing more efficiently how to split the sample
in blocks to compute $\widehat{\Gamma}\left(k\right)$. The pre-test
yields HAR tests with higher power while having the same size as the
corresponding HAR tests with no pre-test. We have not reported the
results with the pre-test for model M1-M6 and M8 because they are
the same as with no pre-test.
Finally, we move to the $t$-test of \citet{giacomini/rossi:09}.
The data-generating process under $H_{1}:\,\mathbb{E}\left(\overline{SL}\right)\neq0$
is given by $y_{t}=1+x_{t-1}+\delta x_{t-1}\mathbf{1}\left\{ t>T_{1}^{0}\right\} +e_{t}$,
where $x_{t-1}\sim\mathrm{i.i.d.}\,\mathscr{N}\left(1.5,\,1\right)$,
$e_{t}=0.3e_{t-1}+u_{t}$, $u_{t}\sim\mathrm{i.i.d.}\mathscr{N}\left(0,\,1\right)$
and $T_{1}^{0}=T\lambda_{1}^{0}$ with $\lambda_{1}^{0}=0.8$. Under
this specification there is a break in the coefficient associated
to the predictor $x_{t-1}$. Thus, there is a forecast failure and
the test of \citet{giacomini/rossi:09} should reject $H_{0}$. From
Table \ref{Table Power GR Test} it appears that all versions of the
classical HAC estimators of \citet{andrews:91} and \citet{newey/west:87},
and KVB's fixed-$b$ lead to $t$-tests that have, essentially, zero
power for all $\delta$. The only exception is Andrews' (1991) HAC
estimator with prewhitening that shows some power but it is not monotonic.
In contrast, the $t$-test based on the DK-HAC estimators have good
power. The failure of existing LRV estimators cannot be attributed
to the sample size because as we raise the sample size to 400 or 800,
the tests still display no power {[}see \citet{casini_diss}{]}.
The failure of the HAR tests based on the existing LRV estimators
occurring in some of the data-generating mechanisms reported here
can be simply reconciled with the fact that in such models the spectrum
of $V_{t}$ is not constant. In other words, the autocovariance of
$V_{t}$ depends not only on the lag order but also on $t$. Existing
LRV estimators estimate an average of a time-varying spectrum. Because
of this instability in the spectrum, they overestimate the extent
of the dependence or variation in $V_{t}$. This is explained analytically
in Casini et al. \citeyearpar{casini/perron_Low_Frequency_Contam_Nonstat:2020}
who showed in a general setting that nonstationarity/misspecification
alters the low frequency components of a time series making the latter
appear as more persistent. Since traditional LRV estimators are a
weighted sum of a large number of low frequency periodogram ordinates,
these estimates turn to be inflated. Similarly, LRV estimators using
long bandwidths (i.e., fixed-$b$) are weighted sum of a large number
of sample autocovariances. Each sample autocovariance is biased upward
so that the latter estimates are even more inflated than the classical
HAC estimators that use a smaller number of sample autocovariances.
This explains why KVB's fixed-$b$ HAR tests are subject to more power
problems, even though the classical HAC estimators are also largely
affected.
The introduction of the smoothing over time in the DK-HAC estimators
avoids the low frequency contamination because observations belonging
to different regimes are not mixed up when computing sample autocovariances.
This guarantees good power properties also under nonstationarity/misspecification
or under nonstationary alternative hypotheses (e.g., HAR tests for
forecast evaluation discussed above). Casini et al. \citeyearpar{casini/perron_Low_Frequency_Contam_Nonstat:2020}
reconciled this issue with some results in the unit root and long
memory literature. Tests for a unit root are known to struggle to
reject the unit root hypotheses if a process is second-order stationary
(i.e., no unit root) but it is contaminated by breaks in the mean
or trend {[}cf. \citeauthor{perron:89} (\citeyear{perron:89}, \citeyear{perron:90}){]}.
Similarly, a short memory sequence contaminated by structural breaks
can approximate a long memory series in that the autocorrelation function
has the same properties as that of a long memory series {[}cf. \citet{diebold/inoue:01},
\citet{hillebrand:05}, \citet{mccloskey/hill:2017} and \citet{mikosh/starica:04}{]}.
\section{\label{Section Conclusions}Conclusions}
Economic time series are highly nonstationary. Methods constructed
under the assumption of stationarity might then have undesirable properties.
This paper developed a theoretical framework for inference in settings
where the data may be nonstationary. A new class of double kernel
heteroskedasticity and autocorrelation consistent (DK-HAC) estimators
was presented. In addition to the usual smoothing procedure over lagged
autocovariances, the estimator applies smoothing over time. This
is important in order to account flexibly for the variation over time
of the structural properties of the economic time series. Optimality
results under MSE criterion concerning bandwidths and kernels have
been established. A data-dependent method based on the ``plug-in''
approach has been proposed. There are empirical relevant circumstances
where HAR tests, either in linear regression models or other contexts,
standardized by existing LRV estimators perform poorly. These may
result in size distortions as well as significant power losses, even
when the sample size is large. In contrast, when the proposed DK-HAC
estimator is used the same HAR tests do not suffer from those issues.
DK-HAC estimators lead to HAR tests that have competitive size control
relative to fixed-$b$ HAR tests, when the latter work well, and have
good power, irrespective of whether there is weak or strong dependence
in the data.
\newpage{}
\bibliographystyle{elsarticle-harv}
\bibliography{References_JoE}
\addcontentsline{toc}{section}{References}
\newpage{}
\newpage{}
\clearpage
\pagenumbering{arabic}