Extracted main text — title through conclusion, appendix excluded. This is what our citation measures are computed over, published so the extraction can be checked by eye.
96,496 characters · 18 sections · 91 citation commands
Heavy Tails and Predictive Ability Testing
\noindentKeywords: Predictive ability, hypothesis testing, heavy tails, stable limit theory, subsampling.
\numberwithin{equation}{section}
The Diebold--Mariano (DM) test of diebold1995paring and the superior predictive ability (SPA) tests of white2000snooping and hansen2005SPA are workhorse tools for forecast comparison in economics and finance. Their asymptotic justification rests on the central limit theorem (CLT) of ibragimov1962clt for strongly mixing processes with 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., gabaix2009power and 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
against $\mathsf{H}_\mathrm{A}:\mathbb{E}[X_t]\neq 0$ with $\mathbb{E}[|X_t|]<\infty$. Depending on how $X_t$ is constructed, (ref) covers hypotheses such as forecast unbiasedness, forecast optimality, forecast encompassing, and, most prominently, the equal predictive ability (EPA) hypothesis of diebold1995paring\footnote{We focus on the unconditional predictive ability hypothesis; the conditional version of GiacominiWhite2006 is beyond the scope of our analysis. We also abstract from parameter estimation error in the forecasts; see west1996asymptotic for asymptotic theory under estimation uncertainty.}. Standard practice is to base inference on
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 ibragimov1962clt CLT provides the foundational justification for this practice, and underlies a large body of work on forecast evaluation, including white2000snooping, hansen2005SPA, hansen_lunde2005, GiacominiWhite2006, hansen_lunde_nason2011MCS, Patton_timmermann_2012rationality, BarendsePatton2022, and odendahl2023forecast.
This paper makes three main contributions. First, we develop a new stable limit theorem for infinite-variance strongly mixing time series. The theorem parallels the structure of the 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 ibragimov1962clt, that the existing forecast-evaluation literature requires, when loss differentials have infinite variance.
Second, we use this result to characterize the limiting distribution of the statistic $T_n$ in (ref) 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.
Third, we consider an inferential procedure that is valid under both finite- and infinite-variance regimes. We construct critical values by subsampling 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 white2000snooping and 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 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 (ref) 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 kim2021forecast, who also consider the DM test under heavy tails. kim2021forecast obtain a non-Gaussian stable limit for the DM statistic by direct application of Davis1983stablelimit and davis1995point and propose subsampling for inference, invoking the validity result of 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 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 kim2021forecast, McElroy_Politis_2002_subsampling, and 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) presents the stable limit theory and the resulting limiting distribution of $T_n$ (Proposition (ref)). Section (ref) presents the subsampling algorithm and establishes its asymptotic validity. Section (ref) contains the extensions to alternatives with undefined mean and to the SPA test. Sections (ref) and (ref) report the simulation experiments and the empirical illustration. The proof of the main Theorem (ref) in the Appendix; additional proofs and simulation results are in the Supplemental Appendix.
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, 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.
This section develops the limit theory underlying our inferential procedure. Section (ref) provides some motivating examples and previews the main result. Sections (ref)--(ref) state the partial-sum limit theory under finite variance and heavy tails, respectively. Section (ref) combines these into the limiting distribution of $T_n$ in (ref), and Section (ref) illustrates the theory in the context of linear autoregressive processes.
Given a random sample $(X_t)_{t=1,\dots,n}$ with a time-invariant mean, we consider testing the hypothesis $\mathsf{H}_0$ in (ref) 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 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.
In order to test $\mathsf{H}_0$, we consider the statistic $T_n$ in (ref) with $\hat{\sigma}_n^2=n^{-1}\sum_{t=1}^nX_t^2$. In contrast to 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
with
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$.
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,
The formal statement is given in Proposition (ref) in Section (ref), which is a direct consequence of the partial-sum limit theorems in Sections (ref)--(ref). Complementing this result, in Proposition (ref) 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.
We assume throughout that the data-generating process (DGP) for $X_t$ is 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
Following 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 (ref) under finite variance of $X_t$, we consider the following assumption.
Assumption (ref) is standard in the literature on forecast evaluations. The summability condition in (ref) 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) 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 (ref) to hold. Consequently, Assumption (ref) 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:
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
where $L$ is a slowly varying function and the tail-balance coefficients $p_+,p_-\geqslant0$ satisfy $p_++p_-=1$. Importantly, under (ref), $\mathbb{E}[{|X|}^{p}]=\infty$ $(<\infty)$ for every $p>\kappa$ $(p<\kappa)$, and, hence, a smaller $\kappa$ implies that $X$ has heavier tails.
Following 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 (ref) and there exists a real-valued stochastic process $(\mathit{\Theta}_t)_{t\in\mathbb{Z}}$ such that for all integers $h\geqslant0$,
where $\overset{w}{\to}$ denotes weak convergence of probability measures. The process $(\mathit{\Theta}_t)_{t\in\mathbb{Z}}$ is called the 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 mikosch2024extreme.
With this notion in hand, we impose the following 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).} Under Assumption (ref) there exists a deterministic sequence $(a_n)_{n\geqslant1}$ satisfying $a_n\to\infty$ and
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), we make the following assumption about the dependence structure of $(X_t)_{t\in\mathbb{Z}}$.
The condition in (ref) is a so-called anti-clustering condition; see e.g. bartkiewicz2011stable. To provide some intuition for this condition, note that (ref) implies
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) play complementary roles in the limit theory, and it is useful to separate them conceptually. Regular variation (Assumption (ref)) together with the anti-clustering condition (ref) 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 (ref) 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) for additional details about this approximation argument.
In addition to playing complementary roles, note that the anti-clustering condition (ref) and the mixing condition (ref) depend on the same deterministic sequence $(r_n)_{n\geqslant1}$. Similar to the role of the constant $\varepsilon>0$ in Assumption (ref), 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 (ref) to hold.
Assumption (ref) 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)), GARCH processes, and other processes given by affine stochastic recurrence equations, as well as stochastic volatility processes; see Section (ref) of the Supplemental Appendix for additional details and technical conditions.
Assumptions (ref)--(ref) ensure that $(X_t)_{t\in\mathbb{Z}}$ obeys the following stable limit 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., davis1995point), our result builds upon recent theoretical results of 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
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 matsui2025self to strongly mixing processes and allow for $\mu\neq 0$. For further details about hybrid characteristic function-Laplace transforms, we refer to mikosch2024extreme, and we refer to the textbook by samorodnitsky1994stable for a treatment of non-Gaussian stable distributions.
The distribution of the limiting vector in (ref) 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), $(X_t)_{t\in\mathbb{Z}}$ is regularly varying and hence, by the definition in (ref), has a spectral tail process $(\mathit{\Theta_t})_{t\in\mathbb{Z}}$. Define the extremal cluster process as
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
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$,
with scale and skewness parameters given, respectively, by
where
and $\Gamma(\cdot)$ the gamma function\footnote{A slightly different expression for the characteristic function for $\xi_\kappa$ can be found in mikosch2024extreme. For the sake of completeness, we provide detailed derivations of (ref) and (ref) in Section (ref) 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 (ref). For certain processes, such as linear autoregressive processes considered in Section (ref) 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.
Combining Theorems (ref) and (ref) yields the limiting behavior of the test statistic $T_n$ in (ref) under both finite and infinite variance:
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)), the limit in case (2) is asymmetric whenever the skewness parameter $\beta$ in (ref) is non-zero, whereas the Gaussian limit in case (1) is always symmetric. This is investigated in more detail in Section (ref) 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 matsui2025moments; nevertheless, as documented in the simulation study in Section (ref), 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.
We end this section by considering additional details about the results in Theorem (ref) and Proposition (ref) 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). Let the real-valued $X_t$ obey the recursion
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 (ref), which is strongly mixing with geometric rate; see, e.g., Proposition 2.2.4 in buraczewski2016stochastic. Consequently, if $\mathbb{E}[{|Z_t|}^{2+\varepsilon}]<\infty$ for some $\varepsilon>0$, then $(X_t)_{t\in\mathbb{Z}}$ satisfies Assumption (ref). 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$:
By Proposition (ref), $(X_t)_{t\in\mathbb{Z}}$ satisfies Assumptions (ref)--(ref) whenever $Z_t$ is regularly varying with tail index $\kappa\in(0,1)\cup(1,2)$, and we have the following result.
Proposition (ref) 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), we consider the properties of the standard diebold1995paring test when $X_t$ is heavy tailed. The Diebold--Mariano test is based on the test statistic in (ref) 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 newey_west_1987, and for the sake of simplicity we here let
for a fixed integer $q>0$ and fixed weights $w_j>0$, $j=1,\dots,q$. We have the following result.
In this section we consider subsampling-based techniques for testing the hypothesis $\mathsf{H}_0$ in (ref) 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., McElroy_Politis_2002_subsampling, kokoszka2004subsampling and 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 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.
The algorithm is tailored to work irrespective of $X_t$ having finite variance (Assumption (ref)) or being regularly varying with $\kappa\in(1,2)$ (Assumptions (ref)--(ref)). In particular, the algorithm does not rely on estimation of the long-run variance $\sigma^2$ in Theorem (ref), and it does not require any knowledge about (or estimation of) the scaling sequence $(a_n)_{n\geqslant1}$ in (ref).
The asymptotic validity of the subsampling-based test in Algorithm (ref) is ensured under the strong mixing conditions in Assumptions (ref) and (ref). In particular, we have the following result.
Theorem (ref) states that the subsampling-based test in Algorithm (ref) is asymptotically valid for any block size $b_n=o(n)$ with $b_n\to\infty$. As pointed out by, e.g., 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 politis2012subsampling. The so-called Minimum Volatility Method politis2012subsampling avoids this assumption, but requires a user-chosen grid of candidate values. kim2021forecast propose a rule depending on estimated tail-index and skewness parameters in the iid stable-distribution setting. We adopt the simple rule
similar to block sizes used in bai2016unified and zhang2013blocksampling, for simplicity and to maintain consistency with existing work on subsampling for heavy-tailed dependent data. The Monte Carlo experiments in Section (ref) indicate that this rule yields reasonable finite-sample rejection frequencies across a range of tail indices and sample sizes.
In this section we consider extensions of the limit theory and related subsampling-based inference considered in Sections (ref)--(ref). In Section (ref) 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) we consider testing a hypothesis about a mean of a vector with particular attention to tests for superior predictive ability (SPA).
A straightforward implication of Theorem (ref) is that when the tail index\footnote{The case $\kappa\in(0,1)$ may indeed be empirically relevant: 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: 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)$,
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), we have that
see Proposition (ref) 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
Under Assumption (ref) or Assumptions (ref)--(ref), 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) we show that, with $a_n$ given by (ref), $\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
We have the following result.
Critical values for a test based on $\widetilde{T}_n$ can be computed by the subsampling Algorithm (ref) in the Supplemental Appendix. Importantly, under $\kappa\in(0,1)$ we show in Section (ref) that the subsample statistics $\widetilde{T}_{i,{b_n}}$ satisfy
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) 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) achieves non-trivial rejection frequencies for sufficiently large sample sizes.
In Section (ref) of the Supplemental Appendix we provide a multivariate version of Theorem (ref). The result, Theorem (ref), 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 hansen2005SPA, see also 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
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 superior predictive ability (SPA) over the other $m$ competing forecast series if
is true, where the inequality holds entry-wise. hansen2005SPA proposes a test for $\mathsf{H}_{\mathrm{SPA}}$ against
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) in Section (ref) 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) in the Supplemental Appendix), we consider a modification of the test statistic proposed by hansen2005SPA. Define
and the test statistic
We have the following result.
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).
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) with block length $b_n$ given by (ref). We refer to Section (ref) 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).
Throughout, we consider a DGP for $X_t$ given in terms of the AR(1) process
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 Stable$(\kappa,0,1,0)$ and 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 (ref) 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) $(X_t)_{t\in\mathbb{Z}}$ is regularly varying with index $\kappa$. Moreover, by Proposition (ref), 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
with $S_n$ given by (ref) and with $\hat{\sigma}^2_{n,\mathrm{NW}}$ the standard newey_west_1987 HAC variance estimator with lag length given by $\lfloor4(n/100)^{2/9}\rfloor$, as considered in 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 (ref). Consequently, from Proposition (ref) 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$.
Table (ref) 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) 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.
Figure (ref) contains rejection frequencies for the subsampling-based test when the null hypothesis is violated. We consider DGPs given by (ref) with $Z_t$ 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). 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.
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 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), 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) 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.
Quantification and hedging of foreign currency risk are important for international portfolio management\footnote{See, e.g., christensen_varneskov_2021_FX and the references therein.}, and we consider in this section an application of Algorithm (ref) to comparing risk forecasts for foreign exchange (FX) rate returns. FX rates are typically volatile, and their returns are heavy-tailed. For instance, 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}$:
In the context of financial risk management, this quantile is the Value-at-Risk (VaR) at risk level $\tau$. Following Example (ref) one may evaluate the time $t$ VaR forecast using the tick loss function,
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
For standard GARCH-type processes, $\sigma_t^2$ can be written as a stochastic recurrence equation (see, e.g., Section (ref) 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 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 patton2019_es_var. The three methods based on score-driven models are labeled "FZ-2F", "FZ-1F", and "Hybrid" and are described in patton2019_es_var.\footnote{Along the lines of 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 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
Following Example (ref), 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) presents Diebold--Mariano $t$-statistics on the loss differentials\footnote{As in 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) with block size given by (ref). 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), 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 davis2018tail_process_inference the tail-balance coefficient $p_+$ in (ref) 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) 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.
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.
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.