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.
96,496 characters
Heavy Tails and Predictive Ability Testing
\title{Heavy Tails and Predictive Ability Testing
\thanks{We thank Mikkel Bennedsen, Bent Jesper Christensen, Michael Jansson, Søren Johansen, Thomas Mikosch, Anders Rahbek, Frederik Vilandt Rasmussen, Mikkel Sølvsten, Peter N. Sørensen, and seminar participants at Aarhus University for comments and suggestions. Frederiksen and Pedersen are partly supported by the Independent Research Fund Denmark (DFF Grant 7015-00028). Support from the Danish Finance Institute (DFI) is gratefully acknowledged by Pedersen. Matsui is partly supported by the JSPS Grant-in-Aid for Scientific Research C (19K11868).} \\
}
\author{Jonas F. Frederiksen\footnote{University of Copenhagen; email: \tt{[email removed]}} \and Muneya Matsui\footnote{University of Osaka; email: \tt{[email removed]}} \and Rasmus S. Pedersen\footnote{University of Copenhagen and Danish Finance Institute; email: \tt{[email removed]} }}
\maketitle
\date{}
\begin{abstract}
We study the asymptotic behavior of widely used tests for evaluating and comparing predictive accuracy when forecast errors exhibit heavy tails. In particular, when loss differentials have infinite variance, the Diebold–Mariano test statistic converges to a nonstandard limit involving non-Gaussian stable random variables. As a consequence, conventional critical values can yield severely distorted inference: a nominal 5\% test may reject a true null as often as 70\% of the time.
To establish these results, we develop a new stable limit theorem for strongly mixing, infinite-variance time series processes. Building on this theory, we consider subsampling-based inference that remains valid irrespective of tail-heaviness and requires no estimation of long-run variances or tail indices. An application to risk forecasts for emerging-market exchange rates shows that accounting for heavy tails can substantially alter conclusions about predictive performance relative to standard procedures.
\end{abstract}
\medskip\noindent\textbf{Keywords:} Predictive ability, hypothesis testing, heavy tails, stable limit theory, subsampling.
\newpage
\numberwithin{equation}{section}
\newpage
\section{Introduction}
The Diebold--Mariano (DM) test of \cite{diebold1995paring} and the
superior predictive ability (SPA) tests of \cite{white2000snooping} and
\cite{hansen2005SPA} are workhorse tools for forecast comparison in
economics and finance. Their asymptotic justification rests on the
central limit theorem (CLT) of \cite{ibragimov1962clt} for strongly mixing processes with \textit{finite} variance.
This finite-variance requirement is at odds with a basic stylized fact
about the data. Many of the underlying variables -- asset returns,
exchange rate fluctuations, and measures of financial volatility --
have power-law (heavy) tails with index $\kappa>0$, and the empirically
relevant range routinely includes $\kappa<2$, corresponding to infinite
(or undefined) variance; see, e.g., \cite{gabaix2009power} and
\citet[Section 1.2]{ibragimov2015heavytails}. If extreme tail
realizations are not fully captured by a forecast series (or if the forecast series itself is heavy-tailed) the
corresponding forecast errors and loss differentials inherit these
heavy tails.
The mismatch has substantial consequences for inference. In a
simulation experiment, a standard DM test at the 5\% nominal level
rejects a true null hypothesis as often as 70\% of the time,
irrespective of the sample size, under empirically plausible
heavy-tailed data-generating processes. In our empirical application to
risk forecasts for emerging-market foreign exchange rates -- where the
loss differentials appear to have infinite variance -- the standard DM
test would lead a researcher to conclude that a
GARCH-based risk model statistically dominates each of nine competing
forecasts. A heavy-tail-robust alternative reaches a substantively
different conclusion: simple non-parametric rolling-window risk
forecasts cannot be rejected as inferior to any of the more
sophisticated GARCH-based or score-driven competitors.
A central problem in the forecast evaluation literature is to conduct
inference about the mean of a time series $(X_t)_{t=1,\dots,n}$, where
$X_t$ is a forecast error, a forecast loss, a loss differential, or a
transformation thereof. We focus first on testing
\begin{equation}\label{eq:zero_mean}
\mathsf{H}_0:\mathbb{E}[X_t]=0,
\end{equation}
against $\mathsf{H}_\mathrm{A}:\mathbb{E}[X_t]\neq 0$ with
$\mathbb{E}[|X_t|]<\infty$. Depending on how $X_t$ is constructed,
\eqref{eq:zero_mean} covers hypotheses such as forecast unbiasedness,
forecast optimality, forecast encompassing, and, most prominently, the
equal predictive ability (EPA) hypothesis of \cite{diebold1995paring}\footnote{We focus on the unconditional predictive ability hypothesis; the conditional version of \cite{GiacominiWhite2006} is beyond the scope of our analysis. We also abstract from parameter estimation error in the forecasts; see \cite{west1996asymptotic} for asymptotic theory under estimation uncertainty.}.
Standard practice is to base inference on
\begin{equation}\label{DM_test}
T_n=\dfrac{\overline{X}_n}{\hat{\sigma}_n/\sqrt{n}},
\end{equation}
where $\overline{X}_n$ is the sample average and $\hat{\sigma}_n^2$ is
typically a heteroskedasticity- and autocorrelation consistent (HAC)
estimator of the long-run variance of $\overline{X}_n$, and to compare
$T_n$ against standard normal critical values. The \cite{ibragimov1962clt} CLT
provides the foundational justification for this practice, and
underlies a large body of work on forecast evaluation, including
\cite{white2000snooping}, \cite{hansen2005SPA}, \cite{hansen_lunde2005}, \cite{GiacominiWhite2006}, \cite{hansen_lunde_nason2011MCS}, \cite{Patton_timmermann_2012rationality},
\cite{BarendsePatton2022}, and \cite{odendahl2023forecast}.
This paper makes three main contributions. \textit{First}, we develop a
new stable limit theorem for infinite-variance strongly mixing time series. The theorem parallels the structure of the
\cite{ibragimov1962clt} CLT: it imposes regular variation ($\kappa<2$), a mild
anti-clustering condition on extremes, and weak rates on mixing
coefficients -- conditions that hold for a range of standard linear and
nonlinear processes, including autoregressive, GARCH, and stochastic
volatility models. The result thus provides the foundational analog to
\cite{ibragimov1962clt}, that the existing forecast-evaluation literature
requires, when loss differentials have infinite variance.
\textit{Second}, we use this result to characterize the limiting
distribution of the statistic $T_n$ in \eqref{DM_test} under
$\mathsf{H}_0$ when $(X_t)_{t=1,\dots,n}$ has tail index $\kappa<2$. We
show that $T_n$ converges in distribution to a ratio of two dependent
non-Gaussian stable random variables, whose joint distribution is fully
characterized by $\kappa$ and the extremal cluster process of
$(X_t)_{t=1,\dots,n}$. The resulting ratio may be skewed, and its
quantiles can differ substantially from those of the standard normal
distribution used in practice.
\textit{Third}, we consider an inferential procedure that is valid under
\textit{both} finite- and infinite-variance regimes. We construct
critical values by subsampling \citep{politis2012subsampling} and
establish their asymptotic validity from primitives. The procedure is
fully data-driven: it does not require knowledge of $\kappa$, the
convergence rate of the sample mean, or consistent estimation of
long-run variances. The same approach yields confidence intervals for
$\mathbb{E}[X_t]$, which may itself be of interest in forecast
evaluation.
We extend the analysis in two further directions. We turn first to the
SPA tests of \cite{white2000snooping} and \cite{hansen2005SPA}, which
address the more demanding problem of deciding whether a benchmark
forecast is dominated by any one of $m$ competing forecasts. We provide
a multivariate version of our limit theory and derive the limiting
distribution of the SPA test statistic under heavy tails, which lays
the groundwork for a heavy-tail-robust SPA test based on subsampling.
This complements the conditional SPA test of \cite{Li2021CSPA}, which,
like the original, assumes finite higher-order moments.
We then consider alternatives under which $\mathbb{E}[|X_t|]$ is infinite, so that $\mathbb{E}[X_t]$ does not exist or equals $\pm\infty$. We show that the statistic $T_n$ in \eqref{DM_test} does not diverge under such alternatives, so that any test based on it -- whether using standard
normal or subsampled critical values -- is inconsistent. We propose a modified test statistic that diverges whenever $\mathbb{E}[X_t]$ differs from zero or is undefined, together with an accompanying subsampling-based critical value construction. To our knowledge, such alternatives have not been addressed in the existing literature on inference about the mean of heavy-tailed time series.
Our work is most closely related to the working paper by
\cite{kim2021forecast}, who also consider the DM test under heavy
tails. \cite{kim2021forecast} obtain a non-Gaussian stable limit for
the DM statistic by direct application of \cite{Davis1983stablelimit}
and \cite{davis1995point} and propose subsampling for inference,
invoking the validity result of \cite{kokoszka2004subsampling}. They
also provide an analytical treatment of how the tail index of the loss
differential depends on the tail indices of the target variable, the
forecasts, and the loss function, and they illustrate empirically the
relevance of subsampling for realized volatility forecast comparisons.
Our analysis differs in three respects: (i) building upon
\cite{matsui2025self}, we provide a new stable limit theorem for
strongly mixing processes that explicitly accounts for the extremal dependence of the underlying time series process. Our assumptions are
strictly more general than those used in \cite{kim2021forecast},
\cite{McElroy_Politis_2002_subsampling}, and
\cite{kokoszka2004subsampling}, and we verify them for a range of
standard linear and nonlinear time series processes including
autoregressive, GARCH, and stochastic volatility models; (ii) our
analysis extends beyond the EPA setup to forecast encompassing,
optimality, and SPA testing, as well as alternatives under which the
mean may fail to exist; and (iii) we provide rigorous arguments for the
asymptotic validity of the subsampling-based tests.
The remainder of the paper is organized as follows.
Section \ref{sec:limit_theory} presents the stable limit theory and the
resulting limiting distribution of $T_n$ (Proposition \ref{Master_Lemma}).
Section \ref{sec:subsampling} presents the subsampling algorithm and
establishes its asymptotic validity. Section \ref{sec:extensions}
contains the extensions to alternatives with undefined mean and to the
SPA test. Sections \ref{sec:Monte_Carlo} and \ref{sec:empirical} report
the simulation experiments and the empirical illustration. The proof of the main Theorem \ref{theo:inf_var_jointCLT} in the Appendix; additional proofs and simulation results are in
the Supplemental Appendix.
\subsection*{Notation}
For a random variable $X$, unless stated otherwise, $\mathbb{E}[X] \neq 0$
if either $\mathbb{E}[|X|]=\infty$ or that $\mathbb{E}[|X|]<\infty$ with
$\mathbb{E}[X]$ non-zero. For two real-valued functions $f$ and $g$ we
write $f(x)\sim g(x)$ as $x\to \infty$ if $\lim_{x\to \infty}f(x)/g(x)=1$.
We generically denote real-valued, \textit{slowly varying} functions by
$L$, that is, the functions satisfy that for any constant $c>0$,
$L(cx)\sim L(x)$ as $x\to\infty$. For short, a strictly stationary process
is referred to as a stationary process.
\newpage
\section{Set-up and limit theory}\label{sec:limit_theory}
This section develops the limit theory underlying our inferential procedure. Section \ref{sec:setup} provides some motivating examples and previews the main result. Sections \ref{sec:limit_partial_sums}--\ref{sec:limit_partial_sums_hv} state the partial-sum limit theory under finite variance and heavy tails, respectively. Section \ref{sec:limit_T_n} combines these into the limiting distribution of $T_n$ in \eqref{DM_test}, and Section \ref{sec:AR_details} illustrates the theory in the context of linear autoregressive processes.
\subsection{Set-up and preview of main result}\label{sec:setup}
Given a random sample $(X_t)_{t=1,\dots,n}$ with a time-invariant mean, we consider testing the hypothesis $\mathsf{H}_0$ in \eqref{eq:zero_mean} against $\mathsf{H}_\textrm{A}:\mathbb{E}[X_t]\neq0$ with $\mathbb{E}[|X_t|]<\infty$. Following much of the existing literature on forecast evaluation and comparisons, we throughout take $X_t$ as a \textit{primitive}. Depending on the construction of $X_t$, the hypothesis test resembles those often encountered in the forecasting literature. We begin with three motivating examples.
\begin{example}[Forecast Optimality]
Let $(f_{t})_{t=1,\dots,n}$ denote a forecast series of a time series $(Y_t)_{t=1,\dots,n}$, and let $e_{t}=Y_t-f_{t}$ denote the forecast error. The forecast is said to be \textit{unbiased}, if $\mathbb{E}[e_t]=0$, and it is said to be \textit{efficient} if $\mathbb{E}[e_{t}f_{t}]=0$; see, e.g., \citet[p.~6]{odendahl2023forecast}. Hence, forecast unbiasedness or forecast efficiency corresponds to the hypothesis $\mathsf{H}_0$ in \eqref{eq:zero_mean} when $X_t$ is given by $e_t$ or $e_{t}f_{t}$, respectively. In both cases, $X_t$ may be heavy-tailed if the forecast errors are heavy-tailed. \hfill $\blacklozenge$
\end{example}
\begin{example}[Forecast Encompassing]
Let $(f_{1,t})_{t=1,\dots,n}$ and $(f_{2,t})_{t=1,\dots,n}$ denote two forecast series of a time series $(Y_t)_{t=1,\dots,n}$, and let $e_{i,t}=Y_t-f_{i,t}$ denote the forecast error of series $i=1,2$. Following \cite{harvey1998encompassing}, forecast series 1 is said to \textit{encompass} forecast series 2 if $\mathbb{E}[(e_{1,t}-e_{2,t})\,e_{1,t}]=0$. Hence, forecast encompassing means that $\mathsf{H}_0$ holds with $X_t=(e_{1,t}-e_{2,t})\,e_{1,t}$,
which may be heavy-tailed whenever at least one of the forecast errors are. \hfill $\blacklozenge$
\end{example}
\begin{example}[Equal Predictive Ability (EPA)]\label{ex:EPA}
Consider two forecast series $(f_{1,t})_{t=1,\dots,n}$ and $(f_{2,t})_{t=1,\dots,n}$ of a time series $(Y_t)_{t=1,\dots,n}$. For a given loss function $\mathcal{L}:\mathbb{R}^2\to\mathbb{R}_+$, the two forecast series are said to have equal predictive ability (EPA) if
\begin{equation*}
\mathbb{E}[\mathcal{L}(Y_t,f_{1,t})-\mathcal{L}(Y_t,f_{2,t})]=0;
\end{equation*}
see \cite{diebold1995paring}. Hence, the EPA hypothesis corresponds to the hypothesis $\mathsf{H}_0$ in \eqref{eq:zero_mean} with $X_t=\mathcal{L}(Y_t,f_{1,t})-\mathcal{L}(Y_t,f_{2,t})$.
A widely used loss function is the squared error, $\mathcal{L}(Y_t,f_{1,t})=(Y_t-f_{1,t})^2=e_{1,t}^2$. Here, if the forecast error $e_{1,t}$ has tail index $\kappa$, then -- by construction -- the squared error $e_{1,t}^2$ has tail index $\kappa/2$, so that $X_t$ may have infinite variance even when $e_{1,t}$ or $e_{2,t}$ themselves have finite variance.
\medskip
Heavy-tailed losses may also occur when evaluating quantile forecasts as considered in the empirical application in Section \ref{sec:empirical}. Let $f_{1,t}$ denote a $\tau$-level (conditional) quantile forecast of $Y_t$. The associated \textit{tick loss} \citep{giacomini_komunjer2005quantile} is
\begin{equation*}
\mathcal{L}_{\mathrm{tick},\tau}(Y_t,f_{1,t})=(\tau-\mathbf{1}\{Y_t-f_{1,t}<0\})(Y_t-f_{1,t})=(\tau-\mathbf{1}\{e_{1,t}<0\})e_{1,t}.
\end{equation*}
A common way to construct quantile forecasts for financial returns (or losses) is by means of GARCH-type models. Under mild conditions, if $Y_t$ follows a GARCH process, then $Y_t$ has a tail index $\kappa>0$ and the tick losses inherit this tail index; see Section \ref{sec:empirical} for further details and an empirical illustration with emerging-market foreign exchange rates, for which the range $\kappa\in(1,2)$ is empirically relevant. \hfill $\blacklozenge$
\end{example}
In order to test $\mathsf{H}_0$, we consider the statistic $T_n$ in \eqref{DM_test} with $\hat{\sigma}_n^2=n^{-1}\sum_{t=1}^nX_t^2$. In contrast to \cite{diebold1995paring}, we do not use a HAC estimator. This is deliberate: under heavy tails the sample mean has infinite (or undefined) long-run variance, so consistent estimation of this quantity is neither possible nor needed for our purposes. The test statistic $T_n$ is, hence, given by the studentized sum
\begin{equation}\label{self_sum}
T_n=\dfrac{S_n}{\gamma_{n}},
\end{equation}
with
\begin{equation} \label{eq:def_S_gamma}
S_n := \sum_{t=1}^nX_t\qquad\text{and}\qquad \gamma_{n}:=\sqrt{\sum_{t=1}^nX_t^2}\;.
\end{equation}
The normalization of $S_n$ by $\gamma_n$ is important, since these quantities are of the same stochastic order irrespective of the tail-heaviness of $X_t$; in particular, the statistic does not require knowledge about the rate of convergence of $S_n/n$.\smallskip
To preview the main results of this section, we briefly summarize the limiting behavior of $T_n$ under $\mathsf{H}_0$. Given appropriate regularity conditions,
\begin{itemize}
\item as predominantly exploited in the existing body of literature on forecast evaluation and comparisons, if $X_t$ has \textit{finite} variance, then $T_n$ has a \textit{Gaussian} limit;
\item if $X_t$ is heavy-tailed with tail index $\kappa\in(1,2)$, so that $X_t$ has infinite variance, then $T_n$ converges in distribution to the ratio $\xi_\kappa/\zeta_{\kappa/2}$, where $\xi_\kappa$ is a \textit{non-Gaussian} $\kappa$-stable random variable and $\zeta_{\kappa/2}^2$ is a $(\kappa/2)$-stable random variable.
\end{itemize}
The formal statement is given in Proposition \ref{Master_Lemma} in Section \ref{sec:limit_T_n}, which is a direct consequence of the partial-sum limit theorems in Sections \ref{sec:limit_partial_sums}--\ref{sec:limit_partial_sums_hv}. Complementing this result, in Proposition \ref{prop:AR_HAC} we show -- in the context of heavy-tailed linear autoregressive processes -- that the same non-Gaussian limit governs the HAC-based variant of $T_n$ that is routinely used in applied forecast comparisons, so that standard normal critical values are potentially unreliable also for this statistic.
\begin{remark}
Throughout, we focus on testing $\mathsf{H}_0$ in \eqref{eq:zero_mean}, but our results also enable the construction of confidence intervals for $\mathbb{E}[X_t]$ when $X_t$ may have infinite variance. This is useful, e.g., for assessing the accuracy of a forecast series: setting $X_t = \mathcal{L}(Y_t, f_t)$, the quantity $\mathbb{E}[X_t]$ represents the expected loss, such as the mean squared forecast error. See Remark \ref{rem:confidence_interval} for details.
\end{remark}
\subsection{Limit theory for partial sums under finite variance}\label{sec:limit_partial_sums}
We assume throughout that the data-generating process (DGP) for $X_t$ is \textit{stationary}, and we write $(X_t)_{t\in\mathbb{Z}}$ to signify that the process is indexed over all integers.
We start by considering the finite variance case. Before stating our assumptions, we recall the definition of strong mixing coefficients, characterizing the time-dependence of a real-valued stationary process $(X_t)_{t\in\mathbb{Z}}$. Let $g_1$ and $g_2$ be real-valued, bounded functions acting on $(X_t)_{t\in\mathbb{Z}}$ and define
\begin{equation}\label{mix_coef}
\alpha_k:=\sup_{|g_1|,|g_2|\leqslant 1}\left|\mathrm{cov}\left[g_1(\dots,X_{-1},X_{0}),g_2(X_{k},X_{k+1},\dots)\right]\right|,\qquad k \geqslant 0.
\end{equation}
Following \citet[p.3]{Doukhan1994mixing} the numbers $(\alpha_k)_{k\geqslant0}$ are the strong mixing coefficients of $(X_t)_{t\in\mathbb{Z}}$, and the process is said to be strongly mixing (or $\alpha$-mixing) if $\alpha_k\to 0$ as $k\to \infty$.
In order to derive the asymptotic properties of $T_n$ in \eqref{self_sum} under finite variance of $X_t$, we consider the following assumption.
\begin{assumption}\label{mix_coef_ass}
There exists $\varepsilon>0$ such that $\mathbb{E}[|X_t|^{2+\varepsilon}]<\infty$ and the $\alpha$-mixing coefficients of $(X_t)_{t\in\mathbb{Z}}$ satisfy
\begin{equation}\label{eq:mix_summability}
\sum_{k=1}^\infty\alpha_k^{\varepsilon/(2+\varepsilon)}<\infty.
\end{equation}
\end{assumption}
Assumption \ref{mix_coef_ass} is standard in the literature on forecast evaluations. The summability condition in \eqref{eq:mix_summability} holds for any $\varepsilon>0$, if the mixing coefficients $(\alpha_k)_{k\geqslant0}$ decay at a geometric rate. Therefore, a wide range of well-known time series processes, including causal stationary solutions to linear autoregressive, GARCH-type and stochastic volatility processes, satisfy the summability condition; see Section \ref{sec:DGPs} in the Supplemental Appendix for details. On the other hand, the moment condition $\mathbb{E}[|X_t|^{2+\varepsilon}]<\infty$ typically imposes additional restrictions on the process, particularly that $X_t$ has finite variance. With $\overline{\varepsilon}:=\inf\{\varepsilon \in(0,\infty]:\mathbb{E}[|X_t|^{2+\varepsilon}]=\infty\}$, we note that a smaller $\overline{\varepsilon}$ imposes a faster decay on $\alpha_k$ for \eqref{eq:mix_summability} to hold. Consequently, Assumption \ref{mix_coef_ass} balances the order of finite moments and the rate of decay of the mixing coefficients. The assumption implies that a CLT and the ergodic theorem apply to $(X_t - \mathbb{E}[X_t])$ and $(X_t - \mathbb{E}[X_t])^2$, respectively:
\begin{theorem}[\cite{ibragimov1962clt}]\label{theo:rio}
Let $(X_t)_{t\in\mathbb{Z}}$ be a stationary process satisfying Assumption \ref{mix_coef_ass}. Then
\begin{equation*}
\dfrac{1}{\sqrt{n}}\left(\sum\nolimits_{t=1}^n(X_t-\mathbb{E}[X_t]),\sqrt{\sum\nolimits_{t=1}^n{(X_t-\mathbb{E}[X_t])}^2}\right)\xrightarrow{d}\left(Z,\sqrt{\mathrm{var}(X_t)}\right),\qquad n\rightarrow\infty,
\end{equation*}
where $Z$ is a $\mathrm{N}(0,\sigma^2)$-distributed random variable with
\begin{equation}\label{eq:def_sigma2}
\sigma^2:= \mathrm{var}(X_t)+2\sum_{t=1}^\infty\;\mathrm{cov}(X_0,X_t),\qquad 0 <\sigma^2<\infty.
\end{equation}
\end{theorem}
\smallskip
\subsection{Limit theory for partial sums under heavy tails}\label{sec:limit_partial_sums_hv}
We now turn to the case where $X_t$ has infinite second moment. To formalize what we mean by heavy tails, we recall the notion of regular variation: a real-valued random variable $X$ is said to be regularly varying with tail index $\kappa>0$ if
\begin{equation}\label{RV_marg}
\mathbb{P}(X>x)\sim p_+\dfrac{L(x)}{x^\kappa}\quad\quad\text{and}\quad\quad \mathbb{P}(X<-x)\sim p_-\dfrac{L(x)}{x^\kappa},\quad\quad x\to\infty,
\end{equation}
where $L$ is a slowly varying function and the tail-balance coefficients $p_+,p_-\geqslant0$ satisfy $p_++p_-=1$. Importantly, under \eqref{RV_marg}, $\mathbb{E}[{|X|}^{p}]=\infty$ $(<\infty)$ for every $p>\kappa$ $(p<\kappa)$, and, hence, a smaller $\kappa$ implies that $X$ has heavier tails.
Following \cite{basrak2009regularly}, a stationary real-valued process $(X_t)_{t\in\mathbb{Z}}$ is said to be regularly varying with index $\kappa>0$, if $X_0$ is regularly varying with index $\kappa>0$ in the sense of \eqref{RV_marg} and there exists a real-valued stochastic process $(\mathit{\Theta}_t)_{t\in\mathbb{Z}}$ such that for all integers $h\geqslant0$,
\begin{equation}\label{RV_TS}
\mathbb{P}(|X_0|^{-1}(X_{-h},\dots,X_{h})\in\cdot\;|\;|X_0|>x)\xrightarrow{w}\mathbb{P}((\mathit{\Theta}_{-h},\dots,\mathit{\Theta}_{h})\in\cdot),\quad\quad x\to\infty,
\end{equation}
where $\overset{w}{\to}$ denotes weak convergence of probability measures. The process $(\mathit{\Theta}_t)_{t\in\mathbb{Z}}$ is called the \emph{spectral tail process} of $(X_t)_{t\in\mathbb{Z}}$ and characterizes the extremal dependence structure of $(X_t)_{t\in\mathbb{Z}}$. Informally, conditional on a very large excursion at time $0$, $(\mathit{\Theta}_t)_{t\in\mathbb{Z}}$ records the shape of the surrounding cluster of large values, normalized by $|X_0|$. An iid\ process has $\mathit{\Theta}_t=0$ almost surely for all $t\neq0$, while persistent processes give rise to non-degenerate $\mathit{\Theta}_t$. By construction, $\mathbb{P}(|\mathit{\Theta}_0|=1)=1$ with $\mathbb{P}(\mathit{\Theta}_0=\pm1)=p_{\pm}\geqslant 0$. For more details on regularly varying time series we refer to \citet[Chapter 4]{mikosch2024extreme}.
With this notion in hand, we impose the following assumption.
\begin{assumption}\label{ass_RV}
The process $(X_t)_{t\in\mathbb{Z}}$ is regularly varying with index $\kappa\in(0,1)\cup(1,2)$ in the sense of \eqref{RV_TS}.
\end{assumption}
We note that for $\kappa\in(1,2)$, $X_t$ has infinite variance, whereas the variance is not defined for $\kappa\in(0,1)$.\footnote{As is standard in the literature on stable limit theory for dependent processes, we rule out the case $\kappa=1$. The case $\kappa\in(0,1)$ is covered both for completeness and in view of the undefined-mean alternatives studied in Section \ref{sec:extensions}.} Under Assumption \ref{ass_RV} there exists a deterministic sequence $(a_n)_{n\geqslant1}$ satisfying $a_n\to\infty$ and
\begin{equation} \label{eq:def_a_n}
n\mathbb{P}(|X_t|>a_n)\rightarrow1,\qquad n\rightarrow\infty.
\end{equation}
In particular, by regular variation of $X_t$, it holds that $a_n\sim L(n)n^{1/\kappa}$ as $n\to\infty$ for some slowly varying function $L$.
Given Assumption \ref{ass_RV}, we make the following assumption about the dependence structure of $(X_t)_{t\in\mathbb{Z}}$.
\begin{assumption}\label{ass:AC+MX}
Let $(a_n)_{n\geqslant1}$ satisfy \eqref{eq:def_a_n}. There exists a deterministic sequence $(r_n)_{n\geqslant1}$ with $r_n\to\infty$ and $r_n=o(n)$ as $n\to\infty$ such that
\begin{equation}\label{AC}
\lim_{k\rightarrow\infty}\limsup_{n\rightarrow\infty}n\sum_{t=k}^{r_n}\mathbb{E}[\min\{|a_n^{-1}X_t|,1\}\;\min\{|a_n^{-1}X_0|,1\}]=0.
\end{equation}
Furthermore, the strong mixing coefficients $(\alpha_k)_{k\geqslant0}$ of $(X_t)_{t\in\mathbb{Z}}$ satisfy
\begin{equation}\label{alpha_coefs_req}
\dfrac{n}{r_n}\alpha_{\ell_n}\rightarrow0,\qquad n\to\infty,
\end{equation}
for some deterministic sequence $(\ell_n)_{n\geqslant1}$ with $\ell_n\rightarrow\infty$ and $\ell_n=o(r_n)$ as $n\to\infty$.
\end{assumption}
The condition in \eqref{AC} is a so-called anti-clustering condition; see e.g. \citet[Section 2.3]{bartkiewicz2011stable}. To provide some intuition for this condition, note that \eqref{AC} implies
\begin{equation}\label{eq:AC_implication}
\lim_{k\rightarrow\infty}\limsup_{n\rightarrow\infty}n\sum_{t=k}^{r_n} \mathbb{P}\left( \vert X_t\vert >a_n \ , \ \vert X_0 \vert > a_n \right) = 0.
\end{equation}
Hence, informally, if we observe a very large $|X_0| > a_n$, then $(X_t)_{t\geqslant 1}$ cannot stay above this threshold $a_n$ for too long with non-negligible probability; the extremes may cluster, but long clusters become increasingly rare.
The two parts of Assumption \ref{ass:AC+MX} play complementary roles in the limit theory, and it is useful to separate them conceptually. Regular variation (Assumption \ref{ass_RV}) together with the anti-clustering condition \eqref{AC} delivers a stable limit for the sum of $k_n = n/r_n$ independent block sums of $X_t$ with each block having length $r_n$. The mixing condition \eqref{alpha_coefs_req} then ensures that this sum of independent block sums is a good approximation of the full partial sum $S_n=\sum_{t=1}^nX_t$; we refer to Appendix \ref{sec:appendix_main_theorem} for additional details about this approximation argument.
In addition to playing complementary roles, note that the anti-clustering condition \eqref{AC} and the mixing condition \eqref{alpha_coefs_req} depend on the same deterministic sequence $(r_n)_{n\geqslant1}$. Similar to the role of the constant $\varepsilon>0$ in Assumption \ref{mix_coef_ass}, balancing the tail-heaviness and the rate of decay of the mixing coefficients in the finite variance case, the sequence $(r_n)_{n\geqslant1}$ balances the clustering of extremes against the mixing coefficients: if $(X_t)_{t\in\mathbb{Z}}$ is likely to exhibit long clusters of extreme observations such that $\mathbb{P}\left( \vert X_t\vert >a_n \ , \ \vert X_0 \vert > a_n \right)>0$ for large $a_n$ and $t$, the anti-clustering condition requires that $r_n$ grows at a sufficiently slow rate. This slow rate must be compensated for by a faster decay of $\alpha_k$ in order for the mixing condition in \eqref{alpha_coefs_req} to hold.
Assumption \ref{ass:AC+MX} can be shown to hold under mild conditions for a wide range of well-known time series processes, including linear autoregressive processes (as considered in Section \ref{sec:AR_details}), GARCH processes, and other processes given by affine stochastic recurrence equations, as well as stochastic volatility processes; see Section \ref{sec:DGPs} of the Supplemental Appendix for additional details and technical conditions.
Assumptions \ref{ass_RV}--\ref{ass:AC+MX} ensure that $(X_t)_{t\in\mathbb{Z}}$ obeys the following stable limit theorem.
\begin{theorem}\label{theo:inf_var_jointCLT}
Let $(X_t)_{t\in\mathbb{Z}}$ satisfy Assumptions \ref{ass_RV}--\ref{ass:AC+MX}, and define $\mu:=\mathbb{E}[X_t]\mathbf{1}(\kappa>1)$. With $(a_n)_{n\geqslant1}$ defined by \eqref{eq:def_a_n}, it holds that
\begin{equation}\label{eq:inf_var_jointCLT}
a_n^{-1}\left( \sum\nolimits_{t=1}^{n}(X_t-\mu), \sqrt{\sum\nolimits_{t=1}^n(X_t-\mu)^2} \right) \overset{d}{\to} (\xi_\kappa,\;\zeta_{\kappa/2}),\qquad n\to\infty,
\end{equation}
where $\xi_\kappa$ is a $\kappa$-stable random variable, $\mathbb{P}(\zeta_{\kappa/2}>0)=1$, and $\zeta_{\kappa/2}^2$ is a $\kappa/2$-stable random variable.
\end{theorem}
To the best of our knowledge, the above result is new. Whereas many existing stable limit results rely on point process convergence (e.g., \cite{davis1995point}), our result builds upon recent theoretical results of \cite{matsui2025self} relying on convergence of so-called hybrid characteristic function--Laplace transforms. For a real-valued random variable $V$ and a non-negative random variable $W$, the hybrid characteristic function--Laplace transform is the map
\begin{equation*}
(u,\lambda)\mapsto \mathbb{E}\big[\mathrm{e}^{iuV-\lambda W}\big],\qquad u\in\mathbb{R},\ \lambda\geq0,
\end{equation*}
which fully characterizes the joint distribution of $(V,W)$ and reduces to the characteristic function of $V$ when $\lambda=0$ and to the Laplace transform of $W$ when $u=0$. This approach typically applies under milder conditions than those imposed to establish point process convergence, and we carefully adapt the results of \cite{matsui2025self} to strongly mixing processes and allow for $\mu\neq 0$. For further details about hybrid characteristic function-Laplace transforms, we refer to \citet[Appendix D.2]{mikosch2024extreme}, and we refer to the textbook by \cite{samorodnitsky1994stable} for a treatment of non-Gaussian stable distributions.
\medskip
The distribution of the limiting vector in \eqref{eq:inf_var_jointCLT} is fully characterized by the tail index and the extremal dependence structure of $(X_t)_{t\in\mathbb{Z}}$. To see this, recall that by Assumption \ref{ass_RV}, $(X_t)_{t\in\mathbb{Z}}$ is regularly varying and hence, by the definition in \eqref{RV_TS}, has a spectral tail process $(\mathit{\Theta_t})_{t\in\mathbb{Z}}$. Define the \textit{extremal cluster process} as
\begin{equation}\label{eq:extremal_cluster_process}
Q_t=\mathit{\Theta}_t/(\sum_{j\in\mathbb{Z}}|\mathit{\Theta}_j|^\kappa)^{1/\kappa},\qquad t\in\mathbb{Z},
\end{equation}
which is the spectral tail process re-normalized so that $\sum_{t\in\mathbb{Z}}|Q_t|^\kappa=1$.
Then it holds that $(\xi_\kappa,\zeta_{\kappa/2}^2)$ has hybrid characteristic function--Laplace transform given by
\begin{align}
&\mathbb{E}\big[\mathrm{e}^{iu\xi_\kappa-\lambda\zeta_{\kappa/2}^2}\big] \nonumber \\
&=\exp\Big(\int_0^\infty \mathbb{E}\Big[\mathrm{e}^{iyu\sum_{t\in\mathbb{Z}}Q_t-y^2\lambda\sum_{t\in\mathbb{Z}}|Q_t|^2}-1-iyu\sum_{t\in\mathbb{Z}}Q_t\boldsymbol{1}\big\{\kappa\in(1,2)\big\}\Big]d(-y^{-\kappa})\Big), \label{eq:hybrid_chf_Laplace}
\end{align}
with $u\in\mathbb{R}$ and $\lambda\geq0$. In particular, it holds that $\xi_\kappa$ and $\zeta_{\kappa/2}$ are dependent. Setting $\lambda=0$ yields the characteristic function of the $\kappa$-stable random variable $\xi_\kappa$,
\begin{equation}\label{eq:chr fct of xi}
\mathbb{E}[\mathrm{e}^{iu\xi_\kappa}]=\exp\big(-|u|^\kappa\tilde{\sigma}^\kappa(1-i\;\beta\;\textrm{sign}(u)\tan(\kappa\pi/2))\big),\qquad u\in\mathbb{R},
\end{equation}
with \textit{scale} and \textit{skewness} parameters given, respectively, by
\begin{equation}\label{eq: scale and skew}
\tilde{\sigma}=\left(c_\kappa\;\mathbb{E}\big[\big|\sum\nolimits_{t\in\mathbb{Z}}Q_t\big|^\kappa\big]\right)^{1/\kappa}\qquad\text{and}\qquad \beta=\dfrac{\mathbb{E}\big[\big(\sum\nolimits_{t\in\mathbb{Z}}Q_t\big)_+^\kappa-\big(\sum\nolimits_{t\in\mathbb{Z}}Q_t\big)_-^\kappa\big]}{\mathbb{E}\big[\big|\sum\nolimits_{t\in\mathbb{Z}}Q_t\big|^\kappa\big]},
\end{equation}
where
\begin{equation}\label{eq:def_c_kappa}
c_\kappa:=\frac{\Gamma(2-\kappa)\cos(\kappa\pi/2)}{1-\kappa},
\end{equation}
and $\Gamma(\cdot)$ the gamma function\footnote{A slightly different expression for the characteristic function for $\xi_\kappa$ can be found in \citet[Section 9.2]{mikosch2024extreme}. For the sake of completeness, we provide detailed derivations of \eqref{eq:chr fct of xi} and \eqref{eq: scale and skew} in Section \ref{sec:chf_scale_skew} of the Supplemental Appendix.}. The distribution of $\xi_\kappa$ is asymmetric whenever the imaginary part of its characteristic function is non-zero, that is, if $\beta\neq0$. Consequently, whether or not the distribution is symmetric depends on the spectral tail process $(\mathit{ \Theta}_t)_{t\in\mathbb{Z}}$ via \eqref{eq:extremal_cluster_process}. For certain processes, such as linear autoregressive processes considered in Section \ref{sec:AR_details} below, $\beta$ is determined entirely by the tail-balance coefficients $p_{\pm} = \mathbb{P}(\mathit{\Theta}_0=\pm1)$. Moreover, note that the imaginary part of the characteristic function tends to zero as $\kappa\uparrow2$. Intuitively, this property is due to the fact that a stable distribution with index of stability $\kappa = 2$ and a location parameter of zero is a normal distribution with zero mean and hence symmetric.
\subsection{Limiting behavior of $T_n$}\label{sec:limit_T_n}
Combining Theorems \ref{theo:rio} and \ref{theo:inf_var_jointCLT} yields the limiting behavior of the test statistic $T_n$ in \eqref{self_sum} under both finite and infinite variance:
\begin{proposition}\label{Master_Lemma}
Let $(X_t)_{t\in\mathbb{Z}}$ be a stationary process of real-valued random variables, and let $T_n$ be given by \eqref{self_sum}.
\begin{itemize}
\item [$(1)$] Suppose that $(X_t)_{t\in\mathbb{Z}}$ satisfies Assumption \ref{mix_coef_ass} and that $\mathbb{E}[X_t]=0$. Then
\begin{equation*}
T_n\xrightarrow{d}\mathrm{N}(0,\sigma^2/\mathbb{E}[X_t^2]),\qquad n\to\infty,
\end{equation*}
where $\sigma^2>0$ is given in \eqref{eq:def_sigma2}.
\item [$(2)$] Suppose that $(X_t)_{t\in\mathbb{Z}}$ satisfies Assumptions \ref{ass_RV}--\ref{ass:AC+MX} with index $\kappa\in(1,2)$ and $\mathbb{E}[X_t]=0$. Then
\begin{equation*}
T_n\xrightarrow{d}\xi_\kappa\;/\zeta_{\kappa/2},\qquad n\to\infty,
\end{equation*}
with $\xi_\kappa$ and $\zeta_{\kappa/2}^2$ non-Gaussian stable random variables whose joint distribution is characterized by the hybrid characteristic function--Laplace transform given in \eqref{eq:hybrid_chf_Laplace}.
\item [$(3)$] Suppose that $(X_t)_{t\in\mathbb{Z}}$ satisfies Assumption \ref{mix_coef_ass} or Assumptions \ref{ass_RV}--\ref{ass:AC+MX} with index $\kappa\in(1,2)$. If $\mathbb{E}[X_t]\neq0$, then $|T_n|\xrightarrow{\mathbb{P}}\infty$ as $n\to\infty$.
\end{itemize}
\end{proposition}
Cases (1) and (2) give the limiting distribution of $T_n$ under $\mathsf{H}_0$ for different tail regimes, and the two limits may differ substantially in shape. Because $\zeta_{\kappa/2}$ is strictly positive (as discussed in Section \ref{sec:limit_partial_sums_hv}), the limit in case (2) is asymmetric whenever the skewness parameter $\beta$ in \eqref{eq: scale and skew} is non-zero, whereas the Gaussian limit in case (1) is always symmetric. This is investigated in more detail in Section \ref{sec:AR_details} in the context of linear autoregressive processes.
Despite being non-Gaussian, the ratio $\xi_\kappa/\zeta_{\kappa/2}$ can be shown to have all moments finite for a wide class of processes, see \cite{matsui2025moments}; nevertheless, as documented in the simulation study in Section \ref{sec:Monte_Carlo}, its quantiles can depart markedly from those of a standard normal distribution, so that comparing $T_n$ against Gaussian critical values can produce a substantially distorted test.
\subsection{Example: Linear autoregressive processes}\label{sec:AR_details}
We end this section by considering additional details about the results in Theorem \ref{theo:inf_var_jointCLT} and Proposition \ref{Master_Lemma} in the case where $(X_t)_{t\in\mathbb{Z}}$ follows an autoregressive process of order one. This example also provides some supporting intuition for the results of the simulation experiment in Section \ref{sec:Monte_Carlo}.
Let the real-valued $X_t$ obey the recursion
\begin{equation}\label{eq:AR1_sec2}
X_t=\delta +\varphi X_{t-1}+Z_t,\qquad t\in\mathbb{Z},
\end{equation}
where $\delta\in\mathbb{R},\;\varphi\in[0,1)$ and $(Z_t)_{t\in\mathbb{Z}}$ is an iid sequence of random variables with strictly positive Lebesgue density on $\mathbb{R}$ and $\mathbb{E}[{|Z_t|}^\varepsilon]<\infty$ for some $\varepsilon>0$.\footnote{The non-negativity restriction on the autoregressive coefficient $\varphi$ is made entirely for the ease of exposition.} We also assume that $\mathbb{E}[Z_t]=0$ whenever $\mathbb{E}[|Z_t|]<\infty$. These conditions ensure that there exists a stationary causal solution $(X_t)_{t\in\mathbb{Z}}$ to \eqref{eq:AR1_sec2}, which is strongly mixing with geometric rate; see, e.g., Proposition 2.2.4 in \cite{buraczewski2016stochastic}. Consequently, if $\mathbb{E}[{|Z_t|}^{2+\varepsilon}]<\infty$ for some $\varepsilon>0$, then $(X_t)_{t\in\mathbb{Z}}$ satisfies Assumption \ref{mix_coef_ass}. On the other hand, if $Z_t$ is regularly varying with some tail index $\kappa>0$, then $(X_t)_{t\in\mathbb{Z}}$ is regularly varying with the same index $\kappa$:
\begin{proposition}\label{prop:AR1 - RV}
Assume that $Z_t$ is regularly varying with tail index $\kappa>0$ and tail-balance coefficients $p_\pm\geqslant0$ with $p_++p_-=1$. Then there exists a stationary and regularly varying solution $(X_t)_{t\in\mathbb{Z}}$ to \eqref{eq:AR1_sec2} with the same index $\kappa$ and has spectral tail process given by
\begin{align}\label{AR_spec_process}
\mathit{\Theta}_t=\mathit{\Theta}_0\;\varphi^t\mathbf{1}(t\geqslant -J),\qquad t\in\mathbb{Z},
\end{align}
where $J$ is a non-negative integer-valued random variable satisfying
\begin{equation}
\mathbb{P}(J=j)=\varphi^{\kappa j}(1-\varphi^\kappa),\qquad j=0,1,\dots,
\end{equation}
and $\mathbb{P}(\mathit{\Theta}_0=\pm1)=p_\pm$. If $\kappa\in(0,2)$ then $(X_t)_{t\in\mathbb{Z}}$ satisfies Assumption \ref{ass:AC+MX}.
\end{proposition}
By Proposition \ref{prop:AR1 - RV}, $(X_t)_{t\in\mathbb{Z}}$ satisfies Assumptions \ref{ass_RV}--\ref{ass:AC+MX} whenever $Z_t$ is regularly varying with tail index $\kappa\in(0,1)\cup(1,2)$, and we have the following result.
\begin{proposition}\label{prop:AR_lim_moments}
Suppose that $(X_t)_{t\in\mathbb{Z}}$ is given by \eqref{eq:AR1_sec2} and satisfies the conditions of Proposition \ref{prop:AR1 - RV} with tail index $\kappa\in(0,1)\cup(1,2)$ and tail-balance coefficients $p_\pm\in[0,1]$. Assume also that the constant term $\delta=0$ if $\kappa\in(1,2)$. Then
\begin{equation}\label{eq:weak_AR1}
T_n\xrightarrow{d}\xi_\kappa/\zeta_{\kappa/2},\qquad n\to\infty.
\end{equation}
The random variable $\xi_\kappa/\zeta_{\kappa/2}$ has all moments finite, and the skewness parameter of $\xi_\kappa$ is given by $\beta=p_+-p_-$. Moreover, if $\kappa\in(1,2)$, then
\begin{align}
\mathbb{E}\left[\dfrac{\xi_\kappa}{\zeta_{\kappa/2}}\right]&=(p_+-p_-)\;\dfrac{\sqrt{1+\varphi}}{\sqrt{1-\varphi}}\;\dfrac{\Gamma((1-\kappa)/2)}{\sqrt{\pi}\;\Gamma(1-\kappa/2)}, \label{eq:limit_mean} \\[0.1in]
\mathbb{E}\left[\Big(\dfrac{\xi_\kappa}{\zeta_{\kappa/2}}\Big)^2\right]&=\dfrac{1+\varphi}{1-\varphi}\left(1+(p_+-p_-)^2\;\dfrac{\kappa}{2}\;\dfrac{\Gamma((1-\kappa)/2)}{\Gamma(1-\kappa/2)}\right), \label{eq:limit_second_moment}
\end{align}
and $\xi_\kappa/\zeta_{\kappa/2}$ has a continuous Lebesgue density.
\end{proposition}
Proposition \ref{prop:AR_lim_moments} states that (under suitable conditions) the distribution of $\xi_\kappa/\zeta_{\kappa/2}$ has all moments finite, suggesting that the distribution is thin-tailed. Note that, in contrast to the standard normal distribution, the mean of $\xi_\kappa/\zeta_{\kappa/2}$ may be non-zero, depending on the tail-balance coefficients.
In the simulation experiment in Section \ref{sec:Monte_Carlo}, we consider the properties of the standard \cite{diebold1995paring} test when $X_t$ is heavy tailed. The Diebold--Mariano test is based on the test statistic in \eqref{DM_test} with $\hat{\sigma}_n^2$ some HAC estimator and with critical values given by quantiles from the standard normal distribution. Typically, $\hat{\sigma}_n^2$ is the estimator by \cite{newey_west_1987}, and for the sake of simplicity we here let
\begin{align}\label{eq:T_HAC}
T_n^{\mathrm{HAC}}=\frac{S_n}{\sqrt{n}\hat{\sigma}_n}, \qquad \hat{\sigma}^2_n=n^{-1}\sum_{t=1}^nX_t^2+2\sum_{j=1}^{q}w_j\left(n^{-1}\sum_{t=1}^{n-j}X_tX_{t+j}\right),
\end{align}
for a fixed integer $q>0$ and fixed weights $w_j>0$, $j=1,\dots,q$.
We have the following result.
\begin{proposition}\label{prop:AR_HAC}
Suppose that $(X_t)_{t\in\mathbb{Z}}$ is given by \eqref{eq:AR1_sec2} and satisfies the conditions of Proposition \ref{prop:AR1 - RV} with tail index $\kappa\in(0,1)\cup(1,2)$, and that $\delta=0$ if $\kappa\in(1,2)$. With $T_n^{\mathrm{HAC}}$ given by \eqref{eq:T_HAC}, it holds that
\begin{equation}\label{eq:sim_limit_ar} T_n^{\mathrm{HAC}}\xrightarrow{d}\left(1+2\sum_{j=1}^qw_j\varphi^j\right)^{-1/2}\dfrac{\xi_\kappa}{\zeta_{\kappa/2}},\qquad n\to\infty,
\end{equation}
with $\xi_\kappa/\zeta_{\kappa/2}$ the limiting distribution in \eqref{eq:weak_AR1}.
\end{proposition}
\section{Critical value construction}\label{sec:subsampling}
In this section we consider subsampling-based techniques for testing the hypothesis $\mathsf{H}_0$ in \eqref{eq:zero_mean} based on the test statistic $T_n$. Subsampling-based inference about the mean of a heavy-tailed time series has been considered in, e.g., \cite{McElroy_Politis_2002_subsampling}, \cite{kokoszka2004subsampling} and \cite{bai2016unified}. While these papers consider the construction of confidence intervals for $\mathbb{E}[X_t]$, our emphasis is on hypothesis testing. In particular, we extend the setting of \citet[Chapter 3.5]{politis2012subsampling} to hypothesis testing for heavy-tailed time series with unknown rate of convergence of the sample mean, and we construct critical values based on subsample versions of the original test statistic. Consider the following algorithm.
\begin{algorithm} \label{algo}
Given a sample $(X_t)_{t=1,\dots,n}$ of size $n > 1$, choose some integer $b_n\in (0,n)$ and some nominal level $\eta\in(0,1)$.
\begin{enumerate}
\item Compute the statistic given in \eqref{self_sum},
\begin{equation*}
T_n=\dfrac{\sum_{t=1}^nX_t}{(\sum_{t=1}^nX_t^2)^{1/2}}.
\end{equation*}
\item With $q_n := n-b_n+ 1$, compute the subsample statistics
\begin{equation}\label{block_stat}
T_{i,b_n} := \dfrac{\sum_{t=i}^{i+b_n-1} X_t}{(\sum_{t=i}^{i+b_n-1} X_t^2)^{1/2}},\qquad i=1,\dots,q_n.
\end{equation}
\item Define
\begin{equation}\label{eq:def_L_n_b}
L_{n,b_n}(x):= q_n^{-1}\sum_{i=1}^{q_n}\mathbf{1}\big(T_{i,b_n}\leqslant x\big),
\end{equation}
and compute the $\eta/2$ and $(1-\eta/2)$ empirical quantiles of $(T_{i,b_n})_{i=1,\dots,q_n}$, i.e.,
\begin{equation}\label{eq:algo_C_n_b_1}
C_{n,b_n}(\eta/2)=\inf\{x:L_{n,b_n}(x)\geqslant \eta/2\}
\end{equation}
and
\begin{equation}\label{eq:algo_C_n_b_2}
C_{n,b_n}(1-\eta/2)=\inf\{x:L_{n,b_n}(x)\geqslant1-\eta/2\}.
\end{equation}
\end{enumerate}
Then a two-sided equal tailed test of level $\eta$ rejects $\mathsf{H}_0$ if $T_n \notin [C_{n,b_n}(\eta/2), C_{n,b_n}(1-\eta/2)]$.
\end{algorithm}
The algorithm is tailored to work irrespective of $X_t$ having finite variance (Assumption \ref{mix_coef_ass}) or being regularly varying with $\kappa\in(1,2)$ (Assumptions \ref{ass_RV}--\ref{ass:AC+MX}). In particular, the algorithm does not rely on estimation of the long-run variance $\sigma^2$ in Theorem \ref{theo:rio}, and it does not require any knowledge about (or estimation of) the scaling sequence $(a_n)_{n\geqslant1}$ in \eqref{eq:def_a_n}.
The asymptotic validity of the subsampling-based test in Algorithm \ref{algo} is ensured under the strong mixing conditions in Assumptions \ref{mix_coef_ass} and \ref{ass:AC+MX}. In particular, we have the following result.
\begin{theorem}\label{cor:validity_subsampling}
With $(b_n)_{n\geqslant 1}$ provided by Algorithm \ref{algo}, suppose that
$b_n\rightarrow\infty$ and $b_n=o(n)$ as $n\to\infty$.
Given a level $\eta>0$, let $C_{n,b_n}(\eta/2)$ and $C_{n,b_n}(1-\eta/2)$ be given by \eqref{eq:algo_C_n_b_1} and \eqref{eq:algo_C_n_b_2}, respectively.
With $(\xi_\kappa,\zeta_{\kappa/2})$ the limiting random variables in Theorem \ref{theo:inf_var_jointCLT}, define $C(y)=\inf\{x:\mathbb{P}(\xi_\kappa/\zeta_{\kappa/2}\leqslant x)\geqslant y)\}$ for $y\in(0,1)$.
Under $\mathsf{H}_0$, suppose either that
\begin{enumerate}
\item Assumption \ref{mix_coef_ass} holds, or
\item Assumptions \ref{ass_RV}--\ref{ass:AC+MX} hold with $\kappa \in (1,2)$, and $\mathbb{P}(\xi_\kappa/\zeta_{\kappa/2}\leqslant x)$ is continuous in neighborhoods around $C(\eta/2)$ and $C(1-\eta/2)$.
\end{enumerate}
Then
\begin{equation*}
\mathbb{P}\big(T_n\notin [C_{n,b_n}(\eta/2), C_{n,b_n}(1-\eta/2)]\big)\rightarrow \eta,\qquad n\to\infty.
\end{equation*}
\noindent Under either Assumption \ref{mix_coef_ass} or Assumptions \ref{ass_RV}--\ref{ass:AC+MX} with $\kappa \in (1,2)$, if $\mathbb{E}[X_t]\neq0$ then
\begin{equation*}
\mathbb{P}\big(T_n \notin [C_{n,b_n}(\eta/2), C_{n,b_n}(1-\eta/2)]\big)\rightarrow 1,\qquad n\to\infty.
\end{equation*}
\end{theorem}
\begin{remark}\label{rem:limit_T_n}
Theorem \ref{cor:validity_subsampling} implies that the subsampling-based test in Algorithm \ref{algo} is asymptotically valid for a range of processes. In the infinite variance case ($\kappa\in(1,2)$), Assumptions \ref{ass_RV}--\ref{ass:AC+MX} cover a more general class of processes than those considered in existing work on subsampling-based inference about the mean of heavy-tailed time series. For instance, \cite{McElroy_Politis_2002_subsampling} focus entirely on heavy-tailed linear processes. \cite{kokoszka2004subsampling} establish \eqref{eq:inf_var_jointCLT}, but their result hinges on a technical assumption about point process convergence.\footnote{That is, Condition (15) in \cite{kokoszka2004subsampling}. To show that this condition holds for a given process, one typically has to argue $(i)$ that $(X_t)_{t\in\mathbb{Z}}$ is reguarly varying (as in our Assumption \ref{ass_RV}), $(ii)$ that an anti-clustering condition similar to our Assumption \ref{ass:AC+MX} holds \citep[Condition (2.8) in Theorem 2.7]{davis1995point}, and $(iii)$ that the limiting point process is non-null.} In addition, they assume that $\mathbb{E}[X_t\boldsymbol{1}(|X_t|\leqslant x)X_s\boldsymbol{1}(|X_s|\leqslant x)]=0$ for all $x>0$, $t\neq s$.\footnote{Their condition (14). The assumption is made in order to prove a so-called small-vanishing-values condition for the case $\kappa\in(1,2)$, c.f. condition (3.2) of Theorem 3.1 in \cite{davis1995point}. This condition is typically hard to verify for a given process, see \citet[Section 2.4]{bartkiewicz2011stable} for a discussion.} This condition rules out a general class of processes, including non-trivial linear processes, and appears unnecessarily restrictive.
\end{remark}
\begin{remark}\label{rem:singularities}
For $\kappa\in(1,2)$ and $\mathbb{E}[X_t]=0$, Theorem \ref{cor:validity_subsampling} imposes that $\mathbb{P}(\xi_\kappa/\zeta_{\kappa/2}\leqslant x)$ is continuous at $C(\eta/2)$ and $C(1-\eta/2)$. This kind of assumption is standard in the literature on subsampling for heavy-tailed time series. Whereas $\xi_\kappa$ is $\kappa$-stable and hence has a continuous distribution, the distribution of the ratio $\xi_\kappa/\zeta_{\kappa/2}$ may have singularities; see, e.g., \cite{Logan1973self-norm}. In the context of heavy-tailed linear AR(1) processes, as considered in Section \ref{sec:AR_details}, it holds by Proposition \ref{prop:AR_lim_moments} that $\xi_\kappa/\zeta_{\kappa/2}$ does not have any singularity points whenever $\kappa\in(1,2)$.
\end{remark}
Theorem \ref{cor:validity_subsampling} states that the subsampling-based test in Algorithm \ref{algo} is asymptotically valid for \textit{any} block size $b_n=o(n)$ with $b_n\to\infty$. As pointed out by, e.g., \cite{kokoszka2004subsampling}, the finite sample properties of a subsampling procedure depend on $b_n$ chosen by the user. Several data-driven approaches have been considered. If one is willing to assume a specific process (e.g., AR(1) or GARCH) for $X_t$, one may determine $b_n$ by means of simulations as in \cite[Section 9.4.1]{politis2012subsampling}. The so-called Minimum Volatility Method \citep[Section 9.4.2]{politis2012subsampling} avoids this assumption, but requires a user-chosen grid of candidate values. \cite{kim2021forecast} propose a
rule depending on estimated tail-index and skewness parameters in the
iid stable-distribution setting. We adopt the simple rule
\begin{equation}\label{eq:b_n_choice}
b_n = \lfloor1.5n^{0.5}\rfloor,
\end{equation}
similar to block sizes used in \citet[Section 5]{bai2016unified} and \citet[Section 3]{zhang2013blocksampling}, for simplicity and to
maintain consistency with existing work on subsampling for heavy-tailed
dependent data. The Monte Carlo experiments in Section \ref{sec:Monte_Carlo} indicate that this rule yields reasonable
finite-sample rejection frequencies across a range of tail indices
and sample sizes.
\begin{remark}
Algorithm \ref{algo} provides critical values for an equal-tailed two-sided test. Alternatively, one could consider a symmetric two-sided test where the critical value is given by the $(1-\eta)$ empirical quantile of $|T_{i,b_n}|$, and the test rejects if $|T_n|$ exceeds this critical value. Likewise, one could alternatively consider one-sided alternatives. The asymptotic validity of both such tests follows directly from the results provided in the Supplemental Appendix.
\end{remark}
\begin{remark}\label{rem:confidence_interval}
As mentioned, most of the existing literature on subsampling-based inference about the mean of a regularly varying time series considers confidence intervals for the mean. In particular, with $\overline{X}_n = n^{-1}\sum_{t=1}^nX_t$ and
\begin{equation*}
T_{i,b_n}^c= \dfrac{\sum_{t=i}^{i+b_n-1} (X_t-\overline{X}_n)}{\big(\sum_{t=i}^{i+b_n-1}(X_t-\overline{X}_n)^2\big)^{1/2}},\qquad i=1,\dots,q_n,
\end{equation*}
define
\begin{equation*}
L_{n,b_n}^c(x)= q_n^{-1}\sum_{i=1}^{q_n}\mathbf{1}\big(T_{i,b_n}^c\leqslant x\big),\qquad x\in\mathbb{R},
\end{equation*}
and
\begin{equation*}
C_{n,b_n}^c(y)=\inf\{x:L_{n,b_n}^c(x)\leqslant y\},\qquad y\in(0,1).
\end{equation*}
Then an equal-tailed two-sided confidence interval for $\mathbb{E}[X_t]$ is given by
\begin{equation}\label{eq:subsamling_CI}
\mathrm{CI}_{b_n,n}=[\overline{X}_n+(\gamma_n/n)C_{n,b_n}^c(\eta/2) \;, \; \overline{X}_n+(\gamma_n/n)C_{n,b_n}^c(1-\eta/2)].
\end{equation}
Under the same assumptions as in Theorem \ref{cor:validity_subsampling}, but with $\mathbb{E}[X_t]$ potentially non-zero,
\begin{equation*}
\mathbb{P}\left(\mathbb{E}[X_t]\in \mathrm{CI}_{b_n,n}\right)\to 1-\eta ,\qquad n\to\infty,
\end{equation*}
that is, the confidence interval has asymptotic correct coverage.
Note that an alternative to the test in Algorithm \ref{algo} is achieved by rejecting $\mathsf{H}_0$ if $0\notin \mathrm{CI}_{b_n,n}$ \citet[p.54]{politis2012subsampling}. The validity of such a test follows directly from the technical proofs in the Supplemental Appendix.
\end{remark}
\section{Extensions}\label{sec:extensions}
In this section we consider extensions of the limit theory and related subsampling-based inference considered in Sections \ref{sec:limit_theory}--\ref{sec:subsampling}. In Section \ref{sec:infinite_mean} we consider testing the hypothesis $\mathsf{H}_0$ allowing for alternatives where $\mathbb{E}[|X_t|]=\infty$, that is, the mean may be infinite or may not exist. In Section \ref{sec:multivariate} we consider testing a hypothesis about a mean of a vector with particular attention to tests for superior predictive ability (SPA).
\subsection{Allowing for alternatives with $\mathbb{E}[|X_t|]=\infty$}\label{sec:infinite_mean}
A straightforward implication of Theorem \ref{theo:inf_var_jointCLT} is that when the tail index\footnote{The case $\kappa\in(0,1)$ may indeed be empirically relevant: \cite{kim2021forecast} find empirical evidence for such values of the tail index when considering loss differentials of realized volatility forecasts under squared loss, but do not address the inferential consequences. Their finding may also be motivated theoretically: \cite{gabaix2006_QJE} argue that stock market returns have tail index of $3$ (the ``cubic law of returns''). Consequently, realized volatility -- constructed by summing squared returns -- should have tail index $3/2$, so that squared realized volatility forecast errors may have tail index $\kappa=3/4<1$.} $\kappa\in(0,1)$,
\begin{equation*}
T_n\xrightarrow{d}\xi_\kappa\;/\zeta_{\kappa/2},\qquad n\to\infty.
\end{equation*}
Consequently, under alternatives where $\mathbb{E}[|X_t|]=\infty$, the test statistic $T_n$ does not diverge. Noting that the limiting distribution of the subsample statistics are identical to the one of $T_n$, with the critical values $C_{n,b_n}(\eta/2)$ and $C_{n,b_n}(1-\eta/2)$ determined according to Algorithm \ref{algo}, we have that
\begin{equation}
\mathbb{P}(T_n\notin[C_{n,b_n}(\eta/2),C_{n,b_n}(1-\eta/2)])\to \eta,\qquad n\to\infty;
\end{equation}
see Proposition \ref{subsample_special} in the Supplemental Appendix.
That is, the two-sided equal-tailed test is inconsistent when $\kappa\in(0,1)$. To the best of our knowledge, this issue has not been considered in the existing body of literature on subsampling techniques for heavy-tailed observations.
To construct a test that delivers non-trivial rejection frequencies under alternatives where $\kappa\in(0,1)$, we consider a modification of the original test statistic $T_n$. Define
\begin{equation}\label{eq:def_gamma_1}
\overline{\gamma}_{n}=\sum_{t=1}^n|X_t|.
\end{equation}
Under Assumption \ref{mix_coef_ass} or Assumptions \ref{ass_RV}--\ref{ass:AC+MX}, if $\mathbb{E}[|X_t|]<\infty$, then by the ergodic theorem $\overline{\gamma}_{n}/n\overset{\mathbb{P}}{\to}\mathbb{E}[|X_t|]$, as $n\to \infty$. In Section \ref{sec:proof_infinite_mean} we show that, with $a_n$ given by \eqref{eq:def_a_n}, $\overline{\gamma}_n/a_n $ has a strictly positive $\kappa$-stable limit if $\kappa\in(0,1)$, while $a_n/n\to\infty$ as $n\to\infty$. This implies that $\overline{\gamma}_{n}/n$ diverges for $\kappa\in(0,1)$ and motivates the modified statistic
\begin{equation}\label{eq:modified_statistic}
\widetilde{T}_n=(\overline{\gamma}_{n}/n)T_n.
\end{equation}
We have the following result.
\begin{proposition}
Let the statistic $\widetilde{T}_n$ be given by \eqref{eq:modified_statistic}.
\begin{enumerate}
\item Under Assumption \ref{mix_coef_ass} and $\mathbb{E}[X_0]=0$,
\begin{equation*}
\widetilde{T}_n\overset{d}{\to} \mathbb{E}[|X_t|]Z,\qquad n\to\infty,
\end{equation*}
with $Z$ a Gaussian random variable given in Theorem \ref{theo:rio}.
\item Under Assumptions \ref{ass_RV}--\ref{ass:AC+MX} with $\kappa\in(1,2)$ and $\mathbb{E}[X_0]=0$,
\begin{equation*}
\widetilde{T}_n\overset{d}{\to} \mathbb{E}[|X_t|] \frac{\xi_\kappa}{\zeta_{\kappa/2}},\qquad n\to\infty,
\end{equation*}
with $(\xi_\kappa,\zeta_{\kappa/2})$ given in Theorem \ref{theo:inf_var_jointCLT}.
\item Under Assumption \ref{mix_coef_ass} or Assumptions \ref{ass_RV}--\ref{ass:AC+MX}, if $\mathbb{E}[X_t]\neq0$,
\begin{equation}
|\widetilde{T}_n|\overset{\mathbb{P}}{\to}\infty,\qquad n\to\infty.
\end{equation}
\end{enumerate}
\end{proposition}
Critical values for a test based on $\widetilde{T}_n$ can be computed by the subsampling Algorithm \ref{algo_infinite_mean} in the Supplemental Appendix. Importantly, under $\kappa\in(0,1)$ we show in Section \ref{sec:proof_infinite_mean} that the subsample statistics $\widetilde{T}_{i,{b_n}}$ satisfy
\begin{equation}\label{eq:modified_subsampling_validity}
|\widetilde{T}_{i,{b_n}}| = o_\mathbb{P}(|\widetilde{T}_n|),
\end{equation}
under $\mathbb{E}[X_t]\neq0$, that is, the original statistic $\widetilde{T}_n$ diverges faster than the subsample statistics, ensuring non-trivial rejection frequencies. The finite-sample properties of a test based on the modified statistic $\widetilde{T}_n$ and subsampling-based critical values, are investigated in a simulation experiment in Section \ref{sim:infinite_mean} in the Supplemental Appendix. Notably, the simulations indicate that for processes with undefined mean ($\kappa <1$), the standard Diebold--Mariano test has rejection frequencies below the nominal level irrespective of the sample size, whereas the subsampling-based test in Algorithm \ref{algo_infinite_mean} achieves non-trivial rejection frequencies for sufficiently large sample sizes.
\subsection{Testing for superior predictive ability} \label{sec:multivariate}
In Section \ref{sec:appendix_multivariate} of the Supplemental Appendix we provide a multivariate version of Theorem \ref{theo:inf_var_jointCLT}. The result, Theorem \ref{theo:mikosch_multivariate}, may be used for testing a hypothesis about the mean vector of a multivariate heavy-tailed stationary time series. We illustrate this with an application to comparing multiple forecast series. Motivated by the work of \cite{hansen2005SPA}, see also \cite{white2000snooping}, suppose that we have $m+1\geqslant 2$ forecast series $(f_{j,t})_{t=1,\dots,n}$, $j=0,1\dots,m$, of a target variable $(Y_t)_{t=1,\dots,n}$ with $m$ fixed. Given a loss function $\mathcal{L}(\cdot , \cdot)$, let $\mathbf{X}_t=(X_{t,1},\dots,X_{t,m})'$ be defined by
\begin{equation}
X_{t,j}:=\mathcal{L}(Y_t,f_{0,t})-\mathcal{L}(Y_t,f_{j,t}),\qquad j=1,\dots,m,
\end{equation}
that is $X_{t,j}$ is the loss differential between forecast series $0$ and $j$ at time $t$. Given that $\mathbf{X}_t$ has time-invariant mean, the forecast series $0$ is said to have \textit{superior predictive ability} (SPA) over the other $m$ competing forecast series if
\begin{equation}\label{eq:def_H_SPA}
\mathsf{H}_{\mathrm{SPA}}:\mathbb{E}[\mathbf{X}_t]\leqslant \mathbf{0}
\end{equation}
is true, where the inequality holds entry-wise. \cite{hansen2005SPA} proposes a test for $\mathsf{H}_{\mathrm{SPA}}$ against
\begin{equation}\label{eq:H_alternative_SPA}
\mathsf{H}_{\mathrm{non-SPA}}: \mathbb{E}[X_{t,j}]>0 \ \text{for some} \ j\in\{1,\dots,m\}
\end{equation}
under the assumptions of $\Vert\mathbf{X}_t\Vert$ having bounded variance and suitable strong mixing coefficients of $(\mathbf{X}_t)_{t\in\mathbb{Z}}$, see Assumption \ref{ass:mixing_multivariate} in Section \ref{sec:appendix_multivariate} of the Supplemental Appendix. Here $\Vert \mathbf{x} \Vert = (\sum_{j=1}^mx_j^2)^{1/2}$ for a vector $\mathbf{x}=(x_1,\dots,x_m)'\in\mathbb{R}^m$. In order to account for the possibility of $\Vert\mathbf{X}_t\Vert$ having infinite variance, that is, being regularly varying with index $\kappa\in(1,2)$ (c.f., Assumption \ref{ass_RV_multivariate} in the Supplemental Appendix), we consider a modification of the test statistic proposed by \cite{hansen2005SPA}. Define
\begin{eqnarray}\label{eq:T_n_multivariate}
\mathbf{S}_n := \sum_{t=1}^n\mathbf{X}_t,\quad \tilde{\gamma}_{n}:=\sqrt{\sum_{t=1}^n\Vert \mathbf{X}_t\Vert^2},\quad \mathbf{T}_n := (T_{n,1},\dots,T_{n,m})'= \tilde{\gamma}_{n}^{-1}\mathbf{S}_n,
\end{eqnarray}
and the test statistic
\begin{equation}\label{eq:def_V_SPA}
V_n^{\mathrm{SPA}}:=\max\left[\max_{j=1,\dots,m}T_{n,j},0\right].
\end{equation}
We have the following result.
\begin{theorem}\label{thm:SPA}
Suppose that $(\mathbf{X}_t)_{t\in\mathbb{Z}}$ is stationary. Let $ V_n^{\mathrm{SPA}}$ be given by \eqref{eq:def_V_SPA}.
\begin{itemize}
\item [1.] Suppose that $\mathsf{H}_{\mathrm{SPA}}$ in \eqref{eq:def_H_SPA} holds.
\begin{itemize}
\item [a.] Under Assumption \ref{ass:mixing_multivariate},
\begin{equation*}
V_n^{\mathrm{SPA}} \overset{d}{\to} \frac{\max\left[\max_{j=1,\dots,m}\tilde{Z}_j,0\right]}{\sqrt{\mathbb{E}[\Vert \mathbf{X}_0 \Vert^2]}},\qquad n\to\infty.
\end{equation*}
where
\begin{equation}
\tilde{Z_{j}}=Z_j\mathbf{1}_{(\mathbb{E}[X_{0,j}]=0)},\qquad j=1,\dots,m,
\end{equation}
with $\mathbf{Z} = (Z_1,\dots,Z_m)'$ a Gaussian vector provided by Theorem \ref{thm:rio_multivariate} in the Appendix.
\item [b.] Under Assumptions \ref{ass_RV_multivariate}--\ref{AC,MX_multivariate} with tail index $\kappa\in(1,2)$,
\begin{equation*}
V_n^{\mathrm{SPA}} \overset{d}{\to} \frac{\max\left[\max_{j=1,\dots,m}\tilde{\xi}_{\kappa,j},0\right]}{\tilde{\zeta}_{\kappa,2}},\qquad n\to\infty.
\end{equation*}
where
\begin{equation}
\tilde{\xi}_{\kappa,j}=\xi_{\kappa,j}\mathbf{1}_{(\mathbb{E}[X_{0,j}]=0)}, \qquad j=1,\dots,m,
\end{equation}
with $\bm{\xi}_\kappa = (\xi_{\kappa,1},\dots,\xi_{\kappa,m})'$ a $\kappa$-stable random vector and $\tilde{\zeta}_{\kappa,2}^2$ a strictly positive $\kappa/2$-stable random variable both provided by Theorem \ref{theo:mikosch_multivariate} in the Appendix.
\end{itemize}
\item [2.] Suppose that $\mathsf{H}_{\mathrm{non-SPA}}$ in \eqref{eq:H_alternative_SPA} holds. Then under either Assumption \ref{ass:mixing_multivariate} or Assumptions \ref{ass_RV_multivariate}--\ref{AC,MX_multivariate} with tail index $\kappa\in(1,2)$,
\begin{equation*}
V_n^{\mathrm{SPA}}\overset{\mathbb{P}}{\to}\infty,\qquad n\to\infty.
\end{equation*}
\end{itemize}
\end{theorem}
As earlier, critical values for a test based on $V_n^{\mathrm{SPA}}$ can be computed by a straightforward multivariate extension of the subsampling procedure in Algorithm \ref{algo}.
\begin{remark}
Under Assumption \ref{ass_RV_multivariate} the entries of $\mathbf{X}_t$ may have different tail heaviness. In this case, the entries of $\mathbf{S}_n$ are of different stochastic orders whereas the normalizing quantity $\tilde{\gamma}_n$ has stochastic order determined by the \textit{smallest} of the tail indices, namely the index of regular variation of $\Vert \mathbf{X}_t \Vert$, $\kappa>0$. Consequently, the limiting distribution of $\mathbf{T}_n$ in \eqref{eq:T_n_multivariate} may be singular with respect to the Lebesgue measure on $\mathbb{R}^m$ when some entries of $\mathbf{X}_t$ have zero mean. This may affect the power properties of the test based on the statistic $V_n^{\mathrm{SPA}}$. An alternative test statistic, similar to the one considered by \cite{hansen2005SPA}, could be based on $\mathbf{S}_n \oslash \bm{\gamma}_n$ with $\bm{\gamma}_n:=(\gamma_{n,1},\dots,\gamma_{n,m})'$, $\gamma_{n,j}=\sqrt{\sum_{t=1}^nX_{t,j}^2}$, and $\oslash$ denoting entry-wise division. It remains an open task to derive stable limit results for $(\mathbf{S}_n,\bm{\gamma}_n)$ under potentially heterogeneous tail heaviness across the entries of $\mathbf{X}_t$ and, consequently, with entry-wise different rates of convergence.
\end{remark}
\section{Monte-Carlo simulations}\label{sec:Monte_Carlo}
In this section, we investigate the finite-sample properties of the standard Diebold--Mariano test as well as the test based on subsampling Algorithm \ref{algo} with block length $b_n$ given by \eqref{eq:b_n_choice}. We refer to Section \ref{sim:infinite_mean} in the Supplemental Appendix for additional results for the Diebold--Mariano test and the subsampling based test allowing for $\mathbb{E}[|X_t|]=\infty$ as considered in Section \ref{sec:infinite_mean}.
Throughout, we consider a DGP for $X_t$ given in terms of the AR(1) process
\begin{equation}\label{eq:sim_AR}
X_t=\delta+0.5\; X_{t-1}+Z_t,\qquad t=1,\dots,n,
\end{equation}
with iid noise $Z_t$.\footnote{The initial $X_0$ is chosen random using a burn-in period of $10^4$ observations initiated at $X_{-10000}=0$.}\footnote{We also considered different values of the autoregressive coefficient ranging from $0.1$ to $0.9$. The simulation results were qualitatively the same as for the ones reported here with $0.5$.} Initially, we consider two different distributions for $Z_t$: a symmetric stable and an asymmetric stable, denoted \textrm{Stable}$(\kappa,0,1,0)$ and \textrm{Stable}$\left(\kappa,4/5,1,4\tan(\kappa\pi/2)/5\right)$, respectively, with $\kappa\in(1,2)$. Both distributions imply that $\mathbb{E}[Z_t]=0$, and, consequently, $\mathbb{E}[X_t]=2\delta$ so that the hypothesis $\mathsf{H}_0$ in \eqref{eq:zero_mean} is true if and only if $\delta=0$. Moreover, $Z_t$ is regularly varying with index $\kappa$ for both choices of distribution. The symmetric distribution has tail balance coefficients $p_+=p_-=1/2$, whereas the asymmetric distribution has $p_+=1-p_-=9/10$. By Proposition \ref{prop:AR1 - RV} $(X_t)_{t\in\mathbb{Z}}$ is regularly varying with index $\kappa$.
Moreover, by Proposition \ref{prop:AR_lim_moments}, the limiting distribution of the statistic $T_n$, $\xi_\kappa/\zeta_{\kappa/2}$, has a continuous density and finite moments of all orders for both distributional specifications of $Z_t$. In the asymmetric case ($p_+ = 0.9$), $\xi_\kappa$ has strictly positive skewness parameter $\beta$, and $\xi_\kappa/\zeta_{\kappa/2}$ has strictly negative mean.
In terms of the Diebold--Mariano test, we make use of the statistic
\begin{align}\label{eq:T_DM}
T_n^{\mathrm{DM}}=\frac{S_n}{\sqrt{n}\hat{\sigma}_{n,\mathrm{NW}}},
\end{align}
with $S_n$ given by \eqref{eq:def_S_gamma} and with $\hat{\sigma}^2_{n,\mathrm{NW}}$ the standard \cite{newey_west_1987} HAC variance estimator with lag length given by $\lfloor4(n/100)^{2/9}\rfloor$, as considered in \citet[Section IV]{newey_west_1994_automatic}. The test rejects the hypothesis $\mathsf{H}_0$ at a 5\% nominal level if $|T_n^{\mathrm{DM}}|>1.96$. Apart from the centering with respect to the sample mean in $\hat{\sigma}^2_{n,\mathrm{NW}}$ and the increasing lag length, the test statistic is identical to $T_n^{\mathrm{HAC}}$ in \eqref{eq:T_HAC}. Consequently, from Proposition \ref{prop:AR_HAC} we expect that, under $\mathsf{H}_0$, $T_n^{\mathrm{DM}}$ is potentially poorly approximated by a standard normal distribution when $X_t$ has tail index $\kappa <2$.
\subsection{Rejection frequencies under $\mathsf{H}_0$}
Table \ref{tab:RFs under null} contains rejection frequencies under $\mathsf{H}_0$ for the two tests for different values of the tail index $\kappa\in(1,2)$. In the case where the innovation $Z_t$ has a symmetric distribution, Panel A of Table \ref{tab:RFs under null} states that the rejection frequencies of the two tests are quite comparable and fairly close to the nominal level of 5\% across tail indices and sample sizes.
Turning to Panel B, where the innovations exhibit extremal skewness, the Diebold--Mariano test is heavily distorted, with rejection frequencies remarkably above the nominal level for $\kappa \leqslant 1.5$. This holds irrespective of the sample size, and the distortion increases with the tail-heaviness ($\kappa \downarrow 1$). This is likely explained by the asymmetry of the limiting distribution $\xi_\kappa/\zeta_{\kappa/2}$, which suggests that critical values based on the (symmetric) standard normal distribution are unreliable. The subsampling-based test does not perform as well as in the symmetric case, but its rejection frequencies are generally substantially closer to the nominal level than those of the Diebold--Mariano test, and converge toward it as the sample size increases. For tail indices $\kappa \geqslant 1.5$, the rejection frequencies of the subsampling-based test appear reasonable across all sample sizes.
\begin{table}[H]
\centering
\caption{Rejection frequencies (in percent) under $\mathsf{H}_0$ for the Diebold--Mariano test and the subsampling-based test in Algorithm \ref{algo}. The Diebold--Mariano test is based on the statistic $T_n^{\mathrm{DM}}$ in \eqref{eq:T_DM}. The DGP $(X_t)_{t=1,\dots,n}$ is given by \eqref{eq:sim_AR} with $\delta=0$. In Panel A the iid noise $(Z_t)_{t=1,\dots,n}$ is $\textrm{Stable}(\kappa,0,1,0)$, and in Panel B it is \textrm{Stable}$\left(\kappa,4/5,1,4\tan(\kappa\pi/2)/5\right)$. All tests are carried out at a 5\% nominal level. The rejection frequencies are based on $M=10^4$ Monte-Carlo replications.}\label{tab:RFs under null}
\resizebox{0.95\textwidth}{!}{
\begin{tabular}{l|ccccc|ccccc}
\hline\hline
& \multicolumn{5}{c|}{Diebold--Mariano} & \multicolumn{5}{c}{Subsampling Algorithm \ref{algo}} \\
\hline
$n$: & 1000 & 2000 & 5000 & 10000 & 100000 & 1000 & 2000 & 5000 & 10000 & 100000 \\
\hline\hline
\multicolumn{11}{c}{Panel A - \textrm{Stable}$(\kappa,0,1,0)$ noise.} \\
\hline
$\kappa=1.1$ & 4.9 & 4.4 & 3.8 & 3.5 & 3.0 & 7.3 & 6.5 & 5.7 & 5.6 & 5.2 \\
$\kappa=1.3$ & 6.1 & 5.5 & 5.0 & 4.2 & 4.1 & 6.7 & 6.2 & 6.0 & 5.2 & 5.3 \\
$\kappa=1.5$ & 6.8 & 5.3 & 5.6 & 5.0 & 4.5 & 5.7 & 4.6 & 4.9 & 4.7 & 4.8 \\
$\kappa=1.7$ & 7.4 & 6.8 & 6.3 & 6.0 & 4.8 & 5.0 & 4.5 & 4.4 & 4.6 & 4.2 \\
$\kappa=1.9$ & 7.6 & 7.2 & 6.7 & 6.3 & 5.5 & 4.5 & 4.0 & 3.8 & 4.0 & 3.7 \\
\hline\hline
\multicolumn{11}{c}{Panel B - \textrm{Stable}$\left(\kappa,4/5,1,4\tan(\kappa\pi/2)/5\right)$ noise.} \\
\hline
$\kappa=1.1$ & 71.2 & 71.8 & 71.4 & 71.0 & 70.6 & 50.4 & 45.1 & 37.8 & 31.2 & 15.1 \\
$\kappa=1.3$ & 35.2 & 36.1 & 34.9 & 34.3 & 33.0 & 15.3 & 12.7 & 9.4 & 8.1 & 6.0 \\
$\kappa=1.5$ & 18.0 & 18.0 & 17.2 & 16.2 & 16.3 & 8.3 & 7.1 & 5.8 & 5.6 & 4.9 \\
$\kappa=1.7$ & 11.1 & 10.1 & 9.4 & 9.3 & 8.1 & 6.1 & 5.0 & 4.7 & 4.8 & 4.3 \\
$\kappa=1.9$ & 8.0 & 6.8 & 7.2 & 6.2 & 6.0 & 4.6 & 3.5 & 3.8 & 3.6 & 4.0 \\
\hline\hline
\end{tabular}
}
\end{table}
\subsection{Rejection frequencies under alternatives}\label{sec:algo_power}
Figure \ref{fig:power_finite_mean} contains rejection frequencies for the subsampling-based test when the null hypothesis is violated. We consider DGPs given by \eqref{eq:sim_AR} with $Z_t$ \textrm{Stable}$(\kappa,0,1,0)$-distributed and $\delta\in(0,2.5]$ such that $\mathbb{E}[X_t]\in(0,5]$. As one would expect, the rejection frequencies are increasing in the sample size. We also note that the rejection frequencies tend to increase quite slowly in the sample size and in $\mathbb{E}[X_t]$ as $\kappa$ tends to one. This is likely explained by the fact that for $\kappa$ close to one, we are near a setting where $X_t$ has undefined mean, as considered in detail in Section \ref{sec:infinite_mean}. When $X_t$ has undefined mean ($\kappa\in(0,1)$) the test statistic $T_n$ does not diverge irrespective of the value of $\delta$, and the subsampling- based test is not expected to provide non-trivial rejection frequencies.
\begin{figure}[H]
\centering
\includegraphics[width=1\linewidth]{Power_Algo1_vs2.png}
\caption{Rejection frequencies under alternatives $\mathbb{E}[X_t]\neq 0$ for the subsampling-based test in Algorithm \ref{algo}. The DGP $(X_t)_{t=1,\dots,n}$ is given by \eqref{eq:sim_AR} with $\delta\in[0,2.5]$ and with iid noise $(Z_t)_{t=1,\dots,n}$ is $\textrm{Stable}(\kappa,0,1,0)$-distributed, $\kappa\in(1,2)$. The test is carried out at a 5\% nominal level. The x-axes indicate the values of $\mathbb{E}[X_t]=2\delta$. The rejection frequencies are based on $M=10^4$ Monte-Carlo replications.}\label{fig:power_finite_mean}
\end{figure}
We end this section by comparing the performances of the Diebold--Mariano and the subsampling-based tests when $X_t$ has finite variance, i.e., when the Diebold--Mariano test is expected to work. To do so, we modify the DGP so that $Z_t$ has a Student's $t$-distribution with $\kappa > 2$ degrees of freedom: $Z_t$ either has a standard (symmetric) Student's $t$-distribution ($\mathrm{Student}(\kappa)$) or a skewed Student's $t$-distribution ($\mathrm{Skt}(\kappa)$).\footnote{We refer to \cite{fernandez1998bayesian} for details about the skewed Student's $t$-distribution.} For the skewed case, $Z_t$ is scaled such that $V(Z_t) = \kappa/(\kappa-2)$, which is also the variance of the symmetric version. Both distributions are regularly varying with index $\kappa$ and have $\mathbb{E}[Z_t] = 0$. For the symmetric case, the tail balancing coefficients are $p_+ = p_- = 0.5$, and the skewed distribution is chosen such that $p_+ = 1 - p_- = 0.9$. Since $\kappa > 2$, the DGPs obey Theorem \ref{theo:rio}, so that $T_n$ has a Gaussian limiting distribution; we expect the same to hold for the Diebold--Mariano statistic $T_n^{\mathrm{DM}}$.
Figure \ref{fig:Algo vs Newey} contains rejection frequencies for the Diebold--Mariano and subsampling-based tests when the null hypothesis is violated. The first row corresponds to symmetric $t$-distribution and the second to the asymmetric one. For the heaviest-tailed symmetric distributions ($\kappa = 2.1, 2.5$), the Diebold--Mariano test yields larger rejection frequencies than the subsampling-based test, whereas the opposite holds under asymmetric distributions. Under lighter tails ($\kappa \geqslant 3$), the rejection frequencies of the two tests are quite similar, especially for large sample sizes ($n = 5000$). Overall, the two tests exhibit comparable power under finite-variance distributions.
\begin{figure}[H]
\centering
\includegraphics[width=1\linewidth]{Power_Algo1_Finite_Variance.png}
\caption{Rejection frequencies under alternatives for tests based on Algorithm \ref{algo} and the Diebold--Mariano test. The tests are carried out at a 5\% nominal level. The DGP is given by \eqref{eq:sim_AR} with $(Z_t)$ iid and $\delta \in [0,0.7]$. In the first row $Z_t$ is $\textrm{Student}(\kappa)$-distributed. In the second row $Z_t$ is \textrm{Skt}$(\kappa)$-distributed. The x-axes indicate the values of $\mathbb{E}[X_t]=2\delta$. The rejection frequencies are based on $M=10^4$ Monte-Carlo replications.}
\label{fig:Algo vs Newey}
\end{figure}
\section{Empirical illustration: Comparisons of risk forecasts for emerging-market foreign exchange rates }\label{sec:empirical}
Quantification and hedging of foreign currency risk are important for international portfolio management\footnote{See, e.g., \cite{christensen_varneskov_2021_FX} and the references therein.}, and we consider in this section an application of Algorithm \ref{algo} to comparing risk forecasts for foreign exchange (FX) rate returns. FX rates are typically volatile, and their returns are heavy-tailed. For instance, \cite{ibragimov2013emerging_markets} find that certain emerging-market FX rate returns are regularly varying, and standard confidence intervals for their tail index include values below two; that is, returns may have infinite variance. Among such FX rates is the South Korean won against U.S. dollars (KRWUSD), which we focus on here.
Let $Y_t$ denote the FX rate log return on day $t$. Given some information set $\mathcal{F}_{t-1}$ available at time $t-1$ denote by $Q_{t-1,\tau}$ the $\tau$-quantile of the conditional distribution of $Y_t$ given $\mathcal{F}_{t-1}$:
\begin{eqnarray}\label{eq:def_VaR}
\mathbb{P}(Y_{t}\leqslant Q_{t-1,\tau}|\mathcal{F}_{t-1})=\tau,\quad Q_{t-1,\tau}\in \mathcal{F}_{t-1},\quad \tau\in(0,1).
\end{eqnarray}
In the context of financial risk management, this quantile is the Value-at-Risk (VaR) at risk level $\tau$. Following Example \ref{ex:EPA} one may evaluate the time $t$ VaR forecast using the tick loss function,
\begin{equation}\label{eq:def_tick_loss}
\mathcal{L}_{\mathrm{tick},\tau}(Y_t,Q_{t-1,\tau})= (\tau - \mathbf{1}\{Y_t-Q_{t-1,\tau}<0\})(Y_t-Q_{t-1,\tau}).
\end{equation}
A common approach to quantifying VaR is by means of a GARCH-type model: suppose that $Y_t = \sigma_t Z_t$, with $\sigma_t > 0$ a function of past returns and $Z_t$ a mean-zero variable independent of $\sigma_t \in \mathcal{F}_{t-1}$. Then the VaR of $Y_t$ is given by $Q_{t-1,\tau}=c_{\tau}\sigma_t$, where the constant $c_\tau$ is the $\tau$-quantile of the distribution of $Z_t$. Consequently, the tick loss is
\begin{equation*}
\mathcal{L}_{\mathrm{tick},\tau}(Y_t,Q_{t-1,\tau}) = (\tau - \mathbf{1}
\{Z_t<c_{\tau}\})(Z_t-c_\tau)\sigma_t.
\end{equation*}
For standard GARCH-type processes, $\sigma_t^2$ can be written as a stochastic recurrence equation (see, e.g., Section \ref{sec:SRE} in the Supplemental Appendix), and under mild conditions, $\sigma_t$ is regularly varying. By Breiman's Lemma, the tick loss $\mathcal{L}_{\mathrm{tick},t}(Y_t,Q_{t-1,\tau})$ is also regularly varying with index $\kappa$, and, likewise, loss differentials may be regularly varying.
Similar to \citet[Section 5]{patton2019_es_var}, we compare the VaR forecasts from ten different methods. The methods can be divided into three groups: rolling window-based estimation, GARCH models, and score-driven models. The rolling window methods, labeled "RW-$H$", compute the risk forecasts using a simple rolling window of $H$ days with $H\in\{125,250,500\}$. Three of the GARCH-based methods, labeled "G-N", "G-Skt", and "G-EDF", compute the forecasts using a standard GARCH(1,1) model with, respectively, the $\mathrm{N}(0,1)$-distribution, a skewed Student's $t$-distribution, and the empirical distribution of the in-sample standardized residuals. The method "G-FZ" computes the risk forecast based on a GARCH(1,1) model estimated under a so-called Fissler-Ziegel loss function; see \citet[Section 2.5.1]{patton2019_es_var}. The three methods based on score-driven models are labeled "FZ-2F", "FZ-1F", and "Hybrid" and are described in \citet[Section 2]{patton2019_es_var}.\footnote{Along the lines of \citet[Section 5]{patton2019_es_var}, the conditional mean of the returns is given by a constant for the GARCH and score-driven models.} Our dataset\footnote{The data were retrieved from the Board of Governors of the Federal Reserve System webpage, \url{https://www.federalreserve.gov/datadownload/}.} contains daily continuously compounded returns on the KRWUSD rate from 4 January 1990 to 26 December 2025. The first ten years of data (4 January 1990 to 31 December 1999) are used for the estimation of the model parameters. These estimates are kept fixed for the entire out-of-sample period (1 January 2000 to 26 December 2025). As in \citet[Section 5]{patton2019_es_var}, we remove all observations with returns exactly zero, yielding a total of 2225 in-sample observations for parameter estimation and $n=6333$ out-of-sample observations for forecast comparisons. Throughout, the risk level $\tau$ is set to 0.05.
For any two competing VaR forecasts $Q_{1, t-1,\tau}$ and $Q_{2, t-1,\tau}$, the time $t$ tick loss differential is
\begin{equation}\label{eq:tickloss_diff}
X_t = \mathcal{L}_{\mathrm{tick},\tau}(Y_t,Q_{1, t-1,\tau}) - \mathcal{L}_{\mathrm{tick},\tau}(Y_t,Q_{2, t-1,\tau}).
\end{equation}
Following Example \ref{ex:EPA}, the two series have equal predictive ability under $\mathsf{H}_0:\mathbb{E}[X_t]=0$, and we test this EPA hypothesis for all pairs of the 10 forecast series. Panel A of Table \ref{tab:EPA_KRWUSD_tickLoss} presents Diebold--Mariano $t$-statistics on the loss differentials\footnote{As in \citet[Section 5]{patton2019_es_var} the statistics are based on the Newey-West estimator with a fixed lag-length of 20 days.}, whereas Panel B presents the statistic $T_n$ on the differentials along with subsampling-based critical values computed by Algorithm \ref{algo} with block size given by \eqref{eq:b_n_choice}. The statistics are computed as "row method minus column method", so that a negative number indicates that the row method has a smaller average loss than the column method. All grey-shaded entries imply a rejection of the EPA hypothesis at a 5\% nominal level. Notably, Panel B has fewer grey-shaded entries than Panel A, that is, the EPA hypothesis is rejected for fewer pairs when relying on the subsampling-based test. For instance, based on the reported Diebold--Mariano statistics and their sign, one would tend to conclude that the forecasts based on the G-FZ method are more accurate than each of the nine competing forecast series. In contrast, applying the subsampling-based test, the EPA hypothesis is not rejected when comparing the G-FZ method with any of the rolling window-based forecasts. A plausible explanation for the discrepancy between the outcomes of the two tests can be found in Figure \ref{fig:hill_KRWUSD_tickloss}, containing the loss differentials between the RW-500 and G-FZ forecasts and a Hill plot of the absolute value of the loss differentials. The loss differentials have an estimated tail index around 1.5, indicating that the series has infinite variance but a well-defined mean. We also note that the loss differential series appears to exhibit extremal skewness. To see this, the solid red lines indicate plus/minus the empirical 99th percentile of the absolute value of the loss differentials. Most of the loss differentials exceeding these thresholds (in magnitude) are positive, suggesting that rolling window-based forecasts are much more prone to exhibiting extreme losses than their counterparts. Following \citet[Section 2.1]{davis2018tail_process_inference} the tail-balance coefficient $p_+$ in \eqref{RV_marg} can be estimated by the proportion of positive loss differentials exceeding a large threshold. Using the aforementioned empirical 99th percentile as the threshold, we obtain an estimate $\hat{p}_+\approx 0.86$. This extremal skewness is also reflected in the reported subsampled equal-tailed critical values, suggesting that the statistic $T_n$, and hence the loss differentials, has a skewed distribution. Recall from the simulation experiments in Section \ref{sec:Monte_Carlo} that in such a setting, even with $n\geqslant5000$, the Diebold--Mariano test tends to over-reject, whereas the subsampling-based test has rejection frequencies closer to its nominal level.
\begin{figure}[H]
\centering
\includegraphics[width=1\linewidth]{final_tick_vs2.png}
\caption{Left: The series of loss differentials when comparing $\mathrm{VaR}$-forecasts based on rolling-window and GARCH approaches, respectively. The red lines indicate plus/minus the empirical 99th percentile of the series in absolute values. Right: Hill plot based on the absolute values, where the dashed red line indicates $\kappa =1.5$.}\label{fig:hill_KRWUSD_tickloss}
\end{figure}
\newpage
\begin{table}[H]
\centering
\caption{Test statistics for EPA tests}\label{tab:EPA_KRWUSD_tickLoss}
\scriptsize
\setlength{\tabcolsep}{3pt}
\resizebox{\linewidth}{!}{
\begin{tabular}{lcccccccccc}
\multicolumn{11}{c}{\small{A: Diebold--Mariano $t$-statistics, $T_n^{\mathrm{DM}}$}}\\
\toprule
& \textbf{RW-125} & \textbf{RW-250} & \textbf{RW-500} & \textbf{G-N} & \textbf{G-Skt} & \textbf{G-EDF} & \textbf{FZ-2F} & \textbf{FZ-1F} & \textbf{G-FZ} & \textbf{Hybrid} \\
\midrule
\textbf{RW-125} & &
\cellcolor{siggray}{-2.177} &
\cellcolor{siggray}{-2.774} &
\cellcolor{siggray}{2.089} &
1.028 &
1.879 &
\cellcolor{siggray}{-3.206} &
-0.475 &
\cellcolor{siggray}{2.962} &
0.879 \\
\textbf{RW-250} &
\cellcolor{siggray}{2.177} & &
\cellcolor{siggray}{-2.470} &
\cellcolor{siggray}{2.525} &
1.769 &
\cellcolor{siggray}{2.382} &
-0.860 &
0.853 &
\cellcolor{siggray}{3.197} &
1.558 \\
\textbf{RW-500} &
\cellcolor{siggray}{2.774} &
\cellcolor{siggray}{2.470} & &
\cellcolor{siggray}{3.212} &
\cellcolor{siggray}{2.691} &
\cellcolor{siggray}{3.118} &
1.156 &
\cellcolor{siggray}{2.181} &
\cellcolor{siggray}{3.764} &
\cellcolor{siggray}{2.437} \\
\textbf{G-N} &
\cellcolor{siggray}{-2.089} &
\cellcolor{siggray}{-2.525} &
\cellcolor{siggray}{-3.212} & &
\cellcolor{siggray}{-4.860} &
\cellcolor{siggray}{-3.113} &
\cellcolor{siggray}{-4.165} &
\cellcolor{siggray}{-3.471} &
\cellcolor{siggray}{2.531} &
\cellcolor{siggray}{-2.474} \\
\textbf{G-Skt} &
-1.028 &
-1.769 &
\cellcolor{siggray}{-2.691} &
\cellcolor{siggray}{4.860} & &
\cellcolor{siggray}{5.629} &
\cellcolor{siggray}{-3.361} &
-1.850 &
\cellcolor{siggray}{4.864} &
-0.060 \\
\textbf{G-EDF} &
-1.879 &
\cellcolor{siggray}{-2.382} &
\cellcolor{siggray}{-3.118} &
\cellcolor{siggray}{3.113} &
\cellcolor{siggray}{-5.629} & &
\cellcolor{siggray}{-4.061} &
\cellcolor{siggray}{-3.103} &
\cellcolor{siggray}{3.203} &
-1.781 \\
\textbf{FZ-2F} &
\cellcolor{siggray}{3.206} &
0.860 &
-1.156 &
\cellcolor{siggray}{4.165} &
\cellcolor{siggray}{3.361} &
\cellcolor{siggray}{4.061} & &
1.799 &
\cellcolor{siggray}{5.093} &
\cellcolor{siggray}{2.555} \\
\textbf{FZ-1F} &
0.475 &
-0.853 &
\cellcolor{siggray}{-2.181} &
\cellcolor{siggray}{3.471} &
1.850 &
\cellcolor{siggray}{3.103} &
-1.799 & &
\cellcolor{siggray}{5.495} &
\cellcolor{siggray}{2.047} \\
\textbf{G-FZ} &
\cellcolor{siggray}{-2.962} &
\cellcolor{siggray}{-3.197} &
\cellcolor{siggray}{-3.764} &
\cellcolor{siggray}{-2.531} &
\cellcolor{siggray}{-4.864} &
\cellcolor{siggray}{-3.203} &
\cellcolor{siggray}{-5.093} &
\cellcolor{siggray}{-5.495} & &
\cellcolor{siggray}{-5.005} \\
\textbf{Hybrid} &
-0.879 &
-1.558 &
\cellcolor{siggray}{-2.437} &
\cellcolor{siggray}{2.474} &
0.060 &
1.781 &
\cellcolor{siggray}{-2.555} &
\cellcolor{siggray}{-2.047} &
\cellcolor{siggray}{5.005} & \\
\midrule
\\[1em]
\multicolumn{11}{c}{\small{B: $T_n$ and equal-tailed critical values}} \\
\toprule
& \textbf{RW-125} & \textbf{RW-250} & \textbf{RW-500} & \textbf{G-N} & \textbf{G-Skt} & \textbf{G-EDF} & \textbf{FZ-2F} & \textbf{FZ-1F} & \textbf{G-FZ} & \textbf{Hybrid} \\
\midrule
\textbf{RW-125} &
&
$\underset{\{-6.06,\,2.76\}}{-3.24}$ &
$\underset{\{-7.47,\,2.39\}}{-4.55}$ &
$\underset{\{-2.05,\,2.49\}}{2.38}$ &
$\underset{\{-1.98,\,2.07\}}{1.18}$ &
$\underset{\{-1.93,\,2.32\}}{2.16}$ &
\cellcolor{siggray}$\underset{\{-2.28,\,1.58\}}{-3.30}$ &
$\underset{\{-3.43,\,2.43\}}{-0.62}$ &
$\underset{\{-1.98,\,3.61\}}{3.42}$ &
$\underset{\{-4.22,\,3.38\}}{1.01}$ \\
\textbf{RW-250} &
$\underset{\{-2.76,\,6.06\}}{3.24}$ &
&
$\underset{\{-8.01,\,4.37\}}{-4.24}$ &
$\underset{\{-2.03,\,4.29\}}{3.17}$ &
$\underset{\{-2.07,\,3.52\}}{2.26}$ &
$\underset{\{-2.01,\,3.96\}}{3.02}$ &
$\underset{\{-2.08,\,2.27\}}{-1.02}$ &
$\underset{\{-3.66,\,3.50\}}{1.07}$ &
$\underset{\{-2.11,\,4.89\}}{3.98}$ &
$\underset{\{-3.42,\,5.02\}}{1.89}$ \\
\textbf{RW-500} &
$\underset{\{-2.39,\,7.47\}}{4.55}$ &
$\underset{\{-4.37,\,8.01\}}{4.24}$ &
&
$\underset{\{-2.04,\,5.14\}}{4.54}$ &
$\underset{\{-1.93,\,4.21\}}{3.91}$ &
$\underset{\{-1.93,\,4.83\}}{4.47}$ &
$\underset{\{-1.96,\,3.02\}}{1.69}$ &
$\underset{\{-2.93,\,5.20\}}{2.85}$ &
$\underset{\{-2.09,\,6.29\}}{5.23}$ &
$\underset{\{-4.60,\,7.02\}}{3.20}$ \\
\textbf{G-N} &
$\underset{\{-2.49,\,2.05\}}{-2.38}$ &
$\underset{\{-4.29,\,2.03\}}{-3.17}$ &
$\underset{\{-5.14,\,2.04\}}{-4.54}$ &
&
\cellcolor{siggray}$\underset{\{-1.75,\,0.98\}}{-4.42}$ &
\cellcolor{siggray}$\underset{\{-1.56,\,1.44\}}{-2.72}$ &
\cellcolor{siggray}$\underset{\{-1.79,\,1.15\}}{-4.52}$ &
\cellcolor{siggray}$\underset{\{-2.77,\,1.39\}}{-3.95}$ &
\cellcolor{siggray}$\underset{\{-2.02,\,1.89\}}{2.24}$ &
$\underset{\{-2.84,\,1.09\}}{-2.37}$ \\
\textbf{G-Skt} &
$\underset{\{-2.07,\,1.98\}}{-1.18}$ &
$\underset{\{-3.52,\,2.07\}}{-2.26}$ &
$\underset{\{-4.21,\,1.93\}}{-3.91}$ &
\cellcolor{siggray}$\underset{\{-0.98,\,1.75\}}{4.42}$ &
&
\cellcolor{siggray}$\underset{\{-0.81,\,1.88\}}{5.25}$ &
\cellcolor{siggray}$\underset{\{-1.97,\,1.91\}}{-3.88}$ &
$\underset{\{-2.68,\,1.38\}}{-1.99}$ &
\cellcolor{siggray}$\underset{\{-1.76,\,1.79\}}{4.14}$ &
$\underset{\{-2.36,\,1.22\}}{-0.05}$ \\
\textbf{G-EDF} &
$\underset{\{-2.32,\,1.93\}}{-2.16}$ &
$\underset{\{-3.96,\,2.01\}}{-3.02}$ &
$\underset{\{-4.83,\,1.93\}}{-4.47}$ &
\cellcolor{siggray}$\underset{\{-1.44,\,1.56\}}{2.72}$ &
\cellcolor{siggray}$\underset{\{-1.88,\,0.81\}}{-5.25}$ &
&
\cellcolor{siggray}$\underset{\{-1.85,\,1.47\}}{-4.51}$ &
\cellcolor{siggray}$\underset{\{-2.75,\,1.31\}}{-3.49}$ &
\cellcolor{siggray}$\underset{\{-2.14,\,1.77\}}{2.75}$ &
$\underset{\{-2.65,\,1.14\}}{-1.67}$ \\
\textbf{FZ-2F} &
\cellcolor{siggray}$\underset{\{-1.58,\,2.28\}}{3.30}$ &
$\underset{\{-2.27,\,2.08\}}{1.02}$ &
$\underset{\{-3.02,\,1.96\}}{-1.69}$ &
\cellcolor{siggray}$\underset{\{-1.15,\,1.79\}}{4.52}$ &
\cellcolor{siggray}$\underset{\{-1.91,\,1.97\}}{3.88}$ &
\cellcolor{siggray}$\underset{\{-1.47,\,1.85\}}{4.51}$ &
&
\cellcolor{siggray}$\underset{\{-1.40,\,1.63\}}{2.02}$ &
\cellcolor{siggray}$\underset{\{-0.42,\,1.76\}}{5.33}$ &
\cellcolor{siggray}$\underset{\{-1.15,\,1.45\}}{2.62}$ \\
\textbf{FZ-1F} &
$\underset{\{-2.43,\,3.43\}}{0.62}$ &
$\underset{\{-3.50,\,3.66\}}{-1.07}$ &
$\underset{\{-5.20,\,2.93\}}{-2.85}$ &
\cellcolor{siggray}$\underset{\{-1.39,\,2.77\}}{3.95}$ &
$\underset{\{-1.38,\,2.68\}}{1.99}$ &
\cellcolor{siggray}$\underset{\{-1.31,\,2.75\}}{3.49}$ &
\cellcolor{siggray}$\underset{\{-1.63,\,1.40\}}{-2.02}$ &
&
\cellcolor{siggray}$\underset{\{-1.05,\,3.35\}}{6.56}$ &
$\underset{\{-2.75,\,2.94\}}{2.30}$ \\
\textbf{G-FZ} &
$\underset{\{-3.61,\,1.98\}}{-3.42}$ &
$\underset{\{-4.89,\,2.11\}}{-3.98}$ &
$\underset{\{-6.29,\,2.09\}}{-5.23}$ &
\cellcolor{siggray}$\underset{\{-1.89,\,2.02\}}{-2.24}$ &
\cellcolor{siggray}$\underset{\{-1.79,\,1.76\}}{-4.14}$ &
\cellcolor{siggray}$\underset{\{-1.77,\,2.14\}}{-2.75}$ &
\cellcolor{siggray}$\underset{\{-1.76,\,0.42\}}{-5.33}$ &
\cellcolor{siggray}$\underset{\{-3.35,\,1.05\}}{-6.56}$ &
&
\cellcolor{siggray}$\underset{\{-4.08,\,0.62\}}{-5.21}$ \\
\textbf{Hybrid} &
$\underset{\{-3.38,\,4.22\}}{-1.01}$ &
$\underset{\{-5.02,\,3.42\}}{-1.89}$ &
$\underset{\{-7.02,\,4.60\}}{-3.20}$ &
$\underset{\{-1.09,\,2.84\}}{2.37}$ &
$\underset{\{-1.22,\,2.36\}}{0.05}$ &
$\underset{\{-1.14,\,2.65\}}{1.67}$ &
\cellcolor{siggray}$\underset{\{-1.45,\,1.15\}}{-2.62}$ &
$\underset{\{-2.94,\,2.75\}}{-2.30}$ &
\cellcolor{siggray}$\underset{\{-0.62,\,4.08\}}{5.21}$ &
\\
\bottomrule
\end{tabular}}
\end{table}
\vspace{-0.5em}
\begin{singlespace}
\noindent {\footnotesize Notes: Panel A contains $t$-statistics, $T_n^{\mathrm{DM}}$, from the Diebold--Mariano test of EPA, using the tick loss function with risk level $\tau=0.05$, over the out-of-sample period from January 2000 to December 2025, for ten different forecasting methods. Method "RW-$H$" computes the risk forecasts using a simple rolling window of $H$ days. Methods "G-N", "G-Skt", and "G-EDF" compute the forecasts using a standard GARCH(1,1) model with, respectively, the $N(0,1)$-distribution, a skewed Student's $t$-distribution, and the empirical distribution of the in-sample standardized residuals. The "FZ-2F", "FZ-1F", "G-FZ", and "Hybrid" methods are described in \citet[Section 2]{patton2019_es_var}. The statistics are based on the Newey-West HAC estimator with 20 lags. A positive value indicates that the row method has higher average loss than the column method. Values exceeding 1.96 in absolute value (grey shaded) imply a rejection of the EPA hypothesis at a 5\% nominal level. Values along the main diagonal are all, by convention, identically zero and are omitted. Panel B contains $t$-statistics, $T_n$, for the null of EPA along with equal-tailed critical values (in brackets) at 5\% nominal level computed according to Algorithm \ref{algo}. Values outside of the interval formed by the critical values (grey shaded) imply a rejection of the EPA hypothesis.}
\end{singlespace}
\section{Conclusion}
Many key variables in economics and finance exhibit heavy-tailed behavior, and the forecast errors and loss differentials derived from them can inherit this feature. We show that standard tools for forecast evaluation, including the widely used Diebold--Mariano test, can be severely unreliable in such settings: under empirically plausible heavy-tailed data-generating processes, a test at the 5\% nominal level can reject a true null of equal predictive ability as often as 70\% of the time. The root cause is that the standard normal critical values are justified by a finite-variance central limit theorem that breaks down when the loss differentials have infinite variance.
To address this, we develop a novel stable limit theorem for strongly mixing time series processes under infinite (or potentially undefined) variance, extending the classical central limit theorem to the heavy-tailed setting. Building on this result, we derive the limiting distribution of the test statistic, which is a ratio of dependent non-Gaussian random variables. The quantiles of this ratio can differ markedly from those of a standard normal distribution. We then consider a subsampling-based inferential procedure that is valid regardless of tail-heaviness and requires no prior knowledge of the tail index, the convergence rate of the sample mean, or estimation of long-run variances.
An empirical application to risk forecasts for emerging-market foreign exchange rates, where the loss differentials appear to have infinite variance, illustrates the practical relevance: the standard Diebold--Mariano test and the subsampling-based test yield substantively different conclusions about equal predictive ability. In particular, simple rolling-window risk forecasts are not found to be inferior to forecasts based on more sophisticated GARCH and score-driven models once the heavy-tailed nature of the loss differentials is properly accounted for.
\subsubsection*{Declaration of generative AI and AI-assisted technologies in the manuscript preparation process}
During the preparation of this manuscript, the authors used Opus 4.6 as a language-editing tool to improve clarity and readability, particularly in the Introduction. All content was subsequently reviewed and edited by the authors, who take full responsibility for the manuscript.
\newpage
\bibliographystyle{ecta}
\bibliography{Mybibliography}