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.
101,238 characters
Illiquidity at Risk
\title{Illiquidity at Risk\thanks{
We thank the participants at the 2nd Italian Conference on Economic Statistics (University of Florence, 2024), the ICEEE conference (University of Palermo, 2025), the SoFiE conference (ESSEC Business School, 2025), the CREST seminar (ENSAE Paris, 2026), and the QFFE conference (Marseille, 2026). Paolo Santucci de Magistris also acknowledges the research support of the European Union’s Next Generation EU program through the Italian PRIN 2022 - M4C2, Investment 1.1 - “Monitoring Risks in Financial Markets” (Codice Cineca: 2022NEL482 - CUP: I53D23003410008).}}
\author{Demetrio Lacava\thanks{University of Messina, Department of Economics, Piazza Pugliatti, 1, 98122 Messina, Italy. E-mail: [email removed]}
\and
Paolo Santucci de Magistris \thanks{Luiss University, Department of Economics and Financial Markets, Viale Romania 32, 00197
Roma, Italy. E-mail: [email removed].}\hspace{1.5mm} \thanks{Corresponding author.} }
\date{}
\maketitle
\begin{abstract}
Market efficiency relies fundamentally on stable liquidity. Consequently, forecasting liquidity dynamics is a priority for both investors and regulators. We introduce a new tail-risk metric, Illiquidity-at-Risk (IlliQaR), designed to quantify the magnitude of extreme liquidity dry-ups. Relying upon the realized Amihud (a precise illiquidity measurement derived from high-frequency data as the ratio of realized volatility to trading volume) we assess the predictive power of various linear and non-linear econometric models, with a specific focus on the impact of discontinuous jump components. Accounting for these jumps is essential for achieving accurate probability coverage and better IlliQaR predictions during periods of systemic stress, where standard continuous models systematically underestimate the severity of liquidity evaporation. Our empirical analysis, encompassing the S\&P 500 index and a cross-section of 25 large U.S. equities, demonstrates that incorporating jumps significantly improves forecasts of illiquidity. Our results suggest that individual stock IlliQaR violations often cluster during periods of S\&P 500 liquidity stress. This indicates that \textit{Illiquidity at Risk} is not just a localized concern but a systemic one, where the main index acts as a leading indicator for extreme dry-ups in individual stock liquidity.
\end{abstract}
\noindent \textbf{Keywords:} Liquidity, Volume, Realized Amihud, Jumps, Forecast, VaR.\\
\vspace{0.3cm}
\emph{J.E.L. classification}: C15, F31, G12, G15
\spacing{1.55}
\newpage
\section{Introduction}
Financial markets rely on two pillars for efficient functioning: price discovery and liquidity. While volatility has long been the primary metric for assessing risk and price discovery, liquidity, i.e. the ability to trade substantial quantities of an asset quickly and at low cost, is equally critical. In many contexts, liquidity is the invisible infrastructure of the market; when it functions well, it is taken for granted, but its disappearance can trigger systemic crises. Conceptually, liquidity risk differs from price risk. While price risk (or volatility) refers to the uncertainty regarding the future value of an asset, liquidity risk refers to the uncertainty regarding the ability to realize that value. It is the risk that an investor cannot convert a financial asset into cash (or vice versa) without incurring a prohibitive cost or significantly moving the price. Conceptually, liquidity is an intrinsic and unobservable market characteristic; its sudden disappearance can transform idiosyncratic shocks into systemic crises, establishing market liquidity as a fundamental risk factor in asset pricing, see among others \cite{pastor2003liquidity} and \cite{acharya2005asset}.
While the existing literature has focused extensively on modeling the level of illiquidity, market participants are often more concerned with the tails, i.e. the sudden ``evaporation'' of liquidity that characterizes flash crashes and financial panics. To address the challenge of quantifying this dimension of liquidity risk, this paper introduces a novel contribution: Illiquidity-at-Risk (IlliQaR). Much like Value-at-Risk (VaR) in market risk management, IlliQaR is designed to capture the tail risk of market liquidity, answering the question: 'What is the maximum level of illiquidity expected over a given horizon at a specific confidence level?'. The reliability of IlliQaR as a risk metric rests on two critical pillars. First, the use of high-frequency data to construct a highly precise measurement of illiquidity allows for the filtration of observation noise, which is shown to mask the true tail behavior when employing lower-frequency proxies. Second, empirical evidence demonstrates that a robust IlliQaR framework must explicitly account for discontinuous illiquidity bursts. As these bursts are the primary drivers of severe liquidity shocks, failure to model them results in a significant underestimation of risk during periods of systemic stress.
We construct IlliQaR by modeling the conditional density of illiquidity, paying particular attention to the role of discontinuous jumps in liquidity dynamics. We argue that failing to account for these jumps leads to a severe underestimation of liquidity risk, leaving portfolios exposed to unexpected transaction costs exactly when trading is most urgent. The inclusion of a jump component is not merely an econometric convenience, but is heavily rooted in market microstructure theory. Sudden evaporation of liquidity—or illiquidity jumps—can be driven by the sudden binding of dealers' capital constraints, see \cite{brunnermeier2009market}, the synchronized withdrawal of algorithmic `phantom liquidity' during periods of stress or flash crashes (\citealp{kirilenko2017flash}), and discrete shocks to adverse selection risk following information arrivals, as in the seminal work of \cite{glosten1985bid}.
The necessity of such a precise risk metric brings us to the fundamental problem of measurement. Before one can model the risk of illiquidity, one must define and measure the underlying latent variable. The academic literature proposes a taxonomy of illiquidity based on three main dimensions: tightness, depth, and resiliency. Tightness refers to the cost of a round-trip transaction and is typically approximated by the bid-ask spread. Following the pioneering work of \cite{Roll1}, numerous studies including \cite{hasbrouck2009trading}, \cite{CorwinSchultz}, and \cite{Abdi2017} have focused on estimating the effective spread from low-frequency data. While bid-ask tightness provides a snapshot of transaction costs, it represents only a single dimension of market liquidity. A market can exhibit narrow spreads yet lack the necessary depth to absorb significant order flow without precipitating adverse price movements. Depth is visible in the limit order book as the volume available at the best bid and ask, but it is harder to reconstruct historically. Finally, and perhaps most importantly for stability, is resiliency: the speed at which prices recover to their equilibrium level after a large trade causes a dislocation.
Measuring resiliency and the price impact of trades requires a more sophisticated approach than simply observing spreads. \cite{kyle1985continuous} famously introduced ``lambda'' to proxy for the elasticity of prices, but it was the contribution of \cite{Amihud2002illiquidity} that provided the most enduring and practical measure: the ratio of absolute returns to trading volume. While the Amihud measure provides an intuitive proxy for price impact by scaling absolute returns by trading volume, daily proxies often introduce substantial noise, necessitating the use of higher-frequency realized measures for greater accuracy. In this paper, we emphasize the importance of high-precision measurement to disentangle genuine liquidity information from measurement error. We therefore rely on the \textit{realized Amihud} measure first proposed by \cite{RanaldoSantucci2020} and refined in \cite{lacava2023realized}. By aggregating high-frequency intraday data, this measure provides a ``noise-reduced" signal of daily illiquidity. Precision here is not merely a statistical luxury; it is a prerequisite for risk management. If the measurement of illiquidity is contaminated by noise, any risk metric derived from it (like IlliQaR) will yield false alarms or, worse, fail to detect stress scenarios.
Building on this high-precision proxy, we evaluate a number of econometric models to forecast illiquidity and estimate IlliQaR. We employ both linear and non-linear specifications, in particular the Multiplicative Error Model (MEM) framework of \cite{Engle:2002}, which is uniquely suited for non-negative processes like illiquidity. To account for the well-known persistence of liquidity (the fact that illiquid days tend to cluster together), we incorporate the Heterogeneous Autoregressive (HAR) structure of \cite{Corsi2009}. Crucially, we augment these models with a ``jump" component, utilizing the MEM-J specification of \cite{Caporin:Rossi:DeMagistris:2017}. This is motivated by the empirical observation that liquidity does not always deteriorate gradually; it often ``jumps" abruptly in response to news or market shocks. We investigate whether these jumps are merely transient noise or if they carry information about future risks. By explicitly modeling the probability of these jumps, we aim to correct the distributional assumptions that underpin standard risk models, particularly addressing the positive skewness and fat tails characteristic of liquidity distributions.
Our empirical analysis yields several robust findings that extend our understanding of illiquidity dynamics. We analyze the daily realized Amihud series for the S\&P 500 and a cross-section of 25 major individual U.S. equities. While the S\&P 500 index serves as a barometer for systemic liquidity risk, the transmission of these dry-ups to individual constituents is non-uniform. By analyzing 25 large-cap equities, we bridge the gap between market-wide IlliQaR and idiosyncratic liquidity failures, examining whether the jump dynamics observed at the index level are a fundamental property of the underlying assets or a result of aggregation.
First, consistent with \cite{lacava2023realized}, we confirm that illiquidity is a long-memory process characterized by strong clustering. This aligns with the findings of \cite{hafner2023dynamic,hafner2025permanent}, who attribute the observed persistence of illiquidity (as measured by the range Amihud proxy introduced by \citealp{lacava2023realized}) to a slowly varying trend component that is potentially subject to structural breaks. Furthermore, we demonstrate that illiquidity dynamics are also characterized by significant jump activity.
Second, and perhaps most consequentially, we find that identifying and modeling these jumps is essential for short-term forecasting. In the context of 1-step-ahead forecasts, the MEM-J specification significantly outperforms models that assume purely continuous dynamics. This superiority is driven by the model's ability to react rapidly to the onset of liquidity stress. Conversely, over medium-to-long horizons, the importance of the jump component fades as illiquidity mean-reverts; in these cases, simpler asymmetric specifications perform adequately.
Third, regarding the estimation of IlliQaR, we find that standard models fail to provide adequate coverage for extreme events, tending to underestimate the frequency of severe illiquidity spikes. Our proposed models, which account for both the discontinuous illiquidity jump component and the specific error distribution, successfully correct this bias and provide accurate probability coverage in the tails. Notably, further improvements in IlliQaR accuracy are achieved by disentangling illiquidity jumps from price jumps—specifically, by utilizing bipower variation rather than realized volatility. Taken together, these results confirm that IlliQaR is a viable and indispensable tool for modern risk management.
Finally, we investigate the economic determinants of IlliQaR violations to understand why liquidity dries up. By regressing the probability of an extreme illiquidity event against various market variables, we find a strong asymmetric effect, analogous to the leverage effect in volatility: liquidity is significantly more likely to evaporate during market downturns than during rallies. Furthermore, we document a strong positive relationship between IlliQaR violations and the VIX index, as well as measures of economic policy uncertainty. Interestingly, risk aversion plays a distinct role; as investors become more risk-averse, the probability of a liquidity crisis increases, likely due to a withdrawal of market-making capital. These relationships hold both in-sample and in out-of-sample forecasting tests, validating the economic rationale behind our risk metric.
The remainder of the paper is organized as follows. Section \ref{sec:IlliQaR} formally defines the concept of Illiquidity-at-Risk, introduces the realized Amihud and the baseline econometric frameworks, specifically the Multiplicative Error Model (MEM) and the Heterogeneous Autoregressive (HAR) specifications. Section \ref{sec:jump_model} extends this framework by incorporating a jump component (MEM-J) and outlines the identification mechanism for conditional jumps. Section \ref{sec:empirical_analysis} is devoted to the empirical application: Section \ref{sec:dataset} provides descriptive statistics for the S\&P 500 index; Section \ref{sec:results} discusses the in-sample estimation results and the significance of the jump parameters; and Section \ref{sec:forecast} evaluates the out-of-sample forecasting performance. Section \ref{sec:illiquidity_determinats} focuses on the evaluation of the IlliQaR metric via coverage tests and analyzes the economic drivers of liquidity stress. Section \ref{sec:individual} extends the analysis by assessing IlliQaR on individual US stocks. Finally, Section \ref{sec:conclusion} concludes. Additional empirical findings and robustness checks are provided in the Supplementary Document.
\section{Illiquidity-at-Risk}\label{sec:IlliQaR}
We introduce the concept of Illiquidity-at-Risk. Analogous to the VaR for large negative returns, IlliQaR is defined as the threshold level, $\text{IlliQaR}_t(p)$, that the realized illiquidity measure, $\text{Illiq}_t$, is expected to exceed with probability $p$. Formally
\begin{equation}\label{eq:IlliQaR}
\text{Pr}(\text{Illiq}_t>\text{IlliQaR}_t(p)|\mathcal{F}_{t-1})=p, \quad t=1,2,\ldots,T,
\end{equation}
where $\mathcal{F}_{t-1}$ is the available information set. The evaluation of the probability in \eqref{eq:IlliQaR} requires several ingredients. In particular, IlliQaR relies on three methodological pillars: the derivation of a precise illiquidity measurement to proxy the unobservable true process, the implementation of a dynamic predictive model to account for time-series dependencies, and the use of a flexible econometric distribution to accurately model tail illiquidity events.
\subsection{Measuring Illiquidity}\label{sec:meas}
As a measurement of illiquidity we consider the \textit{realized Amihud}, introduced by \cite{RanaldoSantucci2020}, which refines the illiquidity proxy proposed by \cite{Amihud2002illiquidity}. Realized illiquidity of a given security is defined as the ratio of two observable quantites referring respectively to realized volatility and volume computed on a given unit interval (e.g. a day, a week, a month). In particular, assuming to split a unit period of time into $M$ subperiods, the realized Amihud on a given interval of unit-length is defined as
\begin{equation}
\text{Illiq}=\frac{RPV}{\nu} \label{eq:realized_Amihud},
\end{equation}
where $RPV=\sum_{i=1}^{M}|r_i|$ denotes the realized power variation \citep{barndorff2003realized} of order one (or the realized absolute variation) with $r_i$ being the log-return on the $i$-th subperiod. Similarly, $\nu=\sum_{i=1}^{M}\nu_i$ denotes the trading volume, which is obtained as the sum of the volume generated on each sub-interval interval $i=1,\ldots M$. \cite{lacava2023realized} provide a comprehensive theoretical derivation of the properties of the realized Amihud in measuring \textit{stochastic illiquidity} by assuming that the market of a given security is populated by traders with different reservation prices. In particular, \cite{lacava2023realized} show that $\text{Illiq}$ converges to the true illiquidity signal aggregated on a given interval of unit length as $M \to \infty$. Following the framework outlined in \cite{RanaldoSantucci2020}, the theoretical quantity measured by the realized Amihud is the \textit{integrated illiquidity}, defined as $\text{Illiq}^* = \frac{1}{\int_{0}^{1} \ell(s) ds}$, where $\ell(s)$ measures the instantaneous sensitivity of the number of trades to new information, representing a dimension of illiquidity related to market depth.
\cite{lacava2023realized} show that the precision of the realized Amihud in measuring daily integrated illiquidity based on intradaily returns sampled at 5-minutes intervals is several times more efficient than that obtained with daily Amihud. In particular, the daily Amihud of \cite{Amihud2002illiquidity} is obtained as a special case of realized Amihud, and it is obtained when $M=1$, i.e. using one observation on returns per day.
\subsection{Linear and non-linear models for illiquidity }\label{sec:lin-nonlin}
We consider several econometric specifications designed to model the dynamic evolution of illiquidity and of its quantiles, with the goal of providing a precise prediction of IlliQaR. By computing an illiquidity measure over disjoint daily intervals, we obtain a time series of illiquidity measurements, $\{\text{Illiq}_t\}_{t=1}^{T}$, where $T$ denotes the sample size. Given that illiquidity is a non-negative, persistent, and heavy-tailed process, its dynamic and distributional features can be effectively analyzed using the econometric framework adopted in the realized variance literature, following the seminal contribution of \cite{andersen2003modeling}.
First, we consider a linear benchmark, the Heterogeneous Autoregressive (HAR) model of \cite{Corsi2009}. The model is defined as
\begin{equation}
\text{Illiq}_t=\omega + \alpha_1 \text{Illiq}_{t-1} +\alpha_2 \text{Illiq}_{t-1:t-5} + \alpha_3 \text{Illiq}_{t-1:t-21} + u_t, \quad u_t|\mathcal{F}_{t-1} \sim N(0,\sigma^2_u),
\label{eq:har_lin}
\end{equation}
where $\text{Illiq}_{t-1:t-5}$ ($\text{Illiq}_{t-1:t-21}$) denote the weekly (monthly) realized Amihud, computed as the rolling average over the last 5 (21) days, and $\mathcal{F}_{t-1}$ is the information set available at time $t-1$. This specification is particularly well-suited for modeling the long-range dependence structure of illiquidity, through a simple linear additive structure of daily, weekly, and monthly components. Recent studies by \cite{hafner2023dynamic, hafner2025permanent} highlight the importance of disentangling slow-moving trends from fast-moving innovations in illiquidity. The HAR specification can be viewed as a reduced-form of this component-based approach; it approximates the multi-scale dynamics of illiquidity through a cascade of autoregressive terms that represent varying speeds of adjustment. It follows that for the Gaussian HAR model, IlliQaR is
\[
\text{IlliQaR}_t^{\text{HAR}}(p)=\omega + \alpha_1 \text{Illiq}_{t-1} +\alpha_2 \text{Illiq}_{t-1:t-5} + \alpha_3 \text{Illiq}_{t-1:t-21}+\Phi^{-1}(1-p)\sigma_u,
\]
where $\Phi^{-1}(1-p)$ denotes the inverse cumulative distribution function (CDF) of the standard Gaussian random variable evaluated at $1-p$ (as we consider exceedance over the right tail).
We also extend the analysis to a class of non-linear specifications specifically designed to accommodate non-negative stochastic processes, namely the class of Multiplicative Error models (MEM), as introduced by \cite{Engle:2002} and revised by \cite{Engle:Gallo:2006} to account for an asymmetric response to the sign of past news (i.e. the Asymmetric MEM - AMEM). The AMEM for $\text{Illiq}_t$ is specified as
\begin{equation}
\begin{array}{l}
\text{Illiq}_{t}=\mu _{t}\epsilon _{t }, ~ \qquad \epsilon _{t}|\mathcal{F}_{t-1}\sim \Gamma(\vartheta ,\frac{1}{\vartheta }),\\
\mu _{t}=\omega +\alpha_1 \text{Illiq}_{t-1}+\beta_1 \mu _{t-1}+\gamma
D_{t-1}\text{Illiq}_{t-1}, \quad t=1,2,\ldots,T,
\end{array}\label{eq:amem}
\end{equation}
where $D_{t-1}$ a dummy variable taking value of 1 if the return of the considered asset is negative, 0 otherwise. Since $\mu_t$ follows the GARCH(1,1) dynamic process of \citep{Bollerslev:1986}, the usual parameter constraints for positiveness ($\omega,\alpha_1,\beta_1,\gamma>0$) and stationarity ($\alpha_1 + \beta_1 +\gamma/2>0$) are imposed. The conditional density of the error term is
\begin{equation*}
\begin{array}{l}
f(\epsilon_t|\mathcal{F}_{t-1})=\frac{1}{\Gamma(\vartheta)}\vartheta^\vartheta\epsilon_t^{\vartheta-1}e^{-\vartheta\epsilon_t},\qquad \epsilon_t>0,
\end{array} \label{eq:density}
\end{equation*}
which is the one of a Gamma distribution with shape parameter $\vartheta>0$ and scale $\frac{1}{\vartheta}$. This means that $\epsilon_t$ is a random variable with unit mean and variance equal to $\frac{1}{\vartheta}$, denoted as $\epsilon_t\sim \Gamma(1,\vartheta)$ in the mean-shape representation. It follows that $\mbox{E}(\text{Illiq}_t|\mathcal{F}_{t-1})=\mu_t$ and $\mbox{Var}(\text{Illiq}_t|\mathcal{F}_{t-1})=\mu^2_t/\vartheta$.
The conditional density of $\text{Illiq}_{t}$ is therefore available in closed form, thus allowing us to estimate the model parameters by maximum likelihood (ML). The asymptotic properties of the ML estimator of the unknown coefficients are discussed in \cite{Engle:2002} and \cite{Engle:Gallo:2006}. By recalling the quasi-ML principle, they show that the ML estimator is consistent and efficient, irrespectively of the appropriateness of the distribution for the error term. We also consider the HAR specification for $\mu_t$ (see \citealp{Gallo:Otranto:2015} and \citealp{Caporin:Rossi:DeMagistris:2017}), that is
\begin{equation}
\mu_t=\omega + \alpha_1 \text{Illiq}_{t-1} +\alpha_2 \text{Illiq}_{t-1:t-5} + \alpha_3 \text{Illiq}_{t-1:t-21} +\beta \mu_{t-1} + \gamma D_{t-1} \text{Illiq}_{t-1}.
\label{eq:har}
\end{equation}
We call this model G-AMEM-HAR. The G-AMEM-HAR nests all the other specifications for $\mu_t$. For instance, the AMEM is obtained by imposing $\alpha_2=\alpha_3=0$, while with $\beta=0$ the model reduces to the AMEM-HAR. Under the various MEM specifications, the IlliQaR is
\begin{equation}\label{eq:illiqar_mem}
\text{IlliQaR}^\text{MEM}_t(p)= \frac{\mu_t}{\vartheta}\gamma^{-1}(\vartheta,\Gamma(\vartheta)(1-p)),
\end{equation}
where $\gamma^{-1}(\cdot)$ is the inverse of the lower incomplete Gamma function.
\subsection{Illiquidity jumps}\label{sec:jump_model}
Although the adoption of the Gamma distribution as in the baseline MEM model is justified in view of its robustness \citep[see][]{Cipollini:Engle:Gallo:2013}, other distributions are widely considered in the current literature: \cite{Engle:2002} suggests to use the exponential distribution, while mixture models are provided by \cite{Lanne:2006}. However, as shown by \cite{Caporin:Rossi:DeMagistris:2017}, often these distributional assumptions do not ensure an adequate coverage of the probability of tail events of volatility. We therefore consider a variant of the MEM specification that features conditional jumps dynamics in illiquidity, as originally proposed by \cite{Caporin:Rossi:DeMagistris:2017}. In such a model, $\text{Illiq}_t$ is the result of the product of three elements
\begin{equation}
\begin{array}{l}
\text{Illiq}_t=\mu_tZ_t\epsilon_t,
\end{array}\label{eq:mem-j}
\end{equation}
where $\mu_t$ and $\epsilon_t$ are defined as in \eqref{eq:amem}, while $Z_t$ is the jump component, assumed to be independent of $\epsilon_t$. The jump term is modeled as a compound Poisson random variable, where the expected number of jump arrivals at time $t$ ($N_t$) is governed by a Poisson random variable with a time-varying intensity $\kappa_t$: when $N_t>0$, the jump size is given by the sum of independent Gamma random variables, $Y_t|\mathcal{F}_{t-1}\sim \Gamma(d_{\kappa_t},\zeta)$ in mean-shape representation, corresponding to a scale given by $\frac{d_{\kappa_t}}{\zeta}$. In particular, $Z_t=d_{\kappa_t}$ if $N_t=0$ and $Z_t=\sum_{j=1}^{N_t}Y_{j,t}$ for $N_t>0$, where $d_{\kappa_t}=(e^{-\kappa_t}+\kappa_t)^{-1}$ ensures $E[Z_t\epsilon_t|\mathcal{F}_{t-1}]=1$.
As a consequence of Poisson distributed jumps, the probability of observing a $j\geq 0$ jumps, conditioning on $\mathcal{F}_{t-1}$, is
\begin{equation*}
\text{Pr}(N_t=j|\mathcal{F}_{t-1})=\frac{e^{-\kappa_t}\kappa^j_t}{j!},\qquad j=0,1,2,\dots.
\label{eq:poisson}
\end{equation*}
Similarly to \cite{Chan:Maheu:2002}, we allow for a time-varying dynamics of the jump size, by specifying an autoregressive structure for $\kappa_t$ as
\begin{equation*}
\kappa_t=\phi_1+\phi_2\kappa_{t-1}+\phi_3\xi_{t-1},
\label{eq:jump_intensity}
\end{equation*}
where the jumps innovations are defined as
\begin{equation*}
\xi_t=E(N_t|\mathcal{F}_t)-E(N_t|\mathcal{F}_{t-1})=\sum_{j=0}^{\infty}j \text{Pr}(N_t=j|\mathcal{F}_t)-\kappa_t.
\end{equation*}
Imposing $\phi_1>0$ and $1>\phi_2>\phi_3>0$ is sufficient to ensure the positiveness and stationarity of $\kappa_t$, see \cite{Chan:Maheu:2002}. In particular, \cite{Maheu:McCurdy:2004} and \cite{maheu2013jumps} successfully applied this specification to study conditional jumps in stock returns, while \cite{Caporin:Rossi:DeMagistris:2016} focused the analysis on realized variance. Furthermore, acknowledging that $N_t$ cannot be directly observed in $\mathcal{F}_t$, $\xi_t$ can be interpreted as the forecast of $N_t$ as new information becomes available. Finally, the filtered probability needed to compute $E(N_t|\mathcal{F}_t)$ is obtained via the Bayes rule as
\begin{equation}
\text{Pr}(N_t=j|\mathcal{F}_t)=\frac{f(\text{Illiq}_t|N_t=j,\mathcal{F}_{t-1})\text{Pr}(N_t=j|\mathcal{F}_{t-1})}{f(\text{Illiq}_t|\mathcal{F}_{t-1})},
\label{eq:bayes_rule}
\end{equation}
where the density of $\text{Illiq}_t$ conditional on $N_t=j$ and $\mathcal{F}_{t-1}$ is
\begin{equation*}
\begin{large}
f_{\text{MEM-J}}(\text{Illiq}_t|N_t=j,\mathcal{F}_{t-1})=
\begin{cases}
\frac{1}{\text{Illiq}_t}\left(\frac{\vartheta \text{Illiq}_t}{d_{\kappa_t} \mu_t}\right)^{\vartheta}~ \frac{e^{\left(\frac{-\vartheta \text{Illiq}_t}{d_{\kappa_t} \mu_t}\right)}}{\Gamma(\vartheta)}, \qquad \qquad \qquad \qquad \qquad \quad ~ N_t=0\\
\frac{2}{\text{Illiq}_t}\left(\frac{\text{Illiq}_t}{\mu_t}\frac{\vartheta \zeta}{d_{\kappa_t}}\right)\left(\frac{j\zeta+\vartheta}{2}\right)\frac{1}{\Gamma(j\zeta)\Gamma(\vartheta)}\mathds{K}_{j\zeta-\vartheta}\left(2\sqrt{\frac{\text{Illiq}_t}{\mu_t}\frac{\vartheta \zeta}{d_{\kappa_t}}}\right),\quad N_t=j>0,
\end{cases}
\end{large}
\end{equation*}
where $\mathds{K}(\cdot)$ is the modified Bessel function of second kind, while the denominator in Eq. \eqref{eq:bayes_rule} is given by
\begin{equation*}
f_{\text{MEM-J}}(\text{Illiq}_t|\mathcal{F}_{t-1})=e^{-\kappa_t}f_{\text{MEM-J}}(\text{Illiq}_t|N_t=0,\mathcal{F}_{t-1})+\sum_{j=1}^{\infty}\frac{e^{-\kappa_t} \kappa_t^j}{j!} f_{\text{MEM-J}}(\text{Illiq}_t|N_t=j,\mathcal{F}_{t-1}).
\end{equation*}
The conditional density, which is therefore a mixture of Gamma and Kappa distributions, is available in closed form, allowing for model parameters estimation by ML. The following proposition derives the analytical expression for the IlliQaR under the MEM-J specification.
\begin{proposition}\label{theo1}
The IlliQaR for the MEM-J model is obtained as the value $\text{IlliQaR}_t(p)$ such that
\begin{equation}
\text{IlliQaR}^\text{MEM-J}_t(p): \qquad p=1-F_{\text{MEM-J}}(\text{IlliQaR}_t^\text{MEM-J}(p)|\mathcal{F}_{t-1}),
\end{equation}
where $F_{\text{MEM-J}}(\cdot|\mathcal{F}_{t-1})$
denotes the conditional CDF of the MEM-J, that is
\begin{equation*}
F_{\text{MEM-J}}(\text{Illiq}_t|\mathcal{F}_{t-1})=e^{-\kappa_t}F_{\text{MEM-J}}(\text{Illiq}_t|N_t=0,\mathcal{F}_{t-1})+\sum_{j=1}^{\infty}\frac{e^{-\kappa_t} \kappa_t^j}{j!} F_{\text{MEM-J}}(\text{Illiq}_t|N_t=j,\mathcal{F}_{t-1}),
\end{equation*}
where $F_{\text{MEM-J}}(\text{Illiq}_t|N_t=0,\mathcal{F}_{t-1})=\frac{1}{\Gamma(\vartheta)}\int_{0}^{\frac{\vartheta \text{Illiq}_t}{d_{\kappa_t}\mu_t}}s^{\vartheta-1}e^{-s}ds$
is the CDF of a Gamma-distributed random variable with shape $\vartheta$ and scale $\frac{d_{\kappa_t} \mu_t}{\vartheta}$ and
\[
F_{\text{MEM-J}}(\text{Illiq}_t|N_t=j,\mathcal{F}_{t-1})=\frac{2^{2-j\zeta-\vartheta}}{\Gamma(j\zeta)\Gamma(\vartheta)}\int_{0}^{2\sqrt{\frac{\zeta\vartheta\text{Illiq}_t}{d_{\kappa_t}\mu_t}}}s^{j\zeta+\vartheta-1}\mathds{K}_{|j\zeta-\vartheta|}\left(s\right)ds,
\]
is the CDF of a Kappa distributed random variable.
\end{proposition}
The CDF of the MEM-J is an essential element for the calculation of IlliQaR as well as for testing the adequacy of this specification as illustrated in Section \ref{sec:testing}.
\subsection{Density Forecast Evaluation: The Berkowitz Test}
\label{sec:testing}
The evaluation of the IlliQaR is based on the Berkowitz test \citep{Berkowitz:2001} for the adequacy of the density function with the realization of the dependent variable. Considering the flexibility of such a test, it can be applied to evaluate the fit of the full density as well as to specific quantile. In our application, we consider three quantiles, corresponding to the 1, 5 and 10\% of the Normal distribution. Importantly, the test is applied on the conditional density function (CDF) of $\text{Illiq}_t$ as expressed in
\begin{equation}
y_t=F(\text{Illiq}_t|\mathcal{F}_{t-1})=\int_{0}^{\text{Illiq}_t} f(x|\mathcal{F}_{t-1})dx.
\label{eq:cdf}
\end{equation}
Taking as an example the linear HAR model, the CDF is that of a Normal distribution. As for the MEM specifications, $F_{\text{MEM}}(\text{Illiq}_t|\mathcal{F}_{t-1})$ is given by the Gamma CDF. Finally, for the jump models in \eqref{eq:cdf}, the CDF is built as a mixture of a Gamma and a Kappa CDFs, as in Proposition \ref{theo1}.
The test is constructed under the null hypothesis that the model is correctly specified: if this is the case, the empirical CDF follows an uniform distribution, i.e. $y_t\sim U(0,1)$. The Berkowitz is therefore derived as a LR test for the comparison of the empirical CDF, $y_t$, and a new truncated variable, computed starting from the inverse of the standard normal distribution of $y_t$ itself. In other words, let the $s_t$ be defined as
\begin{equation*}
s_t=\Phi^{-1}(y_t),
\end{equation*}
we compute the new truncated variable, $s^*_t$ as
\begin{equation}
s^*_t=
\begin{cases}
\text{IlliQaR}_t(p) &\qquad \text{if} \quad s_t\leq \text{IlliQaR}_t(p)\\
s_t \qquad &\qquad \text{if} \quad s_t> \text{IlliQaR}_t(p).
\end{cases}
\end{equation}
The (right) tail coverage test can be derived using the LR principle based on
the censored normal density of $s^*_t$. Under the null of correct tail
coverage the test statistic is distributed as $\chi^2(2)$ as it corresponds
to a test for comparing the mean and variance estimated from
the truncated likelihood of $s^*_t$ to those expected under correct
coverage, i.e. zero mean and unit variance.
\section{Empirical Analysis}\label{sec:empirical_analysis}
We conduct an empirical assessment of IlliQaR by examining a sample of daily illiquidity measures of the S\&P 500 from January 3, 2005 until October 15, 2021, and of 25 US stocks from January 3, 2012 to January 10, 2024.
\subsection{Data}\label{sec:dataset}
Our empirical analysis focuses on the assessment of daily time series of the realized Amihud of the U.S. stock market, as represented by the S\&P 500 index. In analogy with \cite{RanaldoSantucci2020} and \cite{lacava2023realized}, daily illiquidity is computed as the ratio between the daily realized volatility (RV), which is the square root of the realized variance $RV=\sqrt{\sum_{i=1}^{I}r_i^2}$ \citep[computed using returns sampled at 5-minutes intervals, see][among many others]{liu2015does}, and the daily volume.\footnote{In practice, realized power variation (RPV) and realized variance (RV) exhibit an extremely high correlation, with an $R^2$ exceeding 99\%. Consequently, substituting RPV with RV in \eqref{eq:realized_Amihud} has a negligible impact on the resulting illiquidity series. As an alternative to return-based proxies, \cite{lacava2023realized} propose a low-frequency yet efficient estimator based on the daily price range—the difference between the daily high and low log-prices. \cite{hafner2023dynamic, hafner2025permanent} have successfully employed this estimator for illiquidity prediction.} Our data for RV are sourced from the Oxford-Man Institute, while daily trading volume is obtained from Datastream.
\begin{figure}[h!]
\centering
\begin{subfigure}
{\includegraphics[height=9cm,width=16cm]{Plot/realized_vs_daily_spx.pdf}}
\end{subfigure}
\caption{Realized (black solid line) and daily Amihud (red dashed line) series for the S\&P 500. \label{fig:amihud}}
\end{figure}
Figure \ref{fig:amihud} displays the time series of the realized Amihud (in black) for the S\&P 500. The persistence feature of $\text{Illiq}_t$ is evident at first sight, with long periods of low illiquidity followed by long periods of high illiquidity. Within the volatility literature, this is a well-established phenomenon known as {\it volatility clustering} (although in our case it is more appropriate to define it as {\it illiquidity clustering}). Altogether, this gives us a justification for considering well-established volatility models for predicting illiquidity. Furthermore, the series exhibits frequent spikes, particularly during the early portion of the sample which coincides with periods of financial turbulence. Indeed, elevated illiquidity levels have characterized much of the past fifteen years, driven by a succession of systemic shocks following the 2008 Global Financial Crisis, including the European Sovereign Debt Crisis and the COVID-19 pandemic.
Figure \ref{fig:amihud} also reports the daily Amihud (in red). The comparison with the realized Amihud signals two important features of these illiquidity measurements. First, realized Amihud and daily Amihud follow similar dynamic patterns. Second, the latter is much noisier than the former, thus confirming the efficiency gains in measuring illiquidity using high-frequency prices.
Table \ref{tab:stats} reports the sample statistics of realized Amihud, RV, trading volume, daily Amihud, and absolute returns for the S\&P 500. The realized Amihud displays kurtosis in excess compared to the reference value of the Normal distribution, as a consequence of extremely large realizations. These extreme realizations, typically originating from sharp illiquidity episodes, make the realized Amihud positively skewed, thus calling for a proper modeling framework. For instance, the MEM-J is expected to be able to assign the correct probability to these events. Indeed, due to the presence of measurement error, the traditional daily Amihud measure exhibits higher volatility, as indicated by its standard deviation, compared to the realized Amihud. Furthermore, the daily Amihud displays lower kurtosis, suggesting it is less effective at isolating extreme illiquidity realizations due to the noisy measurement of the latter. In the subsequent sections, we evaluate which econometric specifications best capture the dynamic and probabilistic properties of illiquidity to provide robust forecasts of future liquidity conditions.
\begin{table}[h!]
\centering
\setlength{\tabcolsep}{8pt}
\begin{tabular}{lccccc}
\toprule
& RV & $|r|$ & Volume & $\text{Realized Amihud}$ & $\text{Daily Amihud}$ \\
\midrule
Mean & 0.0080 & 0.0086 & 3.8015 & 0.0020 & 0.0018 \\
Median & 0.0061 & 0.0055 & 3.6412 & 0.0018 & 0.0012 \\
St. Deviation & 0.0065 & 0.0105 & 0.0065 & 0.0011 & 0.0018 \\
Skewness & 3.6558 & 3.5439 & 1.1037 & 2.4919 & 2.7596 \\
Kurtosis & 24.151 & 24.142 & 5.628 & 18.183 & 20.632 \\
\midrule
N. obs & 4216 & & & & \\
\bottomrule
\end{tabular}
\caption{S\&P 500 descriptive statistics. Realized Amihud and daily Amihud are scaled by a factor of 1.000E9, absolute returns are scaled by a factor of $\sqrt{\pi/2}$ while volume are divided by a factor of 1.000E8. \label{tab:stats}}
\end{table}
\subsection{Estimation results}\label{sec:results}
The estimation of the various MEM specifications presented in Section \ref{sec:IlliQaR} is performed on the \textit{full-sample} period, which, in the case of the S\&P 500, covers the years between January 3, 2005 to October 15, 2021. Panel a) of Table \ref{tab:results} shows the estimated coefficients for the S\&P 500 index. To model the dynamics of $\mu_t$, we consider a range of specifications including the linear (HAR) and non-linear (MEM) as discussed in Section \ref{sec:lin-nonlin}. Specifically, we consider the MEM and AMEM with both a GARCH(1,1) and HAR structures, alongside several jump-augmented specifications, namely the AMEM-J, AMEM(1,2)-J, G-AMEM-HAR-J, MEM-HAR-J, and AMEM-HAR-J.
\begin{sidewaystable}
\begin{adjustbox}{max width=0.88\linewidth,center}
\begin{tabular}{lccccccccccccc|ccc}
\vspace{0.4cm}
Panel a)& \multicolumn{16}{c}{\Large\textbf{Parameter Estimates}}\\
\vspace{0.4cm}
& \multicolumn{13}{c}{\large{Realized Amihud}} & \multicolumn{3}{c}{\large{Daily Amihud}}\\
& (I) & (II) & (III) & (IV) & (V) & (VI) & (VII) & (VIII) & (IX) & (X) & (XI) & (XII) & (XIII) & (II)& (IV) & (XI) \\
\midrule
$\omega$ & $0.0002^a$ & $0.0002^a$ & $0.0001^a$ & $0.0001^a$ & $0.0001^a$ & $0.0001^a$ & $0.0001^a$ & $0.0001^a$ & $0.0001^a$ & $0.0001^a$ & $0.0001^a$ & $0.0001^a$ & $0.0001^a$ & $0.0004^a$ & $0.0001^a$ & $0.0001^a$ \\
& (0.0000) & (0.0000) & (0.0000) & (0.0000) & (0.0000) & (0.0000) & (0.0000) & (0.0000) & (0.0000) & (0.0000) & (0.0000) & (0.0000) & (0.0000) & (0.0001) & (0.0000) & (0.0000) \\[2mm]
$\alpha_d$ & $0.3814^a$ & $0.2802^a$ & $0.4055^a$ & $0.3077^a$ & $0.3720^a$ & $0.2741^a$ & $0.2909^a$ & $0.4092^a$ & $0.3022^a$ & $0.3183^a$ & $0.2771^a$ & $0.3686^a$ & $0.2603^a$ & $-0.1314^a$ & $0.0056$ & $0.0000$ \\
& (0.0623) & (0.0597) & (0.0234) & (0.0243) & (0.0259) & (0.0288) & (0.0260) & (0.0196) & (0.0186) & (0.0180) & (0.0201) & (0.0221) & (0.0246) & (0.0296) & (0.0106) & (0.0278) \\[2mm]
$\alpha_w$ & $0.3599^a$ & $0.4018^a$ & & & $0.4045^a$ & $0.4462^a$ & $0.0731^b$ & & & & $0.0771^b$ & $0.4174^a$ & $0.4686^a$ & $0.3656^a$ & & $0.0000$ \\
& (0.0848) & (0.0824) & & & (0.0338) & (0.0341) & (0.0373) & & & & (0.0360) & (0.0330) & (0.0344) & (0.0541) & & (0.0258) \\[2mm]
$\alpha_m$ & $0.1772^a$ & $0.1784^a$ & & & $0.1617^a$ & $0.1662^a$ & $0.0692^a$ & & & & $0.0972^a$ & $0.1676^a$ & $0.1732^a$ & $0.4737^a$ & & $0.0432^a$ \\
& (0.0437) & (0.0415) & & & (0.0255) & (0.0254) & (0.0162) & & & & (0.0162) & (0.0275) & (0.0299) & (0.0766) & & (0.0132) \\[2mm]
$\beta_1$ & & & $0.5571^a$ & $0.6112^a$ & & & $0.4481^a$ & $0.5613^a$ & $0.6230^a$ & $0.4724^a$ & $0.4649^a$ & & & & $0.8961^a$ & $0.8411^a$ \\
& & & (0.0260) & (0.0260) & & & (0.0430) & (0.0212) & (0.0199) & (0.0376) & (0.0386) & & & & (0.0100) & (0.0195) \\[2mm]
$\beta_2$ & & & & & & & & & & $0.1318^a$ & & & & & & \\
& & & & & & & & & & (0.0332) & & & & & & \\[2mm]
$\gamma$ & & $0.1134^a$ & & $0.0892^a$ & & $0.0980^a$ & $0.1105^a$ & & $0.0948^a$ & $0.1018^a$ & $0.1175^a$ & & $0.1054^a$ & $0.1479^a$ & $0.1421^a$ & $0.1644^a$ \\
& & (0.0166) & & (0.0084) & & (0.0101) & (0.0101) & & (0.0067) & (0.0071) & (0.0082) & & (0.0092) & (0.0317) & (0.0166) & (0.0223) \\[2mm]
$\vartheta$ & & & $12.4790^a$ & $12.9010^a$ & $12.4745^a$ & $12.8165^a$ & $13.1089^a$ & $17.9103^a$ & $19.1291^a$ & $19.2810^a$ & $19.7753^a$ & $18.1070^a$ & $19.2893^a$ & & $1.2266^a$ & $21.6723^a$ \\
& & & (0.3845) & (0.4033) & (0.3688) & (0.3615) & (0.4239) & (0.5097) & (0.7224) & (0.5722) & (0.5661) & (0.0.5981) & (1.8122) & & (0.0260) & (0.1421) \\[2mm]
$\zeta$ & & & & & & & & $17.3756^a$ & $17.8975^a$ & $17.9585^a$ & $18.4562^a$ & $16.9308^a$ & $18.1812^a$ & & & $0.7757^a$ \\
& & & & & & & & (1.9085) & (0.5399) & (0.4641) & (0.3884) & (0.8604) & (2.8700) & & & (0.0395) \\[2mm]
$\phi_1$ & & & & & & & & $0.0139^b$ & $0.0119^a$ & $0.0111^a$ & $0.0076$ & $0.0099$ & $0.0099$ & & & $0.0344^b$ \\
& & & & & & & & (0.0056) & (0.0040) & (0.0040) & (0.0055) & (0.0078) & (0.0079) & & & (0.0135) \\[2mm]
$\phi_2$ & & & & & & & & $0.9422^a$ & $0.9560^a$ & $0.9588^a$ & $0.9712^a$ & $0.9572^a$ & $0.9617^a$ & & & $0.9896^a$ \\
& & & & & & & & (0.0222) & (0.0148) & (0.0148) & (0.0211) & (0.0346) & (0.0334) & & & (0.0041) \\[2mm]
$\phi_3$ & & & & & & & & $0.2181^a$ & $0.2139^a$ & $0.2011^a$ & $0.1420^a$ & $0.1548^a$ & $0.1491^b$ & & & $0.0224^c$ \\
& & & & & & & & (0.0364) & (0.0329) & (0.0311) & (0.0491) & (0.0542) & (0.0658) & & & (0.0111) \\[2mm]
\midrule
LogLik & 24753.09 & 24830.34 & 25759.31 & 25830.86 & 25758.53 & 25816.78 & 25865.27 & 25922.29 & 26013.86 & 26021.59 & 26052.75 & 25920.96 & 25995.17 & 21160.20 & 22861.33 & 22877.67 \\
\bottomrule
\\
Panel b)& \multicolumn{16}{c}{\Large\textbf{Ljung-Box test}}\\
\multicolumn{16}{l}{Residuals: $\epsilon_t$}\\[2mm]
LB 1 & 0.1015 & 0.0033 & 0.3999 & 0.3658 & 0.0321 & 0.0024 & 0.4420 & 0.4614 & 0.2808 & 0.8426 & 0.7770 & 0.0435 & 0.0045 & 0.3440 & 0.0063 & 0.0007\\[2mm]
LB 5 & 0.0000 & 0.0000 & 0.0163 & 0.0387 & 0.0000 & 0.0000 & 0.4748 & 0.0043 & 0.0129 & 0.0679 & 0.5702 & 0.0000 & 0.0000 & 0.0052 & 0.0805 & 0.0252 \\[2mm]
LB 10 & 0.0000 & 0.0000 & 0.0011 & 0.0054 & 0.0000 & 0.0000 & 0.3228 & 0.0001 & 0.0003 & 0.0004 & 0.3018 & 0.0000 & 0.0000 & 0.0131 & 0.1600 & 0.0526 \\[2mm]
\midrule
\end{tabular}
\end{adjustbox}
\caption{Estimates for S\&P 500. Panel a): Estimated coefficients (robust standard errors in parenthesis); Panel b): $p$-value of the Ljung-Box statistics. Sample period: January 3, 2005 - October 15, 2021. The superscripts a, b and, c denote significant coefficients at $1\%$, $5\%$ and, $10\%$ level, respectively. The estimated models are (I) HAR, (II) AHAR, (III) MEM, (IV) AMEM, (V) MEM-HAR, (VI) AMEM-HAR, (VII) G-AMEM-HAR, (VIII) MEM-J, (IX) AMEM-J, (X) AMEM(2,1)-J, (XI) G-AMEM-HAR-J, (XII) MEM-HAR-J, (XIII) AMEM-HAR-J.\label{tab:results}}
\end{sidewaystable}
The estimation of all model parameters is performed by ML and reported in Panel a) of Table \ref{tab:results}. The coefficients governing the conditional mean $\mu_t$ are all highly significant (at a 1\% significance level), with an average persistence of 95\% among the estimated models. As for $\hat\alpha_d$ coefficient, which summarizes the impact of news on market illiquidty, we find an average value of 0.33, with a peak at 0.41 for the MEM-J (model VII). As for the coefficient $\hat\gamma$, which measures the impact of bad news (in the form of negative returns) on illiquidity, it enters the model with a positive sign, pointing at a more heavily reaction of illiquidity against negative returns rather than against positive returns. This is also confirmed by the fact that, when the asymmetric effect is included in the model, $\hat\alpha_d$ reduces by 22\% on average, meaning that a large part of the impact of the price news on illiquidity is due to negative returns. We interpret this result as evidence in favor of the idea that illiquidity is strongly related to uncertainty among investors and that it causes a larger part of illiquidity persistence.
This translates into estimated values for $\hat\beta$ (0.53, on average) that are typically below the usual value found for volatility \citep[see][Ch. 1]{Bauwens:Hafner:Laurent:2012}. Finally, as expected, $\hat\vartheta$ (which measures the inverse of the variance of the error term) is higher for the more sophisticated MEM-J models, where a significant part of the variability of illiquidity is explained by the conditional jump component. This result emphasizes the need to adopt a model that is able to assign the correct probability to extremely large realizations of the realized Amihud measure. As for the specifications with a jump component, the process $\kappa_t$ for the jump intensity in the model is found to be a very persistent process ($\hat\phi_2$ is above 0.94 in all MEM-J specifications employing the realized Amihud), meaning that, similarly to the continuous part of illiquidity, also the number of jump arrivals strongly depends on its past realizations. The unconditional mean of the jump process, $\frac{\phi_1}{(1-\phi_2)}$, is between 0.23 and 0.27, corresponding to a jump every 4 days.
For comparison, we estimate the AHAR, AMEM, and G-AMEM-HAR-J specifications using the daily Amihud series (see the final three columns of Table \ref{tab:results}). The estimated coefficients for $\mu_t$ indicate that the daily series requires substantially more smoothing than its realized counterpart due to its more erratic nature. Specifically, the AMEM results yield a lower $\hat{\alpha}_d$ of $0.006$ (versus $0.306$ for the realized series) and a higher $\hat{\beta}_1$ of $0.896$ (versus $0.611$). Furthermore, the $\hat{\vartheta}$ estimate of $1.227$ is nearly an order of magnitude smaller than the $12.901$ obtained using realized Amihud. This disparity implies that the error term variance is approximately 11 times larger when the daily Amihud is used as the illiquidity proxy. These results highlight the significant precision gains achieved by utilizing high-frequency price data, which enhances model fit, improves explicability, and allows for a clearer separation between stochastic components of the model. As an additional illustration of the difficulty in disentangling noise from stochastic components, when estimating the G-AMEM-HAR-J on the daily Amihud series (last column of Table \ref{tab:results}), we find that $\hat{\phi}_1 \approx 0.03$ and $\hat{\phi}_2 \approx 0.99$. This leads to an unconditional jump arrival rate of nearly three jumps per day, signaling that liquidity jumps cannot be effectively isolated from the innovation component. Consequently, the log-likelihood remains essentially unchanged relative to a MEM without jumps, while the jump process variance is 25 times higher than that of the model estimated with the realized Amihud. This indicates that the jump component is poorly identified; the surge in variance primarily absorbs measurement noise rather than isolating genuine illiquidity jumps.
Panel b) of Table \ref{tab:results} presents diagnostic tests for all models, reporting $p$-values for the Ljung-Box test at 1, 5, and 10 lags to assess the properties of residuals. As far as models without jump are concerned, i.e. models (I)-(VII), we reject the null of no serial correlation of residuals at a 5\% significance level, with the only exception of the G-AMEM-HAR models (model VII) and for the HAR (I), the MEM (III), and the AMEM (IV), for the first lag. Results do not change if we consider jump models, i.e. models (VIII)-(XIII). In detail, we fail to reject the null hypothesis for the G-AMEM-HAR-J (model XI), which is the only model able to correctly account for the persistent nature of the illiquidity series. A similar performance is observed for the AMEM(1,2)-J (model X), for which we do not reject the null at a 5\% level for lag 1 and 5, while residuals serial correlation is detected for lag 10. Once again, we interpret this result as an evidence in favor of the view that illiquidity needs to be modeled through a model that is able to assign the correct probability to extreme realizations. In Section \ref{sec:forecast}, models are compared by evaluating their out-of-sample forecasting capability by means of the model confidence set procedure.
\subsection{Forecasting Illiquidity}\label{sec:forecast}
This section presents the results of our forecasting analysis, which evaluates the relative performance of the candidate econometric models in predicting realized illiquidity. Accurate illiquidity forecasts are of great importance to both market makers and investors. For market makers, high liquidity facilitates the enforcement of no-arbitrage conditions, thereby ensuring price efficiency. For investors, precise forecasts are essential given the critical role that liquidity risk plays in optimal portfolio construction and risk management.
To evaluate the predictive power of the candidate models, we conduct an out-of-sample forecasting exercise across three distinct horizons: $h=1, 5,$ and $22$ days, corresponding to daily, weekly, and monthly frequencies. The out-of-sample period spans from October 16, 2015, to October 15, 2021. Forecasts are generated using coefficients estimated over the in-sample window (January 3, 2005, to October 15, 2015). The in-sample parameter estimates, as reported in Table \ref*{tab:results_insample} in the Supplementary Document, are highly stable and do not exhibit significant deviations from the full-sample results, as shown in Table \ref{tab:results}.
To statistically evaluate the relative predictive accuracy of the candidate models, we employ the model confidence set (MCS) procedure developed by \cite{Hansen:Lunde:Nason:2011}. This procedure identifies the superior set of models that exhibit statistically indistinguishable forecasting performance at a given significance level. Following \cite{Patton:2011}, we utilize the QLIKE loss function, as it is robust to the choice of the proxy for the latent variable. The MCS is constructed using a 10\% significance level, and our results identify the models that are not outperformed by any others in the set.
\begin{table}[h!]
\centering
\setlength{\tabcolsep}{3pt}
\begin{adjustbox}{max width=0.95\linewidth,center}
\begin{tabular}{lcllcllc}
\multicolumn{2}{c}{1-step ahead} & & \multicolumn{2}{c}{5-step ahead} & & \multicolumn{2}{c}{22-step ahead} \\[0.5mm]
\midrule
\multicolumn{1}{c}{Model} & \multicolumn{1}{c}{p-value} & \multicolumn{1}{c}{} & \multicolumn{1}{c}{Model} & \multicolumn{1}{c}{p-value} & \multicolumn{1}{c}{} & \multicolumn{1}{c}{Model} & \multicolumn{1}{c}{p-value} \\[0.5mm]
MEM & 0.0012 & & AMEM-J & 0.0000 & & AMEM-J & 0.0000 \\
AMEM & 0.0027 & & AMEM & 0.0000 & & AMEM(2,1)-J & 0.0000 \\
MEM-HAR & 0.0036 & & MEM-J & 0.0009 & & MEM-J & 0.0000 \\
HAR & 0.0080 & & AMEM(2,1)-J & 0.0015 & & AMEM & 0.0000 \\
MEM-J & 0.0080 & & MEM & 0.0024 & & G-AMEM-HAR-J & 0.0000 \\
AMEM-J & 0.0080 & & G-AMEM-HAR-J & 0.0073 & & MEM & 0.0000 \\
AMEM-HAR & 0.0145 & & G-AMEM-HAR & 0.0203 & & MEM-HAR-J & 0.0000 \\
MEM-HAR-J & 0.0145 & & HAR & 0.0717 & & AMEM-HAR-J & 0.0000 \\
AMEM(2,1)-J & 0.0145 & & MEM-HAR-J & 0.0717 & & G-AMEM-HAR & 0.0016 \\
AHAR & 0.0488 & & MEM-HAR & 0.0717 & & AMEM-HAR & 0.0144 \\
AMEM-HAR-J & 0.0488 & & {\textbf{AHAR$^*$}} & {\textbf{0.7437$^*$}} & & MEM-HAR & 0.0211 \\
G-AMEM-HAR & 0.0488 & & {\textbf{ AMEM-HAR-J$^*$}} & {\textbf{0.7437$^*$}} & & AHAR & 0.0211 \\
{\textbf{G-AMEM-HAR-J$^*$}} & {$\mathbf{1.0000^*}$ } & & \multicolumn{1}{l}{{\textbf{AMEM-HAR$^*$}}} & $\mathbf{1.0000^*}$& & { \textbf{HAR$^*$}} & $\mathbf{1.0000^*}$ \\
\bottomrule
\end{tabular}
\end{adjustbox}
\caption{Forecast analysis for S\&P 500. Model Confidence Set for the out-of-sample forecasting performance (for 1, 5 and 22 step-ahead). Significance level 10\% (best set of models in bold and identified by an asterisk). Loss function: QLike. Estimation period: January 3, 2005 - October 15, 2015. Out-of-sample period: October 16, 2015 - October 15, 2021.\label{tab:mcs} }
\end{table}
Table \ref{tab:mcs} summarizes the results of the MCS procedure (best set of models in bold and identified by an asterisk). Our results underscore the superior forecasting performance of models incorporating jump dynamics, particularly at the short-term (1-step-ahead) horizon. Specifically, for the S\&P 500, the G-AMEM-HAR-J emerges as the unique best-performing model for the 1-step-ahead horizon, with no other specifications included in the MCS. However, as the forecasting horizon extends to the medium (5-step) and long term (22-step), the importance of the jump component diminishes. At these horizons, the AMEM-HAR and the linear HAR emerge as the top-performing models, respectively. As the forecasting horizon increases, the exclusion of jump models from the superior set of models highlights a shift in dynamics: jumps are critical for identifying `flash' shocks, but less relevant for long-term prediction because illiquidity eventually reverts to its historical average.
\section{Backtesting Illiquidity-at-Risk}\label{sec:illiquidity_determinats}
In this section, we assess the ability of the considered models to provide an accurate characterization of illiquidity risk. This assessment centers on the models' capacity to provide adequate coverage for the probability of extreme illiquidity spikes as measured by the Illiquidity-at-Risk (IlliQaR) metric defined in Section \ref{sec:IlliQaR}.
\begin{figure}[h!]
\includegraphics[height=9cm,width=16cm]{Plot/IaR_0.05_spx.pdf}
\caption{Realized Amihud series (black solid line) and 95\% out-of-sample IlliQaR forecasts from the G-AMEM-HAR-J (red dashed line) and G-AMEM-HAR (blue dashed line) models. \label{fig:IaR}}
\end{figure}
Figure \ref{fig:IaR} illustrates the out-of-sample 5\% IlliQaR for two competing specifications. The first specification, the G-AMEM-HAR-J model (red dashed line), explicitly accounts for discontinuous dynamics through the inclusion of jump components; the expression for the IlliQaR in this framework is derived from Proposition \ref{theo1}. The second is the G-AMEM-HAR baseline specification (blue dashed line), which relies on a standard Gamma MEM framework with the corresponding IlliQaR expression provided in \eqref{eq:illiqar_mem}. While baseline IlliQaR trajectories for both specifications are similar during calm periods, the jump-augmented model is significantly more responsive to new information. Consequently, the G-AMEM-HAR-J specification assigns higher probabilities to tail events, especially during market turmoil, with high illiquidity and increased systemic risk.
\begin{table}[h!]
\centering
\small
\begin{tabular}{lcccclccc}
\toprule
& \multicolumn{3}{c}{\textbf{No-Jump Models}} & & & \multicolumn{3}{c}{\textbf{Jump Models}} \\
\cmidrule(lr){2-4} \cmidrule(lr){7-9}
& 1\% & 5\% & 10\% & & & 1\% & 5\% & 10\% \\
\midrule
\multicolumn{9}{c}{\textit{Panel A: Full Sample}} \\
\midrule
HAR & 0.0000 & 0.0000 & 0.0000 & & MEMJ & 0.1826 & 0.4070 & 0.2392 \\
AHAR & 0.0000 & 0.0000 & 0.0000 & & AMEM-J & 0.4627 & 0.0073 & 0.7580 \\
MEM & 0.0000 & 0.0000 & 0.0000 & & AMEM(2,1)-J & 0.4708 & 0.0017 & 0.4868 \\
AMEM & 0.0000 & 0.0000 & 0.0000 & & G-AMEM-HAR-J & 0.3896 & 0.0183 & 0.4171 \\
MEM-HAR & 0.0000 & 0.0000 & 0.0000 & & MEMJ HAR & 0.1133 & 0.0113 & 0.0494 \\
AMEM-HAR & 0.0000 & 0.0000 & 0.0000 & & MEMJ AHAR & 0.0406 & 0.0038 & 0.1455 \\
G-AMEM-HAR & 0.0000 & 0.0000 & 0.0000 & & & & & \\
\midrule
\multicolumn{9}{c}{\textit{Panel B: Out-of-Sample}} \\
\midrule
HAR & 0.0000 & 0.0000 & 0.0000 & & MEMJ & 0.0095 & 0.0010 & 0.0012 \\
AHAR & 0.0000 & 0.0000 & 0.0000 & & AMEM-J & 0.0155 & 0.0001 & 0.0005 \\
MEM & 0.0462 & 0.0001 & 0.0000 & & AMEM(2,1)-J & 0.0160 & 0.0001 & 0.0008 \\
AMEM & 0.0118 & 0.0000 & 0.0000 & & G-AMEM-HAR-J & 0.0170 & 0.0012 & 0.0111 \\
MEM-HAR & 0.0033 & 0.0008 & 0.0000 & & MEMJ HAR & 0.0090 & 0.0033 & 0.0073 \\
AMEM-HAR & 0.0013 & 0.0000 & 0.0000 & & MEMJ AHAR & 0.0031 & 0.0025 & 0.0199 \\
G-AMEM-HAR & 0.0004 & 0.0000 & 0.0000 & & & & & \\
\bottomrule
\end{tabular}
\caption{Table reports the p-values of the Berkowitz test for the S\&P 500 full sample and out-of-sample IlliQaR($p$) computed at $p=1\%, 5\%, 10\%$.}
\label{tab:berkowitz_full_out}
\end{table}
The top panel of Table \ref{tab:berkowitz_full_out} reports the $p$-value of the Berkowitz test for the full-sample IlliQaR computed at the 1\%, 5\% and 10\% of the Gaussian distribution. The results reveal a clear hierarchy: regardless of the specification chosen for the conditional mean $\mu_t$, jump-augmented models consistently outperform their continuous-path counterparts. Across nearly all quantiles, jump-augmented specifications successfully characterize upper-tail dynamics, yielding $p$-values that consistently exceed the 0.10 significance threshold (marginal exceptions are limited to the 5\% quantile). We interpret this as robust evidence that realized illiquidity is driven by a stochastic process with a significant jump component; failing to account for these discontinuities leads to a systematic mispricing of tail risk.
These limitations are further corroborated by the out-of-sample results presented in the bottom panel of Table \ref{tab:berkowitz_full_out}. While no-jump models are universally rejected across most specifications, largely due to their inability to account for sudden and discontinuous liquidity evaporation, the jump-augmented specifications maintain significantly higher $p$-values, particularly at the 1\% level. This disparity confirms that modeling discrete shocks within the MEM-J framework is essential for achieving accurate probability coverage. By explicitly identifying and incorporating these jump components, the model successfully captures the non-linear dynamics and heavy-tailed distributions characteristic of extreme market stress. However, it is important to note that the superior performance of jump-augmented models is conditional upon the precision of the underlying illiquidity metric. While the MEM-J framework is theoretically superior for IlliQaR forecasting, its predictive gains are highly sensitive to the granularity of the input data; without precise measurement, the identification of genuine liquidity shocks is easily confounded by noise.
\begin{table}[h!]
\centering
\small
\begin{tabular}{lcccclccc}
\toprule
& \multicolumn{3}{c}{\textbf{Full Sample}} & & & \multicolumn{3}{c}{\textbf{Out-of-Sample}} \\
\cmidrule(lr){2-4} \cmidrule(lr){7-9}
& 1\% & 5\% & 10\% & & & 1\% & 5\% & 10\% \\
\midrule
AHAR & 0.000 & 0.000 & 0.000 & & AHAR & 0.000 & 0.000 & 0.000 \\
AMEM & 0.000 & 0.000 & 0.000 & & AMEM & 0.000 & 0.000 & 0.000 \\
JUMPS & 0.000 & 0.000 & 0.000 & & JUMPS & 0.000 & 0.000 & 0.000 \\
\bottomrule
\end{tabular}
\caption{Table reports the p-values of the Berkowitz test for the full sample and out-of-sample IlliQaR($p$) for $p=1\%, 5\%, 10\%$ for the daily Amihud. JUMPS denotes the G-AMEM-HAR-J model. \label{tab:berkowitz_daily}}
\end{table}
The necessity of a precise illiquidity measurement for analyzing tail risk is further reinforced by a comparison with the daily Amihud measure, which proves to be a significantly less efficient proxy for daily illiquidity than the realized Amihud. As reported in Table \ref{tab:berkowitz_daily}, the failure to reject the null hypothesis of the Berkowitz test suggests that the inherent noise in low-frequency measures effectively masks the true tail behavior of the process. Consequently, the daily Amihud is rendered unsuitable for providing reliable risk coverage, as it lacks the granularity required to distinguish between transitory fluctuations and genuine liquidity shocks.
\subsection{Price Jumps}
The results displayed in Tables \ref{tab:berkowitz_full_out} and \ref{tab:berkowitz_daily} point to an improvement in the quality of the IlliQaR coverage when illiquidity jumps are accounted for in the modeling, while simultaneously measuring illiquidity very precisely using the realized Amihud proxy. However, even accounting for liquidity jumps and precisely measuring illiquidity is not completely sufficient to achieve perfect coverage of extreme illiquidity tail events. As shown in Table \ref{tab:berkowitz_full_out}, during the out-of-sample period, even jump models frequently fail to provide correct coverage of the tails.
We argue that this is due to the potential interference of price jumps in the measurement of illiquidity. Specifically, when major public news hits the market (e.g., central bank announcements or earnings releases), reservation prices across all traders shift in the same direction to establish a new equilibrium. This common informational shock generates significant price volatility with minimal corresponding trading volume, as volume is primarily driven by trader disagreement (the investor-specific diffusive component), see among others \cite{boudt2014intraday} or \cite{scaillet2020high} in the context of the Bitcoin market. Consequently, standard realized volatility (RV) aggregates both the continuous diffusive component and discrete price jumps, causing artificial spikes in the realized Amihud measure that do not reflect baseline market illiquidity.
To isolate the pure diffusive component that genuinely drives trading volume and to assess how price jumps affect the adequacy of the \textsf{IlliQaR} framework, we replace RV with jump-robust proxies. Specifically, we replicate the analysis using the bipower variation \citep[BPV,][]{barndorff2004power} to explicitly purge the impact of information jumps (price jumps) on illiquidity measurement.
Since BPV is robust to price jumps by construction, the ratio $\text{Illiq}^C=BPV/\nu$ isolates the continuous component of illiquidity, filtering out the contribution of discontinuous price movements to the realized price impact. This robustness check allows us to disentangle two distinct sources of tail risk: genuine illiquidity jumps driven by trading frictions and the sudden withdrawal of liquidity providers, and price discontinuities that mechanically inflate the realized Amihud through the volatility numerator.
\begin{table}[h!]
\centering
\small
\begin{tabular}{lcccclccc}
\toprule
& \multicolumn{3}{c}{\textbf{No-Jump Models}} & & & \multicolumn{3}{c}{\textbf{Jump Models}} \\
\cmidrule(lr){2-4} \cmidrule(lr){7-9}
& 1\% & 5\% & 10\% & & & 1\% & 5\% & 10\% \\
\midrule
\multicolumn{9}{c}{\textit{Panel A: Full Sample}} \\
\midrule
HAR & 0.0000 & 0.0000 & 0.0000 & & MEM-J & 0.4473 & 0.3585 & 0.9001 \\
AHAR & 0.0000 & 0.0000 & 0.0000 & & AMEM-J & 0.2670 & 0.2158 & 0.7771 \\
MEM & 0.0000 & 0.0000 & 0.0000 & & AMEM(2,1)-J & 0.4659 & 0.5577 & 0.6699 \\
AMEM & 0.0000 & 0.0000 & 0.0000 & & G-AMEM-HAR-J & 0.6554 & 0.4369 & 0.4151 \\
MEM-HAR & 0.0000 & 0.0000 & 0.0000 & & MEM-HAR-J & 0.1530 & 0.3582 & 0.1490 \\
AMEM-HAR & 0.0000 & 0.0000 & 0.0000 & & AMEM-HAR-J & 0.5359 & 0.3383 & 0.1342 \\
G-AMEM-HAR & 0.0000 & 0.0000 & 0.0000 & & & & & \\
\midrule
\multicolumn{9}{c}{\textit{Panel B: Out-of-Sample}} \\
\midrule
HAR & 0.0000 & 0.0000 & 0.0000 & & MEM-J & 0.2303 & 0.0055 & 0.8925 \\
AHAR & 0.0000 & 0.0000 & 0.0000 & & AMEM-J & 0.1392 & 0.0692 & 0.6029 \\
MEM & 0.0000 & 0.0000 & 0.0000 & & AMEM(2,1)-J & 0.1400 & 0.0794 & 0.4231 \\
AMEM & 0.0000 & 0.0000 & 0.0000 & & G-AMEM-HAR-J & 0.2360 & 0.0480 & 0.5584 \\
MEM-HAR & 0.0000 & 0.0000 & 0.0000 & & MEM-HAR-J & 0.3583 & 0.0190 & 0.2455 \\
AMEM-HAR & 0.0000 & 0.0000 & 0.0000 & & AMEM-HAR-J & 0.1402 & 0.1171 & 0.4409 \\
G-AMEM-HAR & 0.0000 & 0.0000 & 0.0000 & & & & & \\
\bottomrule
\end{tabular}
\caption{Table reports the p-values of the Berkowitz test for the S\&P 500 full sample and out-of-sample IlliQaR($p$) computed at $p=1\%, 5\%, 10\%$, based on the jump robust realized Amihud (Illiq$_t^C$).}
\label{tab:berkowitz_full_out_bv}
\end{table}
Table \ref{tab:berkowitz_full_out_bv} reports the $p$-values of the Berkowitz test\footnote{Due to space constraints, the results for model estimation, forecast evaluation, and the determinants of Illiquidity at Risk are available upon request.} for the jump-robust IlliQaR computed at the 1\%, 5\%, and 10\% quantiles for the S\&P 500, for both the full sample and the out-of-sample period. The results confirm the same qualitative hierarchy documented in Table \ref{tab:berkowitz_full_out}: models without a jump component are systematically rejected across all quantiles, while jump-augmented specifications provide adequate coverage. Notably, the $p$-values obtained under the jump-robust measure are generally higher than those obtained with the realized Amihud, in both the full sample and the out-of-sample period. For instance, at the 1\% quantile, the MEM-J $p$-value rises from 0.183 to 0.447, and the AMEM-HAR-J from 0.041 to 0.536. This improvement suggests that discontinuous price movements introduce additional dispersion into the tail of the illiquidity distribution that the MEM-J framework only partially accommodates. When this source of variation is filtered out at the measurement stage, the model's task of accurately characterizing extreme illiquidity events becomes marginally easier. This improvement is also confirmed in the out-of-sample analysis, where a $p$-value lower than 1\% is detected only for the MEM-J at the 5\% quantile.
\subsection{Illiquidity at Risk and Its Determinants}\label{sec:determinants}
Finally, we investigate the economic determinants of IlliQaR violations to identify the factors driving extreme illiquidity episodes. We consider a Logit model for the binary dependent variable, $I_t$, which equals 1 if $\text{Illiq}_t>\text{IlliQaR}_t(p)$ and 0 otherwise. Building on our finding that illiquidity dynamics are closely linked to market uncertainty, we consider the following logit specification
\begin{equation}Pr(\text{I}_t=1 | \mathbf{x}_{t-1}) = \frac{\exp(\delta_0 + \mathbf{x}_{t-1}'\delta)}{1 + \exp(\delta_0 + \mathbf{x}_{t-1}'\delta)},
\label{eq:logit}
\end{equation}
where $\delta$ is a vector of coefficients, and $\mathbf{x}_{t-1}$ is a vector including several macro-financial variables, which we use to explore the determinants of tail liquidity risk. The Cboe Volatility Index (VIX) serves as a proxy for forward-looking market sentiment and implied volatility. The Economic Policy Uncertainty (EPU) index of \cite{Baker:Bloom:Davis:2016} accounts for broader macroeconomic and political risk. The TED spread is a proxy for systemic stress within the interbank lending market, specifically to account for fluctuations in funding liquidity. Furthermore, we introduce a binary indicator, $D$, which equals 1 for negative returns. This allows the model to isolate the asymmetric sensitivity of illiquidity to downward price innovations, a phenomenon often observed during market panics.
Following the methodology of \cite{Bekaert2013}, we also decompose the VIX index into two distinct components: market uncertainty and risk aversion. As a preliminary step, we estimate a linear model by regressing S\&P 500 realized volatility ($\text{RV}_t$) on its lagged values and the lagged VIX
\begin{equation*}
\text{RV}_t = \beta_0 + \beta_1 \text{RV}_{t-1} + \beta_2 \text{VIX}_{t-1} + \epsilon_t.\end{equation*}
The fitted values from this regression serve as our proxy for uncertainty, defined as $\text{Uncertainty}_t := \widehat{E(\text{RV}_t | \mathcal{F}_{t-1})}$, leading to
$$\text{Uncertainty}_t = \hat{\beta}_0 + \hat{\beta}_1 \text{RV}_{t-1} + \hat{\beta}_2 \text{VIX}_{t-1}.$$
The residual difference between the squared VIX and this uncertainty measure represents our proxy for risk aversion
$$\text{Risk Aversion}_t = \text{VIX}^2_t - \text{Uncertainty}_t.$$
We then replace the VIX with these two components, resulting in the following logit specification
\begin{equation}\text{Pr}(\text{I}_t=1 | \mathbf{x}_{t-1}) = \frac{\exp(\delta_0 + \delta_1 \text{Uncertainty}_{t-1} + \delta_2 \text{Risk Aversion}_{t-1} + \delta_3 \text{EPU}_{t-1} + \delta_4 \text{TED}_{t-1} + \delta_5 D_{t-1})}{1 + \exp(\delta_0 + \delta_1 \text{Unc}_{t-1} + \delta_2 \text{Risk Aversion}_{t-1} + \delta_3 \text{EPU}_{t-1} + \delta_4 \text{TED}_{t-1} + \delta_5 D_{t-1})}, \label{eq:logit_2}\end{equation}
Panel (a) of Table \ref{tab:logit_horiz} reports the estimation results for the full sample (columns 1–2) and the out-of-sample analysis (columns 3–4). In the full sample, the only significant predictor is the dummy variable ($D_{t-1}$), with a marginal effect of 0.016. This suggests that the probability of IlliQaR violations increases following negative returns. During the out-of-sample period, macro-financial risk factors, specifically the VIX, EPU, and TED spread, play a significant role in predicting unexpected illiquidity bursts (the LR test suggest a joint significant role of these factors). This shift may be attributed to the 2015–2021 period, which was characterized by increased political uncertainty and potential shifts in market structure and news reactivity. Recognizing that market participants may react to information without delay, Panel (b) presents results using contemporaneous rather than lagged regressors. In this specification, the dummy variable remains significant, with a marginal effect of approximately 0.05. Furthermore, the VIX becomes highly significant, particularly regarding its association with IlliQaR violations in the out-of-sample period. The decomposition of the VIX reveals a notable distinction: while the uncertainty coefficient is negative and large in magnitude, risk aversion is positive and significant across both the full sample and out-of-sample analyses. This provides a valuable insight for investors, suggesting they can better anticipate peaks in tail risk exceeding the threshold IlliQaR by monitoring specific observable signals, namely risk aversion, rather than the aggregate VIX alone. The TED spread is negatively and significantly associated with IlliQaR violations.
\begin{table}[h!]
\centering
\setlength{\tabcolsep}{3pt}
\begin{adjustbox}{max width=\linewidth,center}
\begin{tabular}{lccccccccc}
\toprule
& \multicolumn{4}{c}{Panel a) - Lagged Regressors} & & \multicolumn{4}{c}{Panel b) - Contemporaneous Regressors} \\
\cmidrule{2-5} \cmidrule{7-10}
& \multicolumn{2}{c}{Full sample} & \multicolumn{2}{c}{Out of sample} & & \multicolumn{2}{c}{Full sample} & \multicolumn{2}{c}{Out of sample} \\
\cmidrule{2-3} \cmidrule{4-5} \cmidrule{7-8} \cmidrule{9-10}
& \multicolumn{4}{c}{Parameter Estimates} & & \multicolumn{4}{c}{Parameter Estimates} \\
\midrule
Constant & $-3.2479^a$ & $-3.0703^a$ & $-4.1858^a$ & $-3.6973^a$ & & $-4.0424^a$ & $-3.1838^a$ & $-5.0866^a$ & $-3.1700^a$ \\
& (0.1501) & (0.2506) & (0.2822) & (0.3223) & & (0.1729) & (0.2774) & (0.4397) & (0.5902) \\
VIX & $0.0042$ & ---- & $0.0654^a$ & ---- & & $0.0131^b$ & ---- & $0.1116^a$ & ---- \\
& (0.0067) & ---- & (0.0182) & ---- & & (0.0065) & ---- & (0.0175) & ---- \\
& $[0.0002]$ & ---- & $[0.0021]$ & ---- & & $[0.0005]$ & ---- & $[0.0015]$ & ---- \\
Uncert. & & $-28.4122$ & ---- & $46.4553$ & & ---- & $-152.9560^a$ & ---- & $-99.4808$ \\
& ---- & (47.2839) & ---- & (54.3246) & & ---- & (52.0721) & ---- & (89.7807) \\
& ---- & $[-1.3664]$ & ---- & $[1.5490]$ & & ---- & $[-6.1645]$ & ---- & $[-1.3724]$ \\
Risk Av. & ---- & $0.0001$ & ---- & $0.0001$ & & ---- & $0.0012^a$ & ---- & $0.0021^a$ \\
& ---- & (0.0004) & ---- & (0.0004) & & ---- & (0.0004) & ---- & (0.0008) \\
& ---- & $[0.0000]$ & ---- & $[0.0000]$ & & ---- & $[0.0001]$ & ---- & $[0.0000]$ \\
EPU & $0.0008$ & $0.0013^c$ & $-0.0037^b$ & $-0.0015$ & & $0.0013^c$ & $0.0022^a$ & $-0.0048^b$ & $-0.0029$ \\
& (0.0007) & (0.0007) & (0.0018) & (0.0015) & & (0.0008) & (0.0007) & (0.0021) & (0.0020) \\
& $[0.0000]$ & $[0.0001]$ & $[-0.0001]$ & $[-0.0001]$ & & $[0.0001]$ & $[0.0001]$ & $[-0.0001]$ & $[0.0000]$ \\
TED & $0.0023$ & $0.0026$ & $0.0433^b$ & $0.0420^b$ & & $0.0024$ & $0.0025$ & $-2.9462^a$ & $-3.4556^a$ \\
& (0.0149) & (0.0137) & (0.0184) & (0.0193) & & (0.0170) & (0.0140) & (0.7220) & (0.8579) \\
& $[0.0001]$ & $[0.0001]$ & $[0.0014]$ & $[0.0014]$ & & $[0.0001]$ & $[0.0001]$ & $[-0.0391]$ & $[-0.0477]$ \\
D & $0.3220^b$ & $0.3153^b$ & $0.1928$ & $0.3487$ & & $1.1996^a$ & $1.1104^a$ & $1.5621^a$ & $1.5442^a$ \\
& (0.1402) & (0.1426) & (0.2854) & (0.2799) & & (0.1560) & (0.1551) & (0.3710) & (0.3783) \\
& $[0.0158]$ & $[0.0154]$ & $[0.0062]$ & $[0.0119]$ & & $[0.0544]$ & $[0.0487]$ & $[0.0245]$ & $[0.0251]$ \\
\midrule
& \multicolumn{4}{c}{Model Fit} & & \multicolumn{4}{c}{Model Fit} \\
\midrule
Success rate & 94.80\% & 94.80\% & 96.20\% & 96.30\% & & 94.80\% & 94.90\% & 96.40\% & 96.50\% \\
LR p-value & 0.1062 & 0.1389 & 0.0006 & 0.0507 & & 0.0000 & 0.0000 & 0.0000 & 0.0000 \\
\bottomrule
\end{tabular}
\end{adjustbox}
\caption{Logit model. The table displays logit estimation results for IlliQaR violations. Panel (a) reports coefficients for lagged regressors, while Panel (b) focuses on contemporaneous specifications. Robust standard errors in parenthesis. Average marginal effects in square brackets. The superscripts a, b and c denote significance at 1\%, 5\% and 10\% levels, respectively. Success rate denotes the percentage of correctly classified observations. LR p-value refers to the p-value of the likelihood ratio test for the joint significance of the regressors, where, under the null hypothesis, the model is $\text{Pr}(I_t=1 | \mathbf{x}_{t-1}) = \frac{\exp(\delta_0) }{1 + \exp(\delta_0)}$.\label{tab:logit_horiz}}
\end{table}
Finally, Table \ref*{tab:logit_ma_quantile} in the Supplementary Document reports results from the logit model where the explanatory variables (VIX, Uncertainty, Risk Aversion, EPU, and TED) are defined as deviations from their 22-day moving averages, restricted to observations above the 95th percentile (i.e., focusing on tail events within the regressors also). Intuitively, periods of elevated volatility and macroeconomic uncertainty are likely associated with higher margin requirements, which, in turn, can impair market liquidity. These results largely confirm the findings from Table \ref{tab:logit_horiz}, with one notable exception: once we isolate extreme realizations, volatility uncertainty, which captures the predictable component of the VIX, is no longer significant in the contemporaneous regression. In contrast, the risk aversion component remains highly significant. This suggests that IlliQaR violations are not driven by high volatility per se, but specifically by unexpected spikes in risk aversion. This indicates that market liquidity is most vulnerable when volatility realizations reflect a sudden shift in investor sentiment rather than just heightened fundamental uncertainty.
\section{Individual Stocks}\label{sec:individual}
To assess the robustness of our findings beyond the aggregate market index, we extend the empirical analysis to a cross-section of 25 major individual U.S. equities. This expansion allows us to investigate whether the predictive power of jump dynamics in illiquidity forecasting remains persistent in the presence of idiosyncratic risk, which may interact differently with the market liquidity shocks, as illustrated by \cite{pastor2003liquidity}. The sample, which covers the years between January 3, 2012 until January 10, 2024 is split into an in-sample period from January 3, 2012 to February 4, 2020 and an out-of-sample period from February 5, 2020 to January 10, 2024.
\subsection{Parameter Estimates}
Figure \ref{fig:boxplot_stocks} reports the box plot of the cross-sectional (full-sample) parameter estimates of the AMEM-HAR-J model for the 25 stocks under investigation.\footnote{Due to space constraints, we report the complete estimation results in the Supplementary Document. Specifically, Section \ref*{app:appendix_fullsample_individual_stocks} of the Supplementary Document provides parameter estimates for three representative specifications: the AHAR, AMEM-HAR, and AMEM-HAR-J models, considering a constant jump intensity specification ($\kappa_t=\kappa$). Furthermore, MCS results are reported in Tables \ref*{tab:mcs_1_individual}--\ref*{tab:mcs_22_individual} (Section \ref*{app:individual}), while the results of the Berkowitz test are in Sections \ref*{App:berkowitz_individual_full} and \ref*{App:berkowitz_individual_out} of the Supplementary Document.}
Across all 25 equities, the intercept $\omega$ ranges between 0.07 and 0.015, much larger values than those of the S\&P 500, signaling that the average illiquidity levels of individual stocks are significantly higher than that of the overall US market. As for the HAR components ($\alpha_d$, $\alpha_w$, and $\alpha_m$), they are generally positive and statistically significant at the 1\% level. This confirms that idiosyncratic illiquidity (much like market illiquidity) is characterized by strong persistence and long-memory behavior, with the sum of the three coefficients close to 1. However, there is a slight structural difference compared with the index: $\alpha_d$ is larger for individual stocks, whereas $\alpha_w$ is smaller. Finally, the asymmetry parameter ($\gamma$), whenever significant, is positive, confirming that the ``leverage effect"—where liquidity deteriorates more severely following negative returns—holds at the individual stock level. Compared with the S\&P 500 index, the strength of this leverage effect appears to be weaker. The distributional parameters governing the continuous process are remarkably stable across the cross-section. Under the AMEM-HAR-J specification, the shape parameter of the Gamma distribution ($\vartheta$) consistently falls between 13 and 25, closely mirroring the dispersion levels observed for the S\&P 500.
\begin{landscape}
\begin{figure}[htbp!]
\centering
\subfigure[$\omega$] {\includegraphics[width=0.35\textwidth,height=5cm]{fig_3a.eps}}
\subfigure[$\alpha_1$] {\includegraphics[width=0.35\textwidth,height=5cm]{fig_3b.eps}}
\subfigure[$\alpha_2$] {\includegraphics[width=0.35\textwidth,height=5cm]{fig_3c.eps}}
\subfigure[$\alpha_3$] {\includegraphics[width=0.35\textwidth,height=5cm]{fig_3d.eps}}\\
\subfigure[$\gamma$] {\includegraphics[width=0.35\textwidth,height=5cm]{fig_3e.eps}}
\subfigure[$\vartheta$] {\includegraphics[width=0.35\textwidth,height=5cm]{fig_3f.eps}}
\subfigure[$\zeta$] {\includegraphics[width=0.35\textwidth,height=5cm]{fig_3g.eps}}
\subfigure[$\kappa$] {\includegraphics[width=0.35\textwidth,height=5cm]{fig_3h.eps}}
\caption{Cross-sectional distribution of estimated coefficients across the 25 Dow Jones constituent stocks. Panels report the boxplot of each parameter of the AMEM-HAR-J across the 25 stocks. The black asterisk is the parameter estimate obtained on the S\&P 500 index.}\label{fig:boxplot_stocks}
\end{figure}
\end{landscape}
Instead, the most notable differences relate to the jump intensity and distribution parameters for individual stocks. Specifically, the jump size distribution parameter $\zeta$ across the 25 equities ranges between 4 and 9, substantially lower than the index-level estimate—while the jump intensity coefficients range between 1\% and 5\%, compared to around 22\% for the index. This signals that liquidity jumps occur less frequently for individual stocks than for the broad market, but they carry the potential to trigger larger liquidity dry-ups. Consequently, while market liquidity jumps reflect sustained, systemic stress episodes, single-stock liquidity evaporation is often driven by transient, idiosyncratic shocks that rapidly dissipate despite generating extreme short-term liquidity dry-ups.
\subsection{Forecasting Analysis}
The forecasting results for individual stocks largely corroborate the findings for the S\&P 500, albeit with the degree of heterogeneity expected at a disaggregated level. Figure \ref{fig:mcs_stocks} summarizes the MCS analysis for the 25 stocks under consideration.\footnote{For a complete overview of the MCS results, see Tables \ref*{tab:mcs_1_individual}--\ref*{tab:mcs_22_individual} of the Supplementary Document.}
\begin{figure}[h!]
\centering
\subfigure[1-step ahead]
{\includegraphics[width=0.49\textwidth]{fig_4a.eps}}
\subfigure[5-step ahead] {\includegraphics[width=0.49\textwidth]{fig_4b.eps}}\\
\subfigure[22-step ahead]
{\includegraphics[width=0.49\textwidth]{fig_4c.eps}}
\caption{Forecast analysis for Individual stocks. Model Confidence Set for the 1-, 5-, and 22-step ahead out-of-sample
forecasting performance. Significance level 10\%. The best model is highlighted in yellow; models included in the Model Confidence Set are shown in green; models excluded from the Model Confidence Set are shown in red.
Loss function: QLike. Estimation period: January 3, 2012 - February 4, 2020. Out-of-sample period: February
5, 2020 - January 10, 2024. Stocks are presented in descending order of market capitalization.}\label{fig:mcs_stocks}
\end{figure}
At the short-term horizon ($h=1$) (Panel a), models incorporating jump dynamics consistently populate the MCS. Specifications such as the G-AMEM-HAR-J and the AMEM(2,1)-J are frequently identified as top performers. For instance, for most stocks, the MCS procedure excludes simple linear HAR and standard MEM specifications, favoring instead jump-augmented models. This confirms that for individual assets, daily liquidity is characterized by abrupt, discontinuous changes that require explicit inclusion of jumps to achieve accurate one-step-ahead predictions. Consistent with the index-level analysis, the contribution of the jump component decays as the forecasting horizon extends. At the weekly and monthly horizons ($h=5$ and $22$, Panels b and c), the superior performance of the linear HAR model reaffirms its ability to capture the long-memory behavior of illiquidity through a parsimonious cascade of heterogeneous components. However, across all horizons, the G-AMEM-HAR-J is included in the MCS for most stocks, underscoring the need for accurate modeling of both the long-memory dynamics of illiquidity and its discontinuous, jumpy behavior.
\subsection{IlliQaR Analysis}
The necessity of jump models becomes even more apparent when evaluating the IlliQaR metric. Figure \ref{fig:pval_berkowitz_stocks_1} uses a heatmap to summarize the results of the Berkowitz test for the out-of-sample periods across the 25 stocks under analysis.\footnote{See Tables \ref*{tab:density_forecast_full_1_individual}--\ref*{tab:density_forecast_oos_10_individual} in Sections \ref*{App:berkowitz_individual_full} and \ref*{App:berkowitz_individual_out} of the Supplementary Document for the $p$-values of the Berkowitz test for density forecasts.} Standard linear models (HAR and AHAR) and continuous MEM specifications frequently fail to provide adequate coverage for 1\% tail events, with $p$-values often dropping to 0.0000 and leading to a rejection of the null hypothesis of correct specification. This indicates a systematic underestimation of the probability of extreme liquidity dry-ups. In contrast, jump-diffusion specifications (MEM-J class) outperform baseline models in accurately tracking the tail properties of realized illiquidity, particularly during periods of high market stress. For the majority of the analyzed stocks, jump-augmented models fail to reject the null hypothesis. This result is economically significant: it implies that ignoring the jump component in individual stocks exposes investors to a ``tail risk" that is far greater than what standard Gamma-distributed models predict.
\begin{figure}[h!]
\centering
\subfigure[1\% - Out-of-sample] {\includegraphics[width=0.49\textwidth]{fig_5a.eps}}
\subfigure[5\% - Out-of-sample] {\includegraphics[width=0.49\textwidth]{fig_5b.eps}}
\subfigure[10\% - Out-of-sample] {\includegraphics[width=0.49\textwidth]{fig_5c.eps}}
\caption{Heatmaps of Berkowitz test p-values for the 1\% IlliQaR across 25 individual U.S. equities. Higher p-values (green) indicate adequate tail coverage; lower p-values (red) indicate rejection of the null hypothesis of correct specification.}\label{fig:pval_berkowitz_stocks_1}
\end{figure}
In summary, the analysis of individual stocks confirms that IlliQaR is a jump-driven phenomenon. The ``evaporation" of liquidity is not merely an index-level occurrence; it is inherently discontinuous at the single-stock level as well. Consequently, risk management systems relying solely on continuous approximations of illiquidity are likely ill-equipped to handle periods of market stress. Our results suggest that individual stock IlliQaR violations often cluster during periods of S\&P 500 liquidity stress. This indicates that Illiquidity at Risk is a systemic concern rather than a localized one, with the main index acting as a leading indicator for extreme liquidity dry-ups in individual stocks.
\section{Conclusions}
\label{sec:conclusion}
This paper provides a comprehensive investigation into the dynamics of stock market illiquidity and the extent to which their level and tail behavior can be accurately predicted. Our primary contribution is the introduction of Illiquidity-at-Risk (IlliQaR), a novel risk metric designed to quantify the magnitude of extreme liquidity dry-ups. Our analysis utilizes a refinement of the classic Illiq proxy \citep{Amihud2002illiquidity}, the realized Amihud measure \citep{RanaldoSantucci2020}, which has been demonstrated to provide a precise assessment of illiquidity \citep{lacava2023realized}.
We employ a number of linear and non-linear econometric specifications designed to capture the defining stylized facts of illiquidity, namely positive skewness, excess kurtosis and long-range dependence. Furthermore, we evaluate the specific contribution of illiquidity jumps in shaping the distributional properties of the realized Amihud and the resulting IlliQaR estimates. Several key findings emerge from our empirical analysis of the S\&P 500 and 25 individual US stocks. First, we find that illiquidity is a highly persistent process that is more closely associated with market uncertainty than with risk aversion. Second, jumps constitute a significant component of illiquidity risk. Accurately assigning probabilities to extreme realizations of the realized Amihud is of fundamental importance, particularly for risk management and real-time monitoring of market stress. Regarding the forecasting horizon, our results reveal a clear distinction between short-term and long-term dynamics. For one-step-ahead forecasts, the AMEM-J specification—which explicitly accounts for abrupt, discontinuous jumps—emerges as the superior model. Conversely, at medium-to-long-term horizons, illiquidity tends toward its unconditional mean and the forecast trajectories become smoother. In these instances, the parsimonious HAR-type models, which effectively capture the pseudo long-memory properties of the series, provide better predictive performance. Finally, our investigation into the determinants of IlliQaR violations indicates that the probability of observing an extreme illiquidity event is state-dependent. It responds significantly to both negative lagged returns and shifts in market risk aversion, suggesting that investors can anticipate peaks in tail risk by monitoring these observable variables.
\bibliography{Bibliography}
\begin{appendix}
\section{Proof of Proposition \ref{theo1}}\label{app:proof}
The conditional density of $\text{Illiq}_t$ is defined as an infinite mixture:\begin{equation*}f_{\text{MEM-J}}(\text{Illiq}_t|\mathcal{F}_{t-1}) = e^{-\kappa_t} f_G(\text{Illiq}_t|N_t=0,\mathcal{F}_{t-1}) + \sum_{j=1}^\infty \frac{e^{-\kappa_t}\kappa_t^j}{j!} f_K(\text{Illiq}_t|N_t=j,\mathcal{F}_{t-1}),
\end{equation*}
where $f_G$ and $f_K$ are the conditional Gamma (no jumps) and Kappa (jumps occur) densities, respectively. By interchanging the integral and summation, the conditional CDF simplifies to
\begin{equation}
F_{\text{MEM-J}}(\text{Illiq}_t|\mathcal{F}_{t-1}) = e^{-\kappa_t} F_G(\text{Illiq}_t|N_t=0,\mathcal{F}_{t-1}) + \sum_{j=1}^\infty \frac{e^{-\kappa_t}\kappa_t^j}{j!} F_K(\text{Illiq}_t|N_t=j,\mathcal{F}_{t-1}),
\end{equation}
where $F_G(\cdot)$ and $F_K(\cdot)$ denote the respective Gamma and Kappa CDFs. In particular,
\[
F_G(\text{Illiq}_t|N_t=0,\mathcal{F}_{t-1})=\frac{1}{\Gamma(\vartheta)}\int_{0}^{\frac{\vartheta \text{Illiq}_t}{d_{\kappa_t}\mu_t}}s^{\vartheta-1}e^{-s}ds,
\]
is the cumulative distribution function of a Gamma-distributed random variable with shape $\vartheta$ and scale $\frac{d_{\kappa_t} \mu_t}{\vartheta}$.
As for the the K-distribution CDF, $F_K(y;\mu, \vartheta_1, \vartheta_2)$, expressed as a function of mean and the shape parameters is given by
\begin{equation}
F_K(y;\cdot) = \frac{2^{2-\vartheta_1-\vartheta_2}}{\Gamma(\vartheta_1)\Gamma(\vartheta_2)} \int_0^{2\sqrt{\vartheta_1\vartheta_2 y/\mu}} t^{\vartheta_1+\vartheta_2-1} K_{|\vartheta_1-\vartheta_2|}(t) dt,\end{equation}where $K_a(\cdot)$ is the modified Bessel function of the second kind. For the jump components $K(\text{Illiq}_t|N_t=j)$, we set $y = \text{Illiq}_t$, $\mu = \mu_t j d_{\kappa_t}$, $\vartheta_1 = j\zeta$, and $\vartheta_2 = \vartheta$.
\setcounter{section}{1}
\clearpage
\begin{landscape}
\begin{table}[!h]