EconBase
← Back to paper

The realized copula of volatility

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.

101,511 characters · 9 sections · 61 citation commands

Rendered from LaTeX for readability, not typeset faithfully. Citation keys are highlighted; maths is left as source; figures, tables and equation environments are summarised rather than reproduced; unrecognised commands are greyed out so nothing is silently dropped. Email addresses are removed.

The realized copula of volatility

abstractWe study a new measure of codependency in the second moment of a continuous-time multivariate asset price process, which we name the realized copula of volatility. The statistic is based on local volatility estimates constructed from high-frequency asset returns and affords a nonparametric estimator of the empirical copula of the latent stochastic volatility. We show consistency of our estimator with in-fill asymptotic theory, either with a fixed or increasing time span. In the latter setting, we derive a functional central limit theorem for the empirical process associated with the measurement error of the time-invariant marginal copula of volatility. We also develop a goodness-of-fit test to evaluate hypotheses about the shape of the latter. In a simulation study, we demonstrate that our estimator is a good proxy of both the empirical and marginal copula of volatility, even with a moderate amount of high-frequency data recorded over a relatively short sample. The goodness-of-fit test is found to exhibit size control and excellent power. We implement our framework on high-frequency transaction data from futures contracts that track the U.S. equity and treasury bond market. A Gumbel copula is found to offer a near-perfect bind between the realized variance processes in these data. JEL Classification: C14; C32; C58; C80. Keywords: Copula; empirical process; functional central limit theorem; high-frequency data; nonparametric estimation; stochastic volatility; tail dependence.

\thispagestyle{empty}

Introduction

\setcounter{page}{1}

Modeling the comovement of stochastic volatility across several assets is paramount to risk measurement and management, portfolio allocation, and multi-asset and derivative pricing. It has long been recognized that financial markets are interconnected, causing asset return volatilities---and correlations---to surge in tandem during financial distress, triggering significant volatility spillover effects, posing a potential systemic risk diebold-yilmaz:09a, billio-getmansky-lo-pelizzon:12a. An understanding of volatility codependency is also important for pricing of volatility risk (i.e., the variance risk premium) across asset classes carr-wu:09a, bollerslev-todorov:11a. The importance of the second moment structure of asset returns is further fueled by the so-called leverage effect that captures a negative association between return and volatility black:76a, christie:82a. Thus, from a risk aggregation view the relevant target is often the dependence structure in second moment of returns (namely volatilities and correlations), rather than in the raw returns.

In practice, the assessment of volatility comovements is complicated by the fact that volatility is unobserved, so it requires an empirical route based on an observable proxy. It therefore remains a relatively unexplored area. In recent years, however, the availability of high-frequency data has made it possible to estimate volatility nonparametrically with the realized variance andersen-bollerslev:98a, barndorff-nielsen-shephard:02a. Most of this literature concentrates on estimation of integrated measures of variation, defined as an integral of a smooth function of the latent spot variance process jacod-rosenbaum:13a. However, as li-todorov-tauchen:13a point out, the mapping from the probability distribution of spot volatility to the distribution of integrated volatility is not one-to-one. Thus, the “integrated” approach may be partially informative about some distributional aspects of volatility and codependency, but it faces a limitation due to the inherent “identification failure.” This may weaken, or outright eliminate, power for testing against relevant alternatives. They suggest to work directly with the pathwise distribution of spot volatility via a so-called occupation measure geman-horowitz:80a, which recovers the empirical distribution function of the point-in-time volatility process li-todorov-tauchen:16a, li-todorov-tauchen:16b.\footnote{In addition, christensen-thyrsgaard-veliyev:19a and zhang-li-bollerslev:22a made the univariate framework robust to microstructure noise.} Regardless, the volatility occupation literature is entirely univariate and has concentrated on the marginal distribution, rather than looking at the interrelationship between various volatility processes.

The starting point of our analysis is the remarkable result from probability theory known as Sklar's theorem sklar:59a. It states that for any continuous random vector $(X_{1}, \dots, X_{d})$ with cumulative distribution function $H(x_{1}, \dots, x_{d}) = \mathbb{P}(X_{1} \leq x_{1}, \dots, X_{d} \leq x_{d})$, there exists a unique function $C:[0,1]^{d} \rightarrow [0,1]$ called a copula, such that $H(x_{1}, \dots, x_{d}) = C(H_{1}(x_{1}), \dots, H_{d}(x_{d}))$, where $H_{i}(x_{i}) = \mathbb{P}(X_{i} \leq x_{i})$ is the marginal distribution of $X_{i}$.\footnote{The concise mathematical definition is as follows: A $d$-dimensional copula $C:[0,1]^{d} \rightarrow [0,1]$ is the cumulative distribution function of a random vector $(U_{1}, \dots, U_{d})$, where the marginal distribution of $U_{i}$ is uniform, i.e. $C(u_{1}, \dots, u_{d}) = \mathbb{P}(U_{1} \leq u_{1}, \dots, U_{d} \leq u_{d})$ with $\mathbb{P}(U_{i} \leq u_{i}) = u_{i}$, for $i = 1, \dots, d$ and $0 \leq u_{i} \leq 1$.} Conversely, given $C(u_{1}, \dots, u_{d})$ and $H_{1}(x_{1}), \dots, H_{d}(x_{d})$, $H(x_{1}, \dots, x_{d}) = C(H_{1}(x_{1}), \dots, H_{d}(x_{d}))$ is a probability distribution. Copulas provide a marginal-free representation of dependence in a multivariate setting, and they enable the researcher to model the marginal features of the data separately.\footnote{The term copula originates from the Latin word for “bond” or “tie”.} Indeed, many measures of concordance (and discordance), such as Kendall's tau and Spearman's rho, can be expressed entirely in terms of the copula (see Section (ref)). The latter are often more informative than the Pearson's coefficient of correlation, which captures linear dependence. As emphasized by embrechts-mcneil-straumann:99a, correlation is insufficient for risk management. This is explained by the fact that financial markets typically exhibit nonlinear dependence and, further, may possess a high degree of tail concentration (or corner clustering).\footnote{The application of copula theory had a revival in financial economics after li:00a proposed to model the default correlation in portfolios of bond loans with a Gaussian copula. This enabled tractable analysis of structured credit products (such as collateralized debt obligations). The model became a standard pricing tool on Wall Street. However, since the Gaussian copula has no tail dependence, it tends to understate the probability of concurrent defaults, as seen during the financial crisis salmon:09a. A comprehensive review of the application of copula theory to economic and financial time series is patton:12a.} Therefore, copulas---or descriptive statistics derived from them---are a better representation of dependence. Another advantage of rank-based statistics is their invariance to monotonic transformations (that do not alter the ordering of the data). However, a model-free technique for analyzing the multivariate dependence structure of volatility is still lacking. In particular, tools for estimating and conducting inference on the copula of volatility are not available.

We close that gap and propose to study the dependence structure of several stochastic volatility processes. We make a four-fold contribution to the literature. First, we extend li-todorov-tauchen:13a to the multivariate setting. We operate within a general, albeit standard, continuous-time arbitrage-free It\^{o} semimartingale framework for a (multi-dimensional) log-price process, where the random volatility can exhibit rather unrestricted dynamics. We replace the latent spot variance processes by localized, jump-robust, realized variance measures constructed from discretely observed high-frequency returns over short time intervals. The estimation errors are controllable by the sampling interval and bandwidth selection jacod-protter:12a. We allow for very active jumps in the price process. The latter can be regarded as “nuisance parameters” in our setting, and we handle them with a truncation device mancini:09a. Second, we construct a nonparametric estimator of volatility codependency, called the realized copula of volatility, enabling a meaningful cross-asset comparison irrespective of the marginal distribution of the individual volatility processes. Third, we develop the necessary large-sample theory for our estimator aligned with its target at each horizon. In the in-fill asymptotic limit with fixed time span, we show the realized copula of volatility converges, uniformly, in probability to its empirical counterpart that captures the volatility dependency over the time interval observed so far. Importantly, we allow for nonstationary volatility. In a double asymptotic in-fill and long-span setting---with stationary and weakly dependent volatility---we establish a functional central limit theorem in the form of weak convergence of the measurement error of the realized volatility copula. We show how this can be exploited to construct pointwise or, with a bit more effort, uniform confidence bands. A goodness-of-fit test for assessing the shape of the volatility copula is also proposed. The latter procedure can, for instance, be used to gauge the adequacy of various parametric copulas; an approach we look further at in our application. This builds upon and extends earlier results on the empirical copula process deheuvels:79a, fermanian-radulovic-wegkamp:04a. However, as explained above, our framework is distinctly more complicated than the discrete-time setup with independent and identically distributed data (see fermanian-scaillet:02a for related work in the time series setting). Fourth, we deliver an implementable inference procedure by designing a newey-west:94a-type heteroscedasticity and autocorrelation consistent estimator of the long-run asymptotic variance of our realized copula of volatility.

The small sample properties of our estimator are investigated through a Monte Carlo experiment. We demonstrate its efficacy over the entire support of the copula. In addition, we evaluate the goodness-of-fit test in order to illuminate its size and power. In our empirical work, we consider an application to U.S. equity and treasury futures markets. This leads to forceful evidence that a Gumbel copula---exhibiting upper-tail dependence---offers a near-perfect description of the realized copula of volatility. Thus, our nonparametric procedure complements the existing literature and facilitates a marginal-free analysis of volatility codependency, offering a more refined assessment of systemic risk transmission.

The paper progresses as follows. In Section (ref), we present the theoretical framework and define the empirical copula of the volatility process based on an occupation measure. In Section (ref), we adopt a spot volatility estimator based on discrete high-frequency data, from which we can construct the realized copula of volatility. In Section (ref), we state our assumptions and develop large-sample theory for our estimator. In particular, we show its consistency for the empirical copula of volatility in the in-fill limit with fixed time span. We also derive a functional central limit theorem for the empirical process associated with the estimation error of the time-invariant copula of volatility in a double asymptotic in-fill and long-span setting. We further propose an estimator of the long-run asymptotic variance and show how to exploit the functional analysis to construct uniform confidence bands. At last, we develop a goodness-of-fit test. Section (ref) reports a detailed Monte Carlo study to evaluate the finite sample properties of our realized copula of volatility, along with a parametric version of our goodness-of-fit test. In Section (ref), we present an empirical application. We conclude in Section (ref). The proofs of our statistical analysis are postponed to the \hyperref[app:proofs]{Appendix}.

The setting

We suppose that a filtered probability space (or stochastic basis) $\big( \Omega, \mathcal{F}, (\mathcal{F}_{t})_{t \geq 0}, \mathbb{P} \big)$ describes the evolution of a pair of log-price processes, $X$ and $Y$, that evolve in continuous-time over an interval $[0,T]$.\footnote{Our theoretical framework extends trivially to $d$-dimensional processes, for any fixed $d \geq 2$. We concentrate on the bivariate setting, because it conveys the main idea and avoids the cost of extra notation.} Here, $\mathcal{F}_{s} \subseteq \mathcal{F}_{t} \subseteq \mathcal{F}$ for $s \leq t \leq T$ is a filtration that represents past and current information available to market participants at any time.

As consistent with the “no free lunch with vanishing risk” principle from financial economics, we assume that the motion of $(X,Y)$ is described by an It\^{o} semimartingale delbaen-schachermayer:94a.\footnote{The It\^{o} semimartingale assumption is commonly motivated in economics by the efficient markets hypothesis and rational expectations theory samuelson:65a.} This means that it has the component-wise representation

equation[equation omitted — 306 chars of source]

where the drift $(b_{t}^{Z})_{t \geq 0}$ and the volatility $( \sigma_{t}^{Z})_{t \geq 0}$ are adapted c\`{a}dl\`{a}g processes, while $(W_{t}^{Z})_{t \geq 0}$ is a standard Brownian motion. In addition, $\mu^{Z}$ is a Poisson random measure on $( \mathbb{R}_{+}, \mathbb{R})$ with compensator $\nu^{Z}( \mathrm{d}t, \mathrm{d}z) = dt \otimes \lambda^{Z}( \mathrm{d}z)$, $\lambda^{Z}$ is a $\sigma$-finite measure on $\mathbb{R}$, and $\delta^{Z}: \Omega \times \mathbb{R}_{+} \times \mathbb{R} \rightarrow \mathbb{R}$ is predictable. Here, and in other places, we employ a generic process $Z$ (typically with $Z = X,Y$) to avoid repeating definitions.

In general, a semimartingale can be decomposed into a finite variation component and a local martingale. The “It\^{o}” classifier imposes an absolute continuity condition on these, so that we can express them as integrals of “spot” processes. This restriction is common, because the latter are more amenable to statistical analysis from high-frequency data. Apart from that, our framework is nonparametric, model-free, and suffices to capture the most dominant features of empirical asset price processes, such as time-varying expected returns, stochastic volatility, and leverage effect. Moreover, it encompasses the presence of large and small price jumps, where the latter can be infinitely active and of infinite variation. We also permit complex within- and cross-asset dependencies between the various model components. For example, the standard Brownian motions can be correlated (i.e. $\rho_{t} \mathrm{d}t = \mathrm{d} \langle W^{X}, W^{Y} \rangle_{t}$, where $\langle W^{X}, W^{Y} \rangle_{t}$ is the quadratic covariation process). Moreover, the jumps in volatility can be related to the occurrence of jumps in the log-price, and there can be common jumps in both todorov-tauchen:11a, bibinger-winkelmann:18a. We are, however, going to impose additional structure and smoothness conditions for our econometric procedure.

Our aim is to study the empirical copula of volatility, which begins with a notion of the empirical distribution function of volatility, defined as

equation[equation omitted — 155 chars of source]

where $0 < x,y < \infty$, while $V_{t}^{X} = ( \sigma_{t}^{X})^{2}$ and $V_{t}^{Y} = ( \sigma_{t}^Y)^{2}$ are point-in-time variances of $X$ and $Y$, and $\mathbbm{1}_{ \left\{ \cdot \right\}}$ is the indicator function. $H_{T}(x,y)$ measures the relative amount of time in $[0,T]$ that the stochastic volatility processes spent at different levels over their support. Hence, (ref) is the pathwise counterpart of the cumulative distribution function of volatility. Without the $T$ normalization, $H_{T}(x,y)$ is a multivariate extension of the volatility occupation measure studied in li-todorov-tauchen:13a, which corresponds to the marginal empirical distribution function $F_{T}(x) = H_{T}(x, \infty)$ and $G_{T}(y) = H_{T}( \infty, y)$.

The empirical copula of volatility is motivated by Sklar's theorem. In particular, since $H_{T}(x,y)$ in (ref) is the pathwise analogue of the distribution function of the latent variance process $(V_{t}^{X}, V_{t}^{Y})$ over $[0, T]$, it is natural to separate its marginal features from its dependence structure. We therefore define the empirical copula of volatility as

equation[equation omitted — 182 chars of source]

where $0 < u, v < 1$.

If $F_{T}(x)$ and $G_{T}(y)$ are invertible (e.g., continuous and strictly increasing), (ref) can be written as $C_{T}(u,v) = H_{T}(F_{T}^{-1}(u),G_{T}^{-1}(v))$, which leads to the “inversion formula” $H_{T}(x,y) = C_{T}(F_{T}(x),G_{T}(y))$. Thus, the empirical distribution function of volatility can be expressed as the marginal empirical distribution of the individual volatility processes and a copula that captures their dependence structure up to $T$. $C_{T}(u,v)$ delivers a margin-free characterization of volatility codependency that is invariant to strictly monotone transformations of the marginal volatility processes. This is useful in financial applications, where the marginal behavior of volatility can differ markedly across assets, while their dependence structure can be decided separately. In addition, and in contrast to Pearson's measure of linear correlation, the copula can accommodate nonlinear dependence and tail concentration.

Without prior knowledge of the volatility processes, however, we define the empirical quantile function of $F_{T}(x)$ and $G_{T}(y)$ as the generalized inverse:

equation[equation omitted — 226 chars of source]

with the convention $\inf \emptyset = \infty$.\footnote{The generalized inverse of a distribution function is non-decreasing, left-continuous, and admits a limit from the right. If the distribution function is discontinuous, the quantile function has flat spots, and vice versa.} In the general setting, we still set $C_{T}(u,v) = H_{T}(Q_{T}^{X}(u),Q_{T}^{Y}(v))$, but it should now be given an infimum-based interpretation. Hence, $C_{T}(u,v)$ measures the proportion of time that the volatility processes spend below the smallest coordinate $x$ and $y$ for which $F_{T}(x)$ and $G_{T}(y)$ exceed $u$ and $v$.

To the best of our knowledge, both the multivariate extension of the empirical distribution function of volatility, and the associated copula, are novel to the literature.

High-frequency estimation

In practice, $V^{X}$ and $V^{Y}$ are not directly observed, so neither is the empirical distribution function of volatility nor the copula transformation. Our empirical strategy is based on the premise that we can approximate the spot variance processes with a realized measure constructed from discretely observed high-frequency data of $X$ and $Y$ over small time intervals. In the asymptotic theory, we are going to both increase the number of blocks and shrink the timespan of each block, while padding it with an increasing number of observations, such that we attain a better and better assessment of the volatility's path. In doing so, we control for the price jump component. We also deal with the drift, but this is easier since it is asymptotically negligible in our setting.\footnote{christensen-oomen-reno:22a propose a drift burst model, in which the drift term is locally of larger order than the volatility component.}

We assume that $(X,Y)$ are sampled at equidistant time points $(i \Delta_{n})_{i=0}^{n}$, where $\Delta_{n}$ is the time gap and $n = \lfloor T/ \Delta_{n} \rfloor$ denotes the number of observations over $[0, T]$. We define an increment of the process $Z$ between $(i-1) \Delta_{n}$ and $i \Delta_{n}$ as $\Delta_{i}^{n}Z = Z_{i \Delta_{n}} - Z_{(i-1) \Delta_{n}}$, for $i = 1, \dots, n$. Then, as in li-todorov-tauchen:13a, we adopt the following estimator of $V_{t}^{Z}$ (at time $t$):

equation[equation omitted — 744 chars of source]

$\hat{V}_{t}^{Z}$ is the spot realized variance jacod-protter:12a. $k_{n}$ is the number of log-price increments included in the estimation on each subinterval and tends to infinity in the asymptotic analysis, whereas $\alpha > 0$ and $\varpi \in (0, 1/2)$ relate to the jump truncation device of mancini:09a. The threshold $\alpha \Delta_{n}^{ \varpi}$ is designed to remove log-price increments perturbed by the price jump component and is asymptotically shrinking, slowly enough, to zero. The rate conditions of the tuning parameters are made explicit below.\footnote{Note that in much of the following, we suppress the explicit dependence of various statistics on $n$, when it can be avoided without causing confusion.}

We define the realized distribution function of volatility, an estimator of the empirical distribution function of volatility, as follows:

equation[equation omitted — 159 chars of source]

The univariate version can be retrieved as $\widehat{F}_{n,T}(x) = \widehat{H}_{n,T}(x, \infty) = \frac{1}{T} \int_{0}^{T} \mathbbm{1}_{ \left \{ \hat{V}_{t}^{X} \leq x \right\}} \mathrm{d}t$ and $\widehat{G}_{n,T}(y) = \widehat{H}_{n,T}(\infty,y) = \frac{1}{T} \int_{0}^{T} \mathbbm{1}_{ \left\{ \hat{V}_{t}^{Y} \leq y \right\}} \mathrm{d}t$. The latter are again equivalent to those proposed in li-todorov-tauchen:13a.

Then, the realized copula of volatility is the statistic:

equation[equation omitted — 237 chars of source]

and the realized quantile functions of volatility are extended as follows:

equation[equation omitted — 198 chars of source]

such that $\widehat{C}_{n,T}(u,v) = \widehat{H}_{n,T}( \widehat{Q}_{n,T}^{X}(u), \widehat{Q}_{n,T}^{Y}(v))$.

Asymptotic theory

To derive the asymptotic theory for the realized distribution function of volatility and realized copula of volatility, we prove a consistency result in distinct settings with i) $\Delta_{n} \rightarrow 0$ and $T$ fixed, and ii) $\Delta_{n} \rightarrow 0$ and $T \rightarrow \infty$. In the first, the realized statistic is a natural estimator of the empirical counterpart. In the second, the realized statistic converges to the invariant distribution, assuming it exists. Under further smoothness conditions and weak dependence of the volatility process, we establish a central limit theorem. The latter is done in the functional sense of weak convergence of the probability measures associated with the empirical process of the normalized measurement error of volatility.

To begin with, we introduce a number of additional conditions. The first assumption concerns the regularity of the price jump component.

asuThe process $Z$ follows (ref) with $b^{Z}$ locally bounded and $\sigma^{Z}$ c\`{a}dl\`{a}g. Moreover, for a constant $r \in [0,2]$, and a sequence of stopping times $( \tau_{m})_{m=1}^{ \infty}$, such that $\tau_{m} \rightarrow \infty$ as $m \rightarrow \infty$, there exists a sequence of functions $( \Gamma_{m})_{m=1}^{ \infty}$ such that for each $m: \min(| \delta^{Z}( \omega,t,z)|,1) \leq \Gamma_{m}(z)$ for $( \omega,t,z)$ with $t \leq \tau_{m}( \omega)$, and $\int_{ \mathbb{R}} \Gamma_{m}(z)^{r} \lambda^{Z}( \mathrm{d}z) < \infty$.

The exponent $r$ in Assumption (ref) is an upper bound on the generalized Blumenthal-Getoor index, an adaptation from L\'{e}vy processes to the semimartingale setting, see jacod-protter:12a.\footnote{The Blumenthal-Getoor index of a L\'{e}vy process is defined as $\beta = \inf \left\{ r \geq 0 \;:\; \int_{|x| < 1} |x|^{r} \nu( \mathrm{d}x) < \infty \right\}$, where $\nu$ is the L\'{e}vy measure. $\beta$ can be interpreted as the smallest power at which the distribution of the small jumps (arbitrarily chosen to be of size $|x| < 1$) has finite $r$th moment.} It is therefore related to the jump activity of the log-price process with lower values of $r$ being more binding. The point $r = 1$ separates jump processes with sample paths of finite and infinite variation, i.e., it concerns their absolute summability, while $r = 0$ implies that the jumps are finitely active (on finite time intervals, almost surely). $r$ plays a crucial role in determining rate conditions for the tuning parameters in our estimation procedure, and also in establishing the rate of convergence of our estimator.

asu$H_{T}(x,y)$ is continuous, almost surely.

The continuity of the empirical distribution function imposed by Assumption (ref) is used to control the approximation error in the uniform metric.

thmWe suppose that Assumption (ref) with $r = 2$ and Assumption (ref) hold. As $\Delta_{n} \rightarrow 0$ and $k_{n} \rightarrow \infty$, such that $k_{n} \Delta_{n} \rightarrow 0$, and with $T$ fixed, it follows that \begin{equation} \sup_{(x,y) \in \mathbb{R}_{+}^{2}}| \widehat{H}_{n,T}(x,y) - H_{T}(x,y)| \overset{ \mathbb{P}}{ \longrightarrow} 0. \end{equation} Furthermore, let $\mathcal{Q} = \{(u,v) \in [0,1]^{2} : Q_{T}^{X}(u)$ is continuous at $u$ and $Q_{T}^{Y}(v)$ is continuous at $v$, almost surely$\}$. Then, for each $(u,v) \in \mathcal{Q}$, it further holds that \begin{equation} \widehat{C}_{n,T}(u,v) \overset{ \mathbb{P}}{ \longrightarrow} C_{T}(u,v). \end{equation} Moreover, if $Q_{T}^{X}(u)$ and $Q_{T}^{Y}(v)$ are continuous, almost surely, on $[0,1]^{2}$: \begin{equation} \sup_{(u,v) \in [0,1]^{2}}| \widehat{C}_{n,T}(u,v) - C_{T}(u,v)| \overset{ \mathbb{P}}{ \longrightarrow} 0. \end{equation}

Theorem 1 establishes uniform convergence in probability of $\widehat{H}_{n,T}(x,y)$ as an estimator of $H_{T}(x,y)$. As $\Delta_{n} \rightarrow 0$ and $k_{n} \rightarrow \infty$, such that $\Delta_{n} k_{n} \rightarrow 0$, we get an increasing number of error-free estimates of the spot variance, forming a dense subset on $[0,T]$, from which we can construct an arbitrarily accurate estimate of $H_{T}(x,y)$. The exact rate of $k_{n}$ does not influence this result, so long as $\Delta_{n} k_{n} \rightarrow 0$. We can also show pointwise convergence in probability of $\hat{C}_{n,T}(x,y)$ on the continuity set $\mathcal{Q}$. However, in the discontinuous case at a point $(u,v)$ where either $Q_{T}^{X}(u)$ or $Q_{T}^{Y}(v)$ are not continuous, even uniform convergence of $\widehat{H}_{n,T}(x,y)$ does not imply that $Q_{n,T}^{X}(u) \overset{ \mathbb{P}}{ \longrightarrow} Q_{T}^X(u)$ (or $Q_{n,T}^Y(u) \overset{ \mathbb{P}}{ \longrightarrow} Q_{T}^{Y}(v)$), so $\widehat{C}_{n,T}(u,v)$ is generally inconsistent for $\widehat{C}_{T}(u,v)$ at that point. We can strengthen the convergence of $\widehat{C}_{n,T}(x,y)$ to a uniform statement if the limiting empirical quantile function is continuous, almost surely.

The smoothness imposed on $H_{T}(x,y)$ does not require the volatility path itself to be continuous. Even if $V_{t}^{X}$ or $V_{t}^{Y}$ exhibit jumps, $H_{T}(x,y)$ can be continuous provided the sample path variation is rich enough to fill out the gaps. The intuition is that $H_{T}(x,y)$ aggregates the amount of time that the processes spend below a given threshold, and it can thus “iron out” local irregularities through temporal aggregation.

To establish the asymptotic theory for the combined in-fill and long-span setting, we need stricter control of the error embedded in the recovery of the volatility path. First, we need to restrict the c\`{a}dl\`{a}g dynamic of the volatility processes.

asuThe volatility process $\sigma^{Z}$ is of the form: \begin{equation} \sigma_{t}^{Z} = \sigma_{0}^{Z} + \int_{0}^{t} \tilde{b}_{s}^{Z} \mathrm{d}s + \int_{0}^{t} \tilde{ \sigma}_{s}^{Z} \mathrm{d}\tilde{W}_{s}^{Z} + \int_{0}^{t} \int_{ \mathbb{R}} \tilde{ \delta}^{Z}(s,z)( \tilde{ \mu}^{Z} - \tilde{ \nu}^{Z})( \mathrm{d}s, \mathrm{d}z). \end{equation} Here, $\tilde{b}^{Z}$ and $\tilde{ \sigma}^{Z}$ are adapted and locally bounded processes, $\tilde{W}^{Z}$ is a standard Brownian motion (that may depend on $W^{Z}$), while $\tilde{ \delta}^{Z}$ is a predictable function. Moreover, for a constant $\tilde{r} \in [0,2]$, and a sequence of stopping times $(\tilde{ \tau}_{m})_{m=1}^{ \infty}$ such that $\tilde{ \tau}_{m} \rightarrow \infty$ as $m \rightarrow \infty$, there exists a sequence of functions $( \tilde{ \Gamma}_{m})_{m=1}^{ \infty} : \min \{| \tilde{ \delta}^{Z}( \omega,t,z)|,1 \} \leq \tilde{ \Gamma}_{m}(z)$ for $( \omega,t,z)$ with $t \leq \tilde{ \tau}_{m}( \omega)$, and $\int_{ \mathbb{R}} \tilde{ \Gamma}_{m}(z)^{ \tilde{r}} \tilde{ \lambda}^{Z}( \mathrm{d}z) < \infty$.

Assumption (ref) states that the volatility is an It\^{o} semimartingale. In contrast to the log-price process, this is not a consequence of the no-arbitrage principle. However, it is commonly assumed, because it restricts the local behavior of volatility that is necessary to apply standard estimates for semimartingales in the high-frequency paradigm jacod-protter:12a. For example, if $\sigma^{Z}$ is a diffusion, the assumption follows by an application of It\^{o}'s Lemma under certain smoothness conditions. Either way, it is fulfilled for most stochastic volatility models applied in practice, even when volatility has sources of randomness not captured by $Z$. The ability to correlate the Brownian motions of $Z$ and $\sigma^{Z}$ is particularly relevant for financial applications that involve equity data, in order to capture the so-called leverage effect black:76a, christie:82a. Note that (ref) does rule out volatility models driven by a fractional Brownian motion, which feature prominently in some strands of research comte-renault:98a, gatheral-jaisson-rosenbaum:18a. The main difficulty is that in the “rough” regime, the volatility behaves very erratically over short time intervals, and this makes it harder to control the discretization error chong-todorov:23a. christensen-thyrsgaard-veliyev:19a allow $\sigma^{Z}$ to satisfy a H\"{o}lder-type condition, in the expectation sense. It may be possible to extend our theoretical framework in this direction as well, but we leave it for a future endeavor.

Second, we restrict the memory of volatility processes. To accomplish this, we introduce the additional notation for a generic process $Z$:

equation[equation omitted — 163 chars of source]

Here, $\mathcal{F}_{s} = \sigma(Z_{u}; u \leq s)$ and $\mathcal{F}^{t} = \sigma(Z_{u}; u \geq t)$ are the backward- and forward-looking $\sigma$-algebras generated by $Z$ until time $s$ and from time $t$, respectively. Moreover, $\alpha_{t}$ is the mixing coefficient function that measures the degree of stochastic (in)dependence between them. This is used in the definition of strong mixing introduced by rosenblatt:56a.

definition$Z$ is strong mixing (or $\alpha$-mixing) if $\alpha_{t} \rightarrow 0$ as $t \rightarrow \infty$.
asuThe volatility process $(V_{t}^{X}, V_{t}^{Y})_{t \geq 0}$ is stationary and $\alpha$-mixing in the sense of Definition (ref) with $\alpha_{t} = O(t^{-(1+ \tau)})$ for some $\tau > 0$.

Assumption (ref) is essential in deriving the asymptotic distribution theory. It corresponds to Assumption 5 in christensen-thyrsgaard-veliyev:19a, but for the two-dimensional case, and is again a weak regularity condition that can be verified for most stochastic volatility models. The stationary condition makes it meaningful to speak of a time-invariant distribution function of volatility, $H(x,y) = \mathbb{P}(V_{0}^{X} \leq x, V_{0}^{Y} \leq y)$, from which $F(x) = \mathbb{P}(V_{0}^{X} \leq x)$ and $G(y) = \mathbb{P}(V_{0}^{Y} \leq y)$ can be extracted. Sklar's theorem then guarantees existence, and uniqueness under further conditions, of the copula $C(u,v) = \mathbb{P}(F(V_{0}^{X}) \leq u, G(V_{0}^{Y}) \leq v)$. The question is then to what extent $\widehat{H}_{n,T}(x,y)$, consistently estimating $H_{T}(x,y)$ as $\Delta_{n} \rightarrow 0$, also provides a good proxy of the stationary distribution $H(x,y)$, since there is only a single realization to work with. This, at a bare minimum, additionally requires $T \rightarrow \infty$, so sufficient conditions to apply the ergodic theorem on the volatility path are needed to prompt the law of large numbers for consistency. This is implied by the strong mixing condition, and the polynomial decay of $\alpha_{t}$ ensures that volatility persistence decays suitably fast to derive the limiting distribution of the statistic. In this regard, we should note that since our CLT is derived for a sequence of bounded random variables, all moments exist, so the usual moment condition “trade-off” in Davydov’s inequality is irrelevant, and hence the requirement for the mixing coefficient function is the weakest possible and the autocorrelation function merely needs to be integrable.

A prominent counterexample (ruled out by Assumption (ref)), where Assumption (ref) typically fails, is long-memory processes, such as the fractional Gaussian noise with Hurst exponent $H \in (0.5,1.0)$, which is stationary and weak mixing, but not strong mixing.

The next assumption imposes a uniform boundedness on the coefficient processes.

asuFor all $p \geq 1$, \begin{equation} \sup_{t \geq 0} \left\{ \mathbb{E} \left[|b_{t}^{Z}|^{p} + | \sigma_{t}^{Z}|^{p} + | \tilde{b}_{t}^{Z}|^{p} + | \tilde{ \sigma}_{t}^{Z}|^{p} \right] + \int_{ \mathbb{R}}(1 \wedge \Gamma(z))^{p} \lambda^{Z}( \mathrm{d}z) \right \} \leq C, \end{equation} for a constant $C$.

This follows Assumption 1 in christensen-thyrsgaard-veliyev:19a and andersen-thyrsgaard-todorov:19a. It imposes a light-tailedness condition on the distribution of the various coefficient processes driving the continuous part of $X$ and $Y$, and the large price jumps, but not for the jumps of volatility. It is a common assumption in the long-span setting. The existence of moments of any order is more than we require, but it streamlines the exposition. The technical reason for this condition is that the standard estimates for the increments of semimartingales based on the localization procedure of jacod-protter:12a presumes $T$ is fixed. The above assumption facilitates the derivation of such estimates over $[0, \infty)$.

asu$H_{T}(x,y)$ is continuously differentiable on $\mathbb{R}_{+}^{2}$, almost surely. The partial derivatives (that exist, almost surely) are denoted by $\partial_{x} H_{T}(x,y)$ and $\partial_{y} H_{T}(x,y)$, respectively, and we assume that $\partial_{ \bullet}H_{T}(x,y) \neq 0$ such that $\mathbb{E} \left[ \partial_{ \bullet} H_{T}(x,y) \right] \leq C$ for almost every $(x,y) \in \mathbb{R}_{+}^{2}$, where $C$ is a bounding constant that does not depend on $(x,y)$.

This combines Assumption 2 of christensen-thyrsgaard-veliyev:19a with Assumption B of li-todorov-tauchen:13a. The smoothness condition is used to derive a uniform upper bound for estimation error of the empirical copula function. The requirement $\partial_{ \bullet}H_{T}(x,y) \neq 0$ rules out “flat spots” in the marginal distribution to ensure that $H_{T}(x,y) \mapsto {C}_T(u,v)$ is smooth and non-singular, so that we can apply the functional delta method without getting a degenerate limit process with zero variance. We note that the boundedness in expectation is weaker than asking $\mathbb{E}[ \sup_{(x,y) \in \mathbb{R}_{+}^{2}} \partial_{ \bullet}H_{T}(x,y)] \leq C$, since the latter requires almost sure control of the pathwise behavior of volatility.

At last, we state an extra condition that is required to establish the central limit theorem. To accomplish this, we introduce the notation $k_{n} = \Delta_{n}^{- \gamma}$ and for any $\iota >0$

eqnarray[eqnarray omitted — 385 chars of source]
asuWe assume that $\frac{r-1}{r} < \varpi < \frac{1}{2}$ and $r(1- \varpi) < \gamma <1$ for $r>1$. In addition, $0 < \varpi < \frac{1}{2}$ and $1-(2-r) \varpi < \gamma <1$ for $r \leq 1$. Moreover, $\frac{1- \gamma}{1+ \tilde{r}} > 0$ for $\tilde{r}>0$. Finally, we assume that either of the following conditions holds for any $\iota > 0$: \begin{equation} 1). \:\: r\leq 1 and \sqrt{T}d_{n} \rightarrow 0 \quad or \quad 2). \:\: r > 1 and \sqrt{T}d'_{n} \rightarrow 0. \end{equation}

This is identical to the conditions in Theorem 4 of li-todorov-tauchen:13a, except for the appearance of $\sqrt{T}$ since we operate in the long-span setting. It implies that the volatility error, $\hat{V}_{t}^{Z} - V_{t}^{Z}$, is dominated by $T$, so the high-frequency discretization part is asymptotically negligible. The four terms in $d_n$ and $d_n'$ represent the various sources of the estimation error in $H_{n,T}(x,y)$. Specifically, the first term captures the order of the truncation procedure, induced by replacing truncated increments by those of the continuous part. It depends on whether the jump component has finite or infinite variation. The second and third terms capture the sampling variability and the discretization bias in approximating the spot volatilities. The last order is for the collection of summands over the sampling intervals that contain jumps (with a size larger than some level) in the volatility. The choices of $\varpi, r$, and $\gamma$ imply $d_{n} \rightarrow 0$ and $d_{n}' \rightarrow 0$. Thus, we can assume $T \rightarrow \infty$ such that Assumption (ref) holds. In particular, for $r \in(0, \frac{1}{2})$ and $\tilde{r}=0$, we can take $\varpi \in( \frac{3}{8-4r}, \frac{1}{2})$ and $\gamma=1/2$, from which it follows that $d_{n}=O( \Delta_n^{ \frac{1}{4}- \iota})$. In this case, $T = O( \Delta_{n}^{-1/4+ \epsilon})$ with $\epsilon > 0$ verifies the requirement.

The next theorem delivers the functional CLT of the estimation procedure for the stationary distribution function of volatility and the associated copula. It is grounded in a double asymptotic setting with $\Delta_{n} \rightarrow 0$ and $T \rightarrow \infty$ and extends Theorem 1 in li-todorov-tauchen:13a and Theorem 3.1 in christensen-thyrsgaard-veliyev:19a.

thmSuppose that Assumptions (ref) - (ref) hold (with $r = 2$ in Assumption (ref) and $\tilde{r} = 2$ in Assumption (ref)). As $\Delta_{n} \rightarrow 0$ and $T \rightarrow \infty$, it follows that \begin{equation} \sqrt{T} \left( \widehat{H}_{n,T}(x,y) - H(x,y) \right) \Rightarrow \mathcal{G}. \end{equation} Here, “$\Rightarrow$” denotes weak convergence in the space $\mathbb{D}( \mathbb{R}_{+}^{2})$ of c\`{a}dl\`{a}g functions equipped with the uniform topology. Also, $\mathcal{G}$ is a Brownian bridge on $\mathbb{R}_{+}^{2}$ with covariance function \begin{equation} \mathrm{cov} \big( \mathcal{G}(x,y), \mathcal{G}(x',y') \big) \equiv \mathrm{avar}_{ \mathcal{G}}(x,y,x',y') = 2 \int_{0}^{ \infty} \Big( H_{t}(x,y,x',y') - H(x,y)H(x',y') \Big) \mathrm{d}t, \end{equation} where $H_{t}(x,y,x',y') = P(V_{0}^{X} \leq x, V_{0}^{Y} \leq y, V_{t}^{X} \leq x', V_{t}^{Y} \leq y')$. Moreover, it follows that \begin{equation} \sqrt{T} \left( \widehat{C}_{n,T}(u,v) - C(u,v) \right) \Rightarrow \mathcal{C}, \end{equation} where $\mathcal{C}$ is a stochastic process on $[0,1]^{2}$, derived from $\mathcal{G}$, that is defined as: \begin{equation} \mathcal{C}(u,v) = \mathbf{g}(u,v)^{ \top} \mathbf{Z}(u,v). \end{equation} Here, \begin{equation} \begin{aligned} \mathbf{Z}(u,v) &= \left[ \mathcal{G} \big(Q^{X}(u),Q^{Y}(v) \big), \mathcal{G} \big(Q^{X}(u), \infty \big), \mathcal{G} \big( \infty,Q^{Y}(v) \big) \right]^{ \top}, \\[0.25cm] \mathbf{g}(u,v) &= \left[1, -\partial_{u}C(u,v), -\partial_{v}C(u,v) \right]^{ \top}, \end{aligned} \end{equation} while $Q^{X}(u) = \lim_{T \rightarrow \infty}Q_{T}^{X}(u)$ and $Q^{Y}(v) = \lim_{T \rightarrow \infty}Q_{T}^{Y}(v)$.

The proof of Theorem (ref) draws on the functional CLT for stationary mixing sequences of bounded random variables (see Theorem 18.5.4 in ibragimov:75a and Theorem 5.2 in dehay:05a) building on rosenblatt:56a. We reiterate that the rate of convergence and the asymptotic distribution of $\widehat{H}_{n,T}(x,y)$ and $\widehat{C}_{n,T}(x,y)$ are unaffected by the error from recovering volatility, $V_{t}^{Z} - \widehat{V}_{t}^{Z}$, which is asymptotically negligible under the rate conditions in (ref).

The covariance function of $\mathcal{C}(u,v)$ follows from Theorem (ref) and the functional delta rule, i.e.

equation[equation omitted — 208 chars of source]

where

equation[equation omitted — 117 chars of source]

As a consequence,

equation[equation omitted — 349 chars of source]

where $\overset{d}{ \longrightarrow}$ is regular convergence in law. This result can be used to construct pointwise confidence intervals, or conduct hypothesis tests, once an estimator of the asymptotic covariance function has been designed.

To construct such an estimator, we set:

equation[equation omitted — 215 chars of source]

with $\widehat{H}_{t,n,T}(x,y,x',y') = \frac{1}{T-T^{ \xi}} \int_{0}^{T-T^{ \xi}} \mathbbm{1}_{ \left\{ \hat{V}_{s}^{X} \leq x, \hat{V}_{s}^{Y} \leq y, \hat{V}_{s+t}^{X} \leq x', \hat{V}_{s+t}^{Y} \leq y' \right\}} \mathrm{d}s$, where $\xi$ is a small positive constant.

corollaryWe suppose that Assumptions (ref) - (ref) hold, and that $\xi \in (0, 1/3)$. As $\Delta_{n} \rightarrow 0$ and $T \rightarrow \infty$, we deduce that \begin{equation} \widehat{ \mathrm{avar}}_{ \mathcal{G}}(x,y,x',y') \overset{ \mathbb{P}}{ \longrightarrow} \mathrm{avar}_{ \mathcal{G}}(x,y,x',y'). \end{equation}

In the above, $T^{ \xi}$ plays the role of the “lag length” in the estimation of the long-run variance of a stationary time series, which should not grow too fast for a consistent estimation. In fact, the optimal order of the bandwidth is often $T^{1/3}$ newey-west:94a. The choice of $\xi$ relates to this observation.

By Slutsky's theorem,

equation[equation omitted — 176 chars of source]

Conducting inference about $C(u,v)$ is a more complicated endeavor, since it also involves partial derivatives of $C(u,v)$ on top of the covariance function $\mathrm{avar}_{ \mathcal{G}}(x,y,x',y')$. To construct estimators of the former, we propose a nonparametric kernel-based smoothing approach:\footnote{This idea was introduced in rosenblatt:56b.}

equation[equation omitted — 466 chars of source]

Here, $K$ is defined on a compact subset of $\mathbb{R}^{2}$, Lipschitz continuous, and such that

equation[equation omitted — 171 chars of source]

This is a common regularity condition for kernel density estimation. Note that $K$ is not required to be symmetric.

To show consistency of $\widehat{ \partial_{u}C}(u,v)$ and $\widehat{ \partial_{v}C}(u,v)$, we take a slight detour by forming a preliminary nonparametric estimator of the copula density, i.e. $c(u,v) = \partial_{u} \partial_{v} C(u,v)$,

equation[equation omitted — 232 chars of source]

The following condition is sufficient to get consistency of $\hat{c}_{n,T}(u,v)$.

asu$H_{T}(x,y)$ is three times continuously differentiable on $\mathbb{R}_{+}^{2}$.

In other words, we implicitly assume that the empirical copula density $c_{T}(u,v)$ exists and is itself continuously differentiable on $[0,1]^{2}$.

corollaryWe suppose that Assumptions (ref) - (ref) hold. As $\Delta_{n} \rightarrow 0$, $T \rightarrow \infty$, $h \rightarrow 0$, such that \begin{equation} \quad \frac{d_{n}}{h^{3}} \vee \frac{1}{ \sqrt{T}h^{3}} \rightarrow 0, \end{equation} where $d_{n}$ is defined in (ref), for any $(u, v) \in [0,1]^{2}$ it holds that \begin{eqnarray} T} \hat{c}_{n,T}(u, v) \overset{ \mathbb{P}}{ \longrightarrow} c(u,v). \end{eqnarray}

The partial derivatives of $C(u,v)$ can be written as $\partial_{u}C(u,v) = \int_{0}^{v} c(u,w) \mathrm{d}w$ and $\partial_{v} C(u,v) = \int_{0}^{u} c(w,v) \mathrm{d}w$ by Assumption (ref). Hence, $\widehat{ \partial_{u} C}(u,v)$ and $\widehat{ \partial_{v} C}(u,v)$ are Riemann approximations of these integrals based on the kernel density estimator $\hat{c}_{n,T}(u,v)$. Therefore, the consistency of the partial derivatives in (ref) follows directly, namely

equation[equation omitted — 246 chars of source]

In view of (ref) and (ref), we define the following plug-in estimator for the covariance functional of the realized copula of volatility, $\mathrm{avar}_{ \mathcal{C}}(u,v,u',v')$:

equation[equation omitted — 166 chars of source]

where

equation[equation omitted — 269 chars of source]

and

equation[equation omitted — 211 chars of source]

is a three-dimensional, mean zero Gaussian vector whose covariance is estimated via the long-run covariance functional in (ref). Then, by dominance of convergence in probability,

equation[equation omitted — 176 chars of source]

Uniform confidence bands

The pointwise inference derived above ensures asymptotically correct coverage and testing size for each fixed value of $(u,v)$, or $(x,y)$, but it conceals that the functional CLT established in Theorem (ref) offers stronger control of the measurement errors. In this subsection, we exploit this observation to construct uniform confidence bands for the volatility copula.

corollaryUnder the assumptions of Theorem (ref): \begin{equation} \sup_{(u,v) \in[0,1]^{2}} \sqrt{T} \left( \widehat{C}_{n,T}(u,v) - C(u,v) \right) \overset{d}{ \longrightarrow }\sup_{(u,v) \in[0,1]^{2}}{\mathcal{C}(u,v)}, \end{equation} where $\mathcal{C}$ is the Gaussian process defined in Theorem (ref).

Corollary (ref) employs the continuous mapping theorem for functional convergence in law, since the supremum functional is continuous in the uniform metric. However, the limiting distribution of the supremum point is non-standard, so we propose a simulation procedure to construct valid confidence bands.

We define the $(1- \alpha)$-quantile of the right-hand side in (ref) as follows:

equation[equation omitted — 182 chars of source]

Thus, it follows that as $T \rightarrow \infty$,

equation[equation omitted — 177 chars of source]

so an asymptotic $(1- \alpha)$ uniform confidence band for $C(u,v)$ is given by

equation[equation omitted — 154 chars of source]

To compute the constant $q_{1- \alpha}^{CB}$, we follow a Monte Carlo approach, which builds on a consistent estimator of the covariance function of $\mathcal{C}(u,v)$. To implement it, we discretize $[0,1]^{2}$ into an $m \times m$ grid $\left \{(u_{i}, v_{j}) \right\}_{i,j=1}^{m}$ and form the corresponding $m^{2} \times m^{2}$ covariance matrix, $\widehat{ \mathrm{avar}}_{ \mathcal{C}}$. Based on the Cholesky decomposition of $\widehat{ \mathrm{avar}}_{ \mathcal{C}} = LL^{ \top}$, we repeat the following procedure for $b = 1, \ldots, B$:

enumerate• Draw $z^{(b)} \overset{\text{i.i.d.}}{ \sim} N(0, I_{m^{2}})$. • Set $\mathcal{Z}^{(b)} = Lz^{(b)}$, such that $\mathcal{Z}^{(b)} \sim N(0, \widehat{ \mathrm{avar}}_{ \mathcal{C}})$. • Compute $S_{(b)} = \max_{1 \leq k \leq m^{2}} \bigl| \mathcal{Z}^{(b)}_{k} \bigr|$. • Retrieve the $(1- \alpha)$-quantile of $\left \{S_{(1)}, \ldots, S_{(B)} \right\}$: \begin{equation} \widehat{q}_{1- \alpha}^{CB} = S_{(\lceil(1- \alpha)B \rceil)}. \end{equation} • Construct the feasible uniform confidence band: \begin{equation} \left\{ \left[ \widehat{C}_{n,T}(u,v) \pm \frac{ \widehat{q}_{1- \alpha}^{CB}}{ \sqrt{T}} \right], \quad for all (u,v) \in [0,1]^{2} \right\}. \end{equation}

The validity of this resampling procedure is formally established by the next proposition, which is a functional version of Slutsky's theorem.

proWe assume the conditions of Theorem (ref) and Corollary (ref) - (ref) hold, such that \begin{equation} \sqrt{T} \left( \widehat{C}_{n,T}(u,v) - C(u,v) \right) \Rightarrow \mathcal{C}, \end{equation} where $\mathcal{C}$ is a centered Gaussian process with continuous covariance function $\mathrm{avar}_{ \mathcal{C}}(u,v,u,v)$ bounded away from zero, i.e. $\inf_{(u,v) \in [0,1]^{2}} \mathrm{avar}_{ \mathcal{C}}(u,v,u,v) > 0$. Also, let $\widehat{ \mathrm{avar}}_{ \mathcal{C}}$ be the plug-in estimator defined in (ref), which consistently estimates the covariance function of $\mathcal{C}$ at the grid points. Conditionally on $\widehat{ \mathrm{avar}}_{ \mathcal{C}}$, let $( \widehat{ \mathcal{C}}_{n,T})$ be a centered Gaussian processes with covariance matrix $\widehat{ \mathrm{avar}}_{ \mathcal{C}}$. Write \begin{equation} \widehat{q}_{1- \alpha}^{CB} = \inf \left\{c>0: \mathbb{P} \left( \sup_{i,j} \big| \widehat{ \mathcal{C}}_{n,T}(u_{i},v_{j}) \big| \leq c \mid \widehat{ \mathrm{avar}}_{ \mathcal{C}} \right) \geq 1-\alpha \right\}, \end{equation} Then, as the grid is made arbitrarily fine, i.e. $m \rightarrow \infty$, with $\Delta_{n} \rightarrow 0$ and $T \rightarrow \infty$, for any constant $c \geq 0$, we deduce that \begin{equation} \mathbb{P} \left( \sup_{(u, v) \in [0,1]^{2}} \big| \widehat{ \mathcal{C}}_{n,T}(u, v) \big| \leq c \mid \widehat{ \mathrm{avar}}_{ \mathcal{C}} \right) \overset{ \mathbb{P}}{ \longrightarrow} \mathbb{P} \left( \sup_{(u,v) \in [0,1]^{2}} \left| \mathcal{C}(u,v) \right| \leq c \right), \end{equation}

It follows from the proposition that $\widehat{q}_{1-\alpha}^{CB} \overset{ \mathbb{P}}{ \longrightarrow} q_{1-\alpha}^{CB}$. Therefore, the uniform confidence band constructed from $\widehat{q}_{1- \alpha}^{CB}$ has the correct asymptotic coverage:

equation[equation omitted — 189 chars of source]

Goodness-of-fit testing

The previous analysis can also be extended to construct a goodness-of-fit test of the volatility copula function, which is another prominent application of our functional convergence results. In particular, we introduce a null hypothesis and an alternative hypothesis:

equation[equation omitted — 124 chars of source]

Here, equality is almost everywhere on $[0,1]^{2}$.

To construct a testing procedure, we let

equation[equation omitted — 96 chars of source]

We now arrive at the following statement, which again follows from the continuous mapping theorem for functional convergence applied to Theorem (ref).

corollaryUnder the assumptions of Theorem (ref) and conditional on $\mathcal{H}_{0}$: \begin{equation} \lVert G_{n,T} \rVert_{L_{2}}^{2} \overset{d}{ \longrightarrow} \lVert \mathcal{C} \rVert_{L_2}^{2}. \end{equation}

The limiting distribution in Corollary (ref) (a sort of weighted infinite mixture of $\chi^{2}$ random variables) is again not standard, so we employ a Monte Carlo procedure, closely related to the one outlined above, to form a test of $\mathcal{H}_{0}$ versus $\mathcal{H}_{1}$. In particular, we let $q_{1- \alpha}^{GF}$ be a consistent estimator of the $(1- \alpha)$-quantile of the distribution of $\lVert \mathcal{C} \rVert_{L_{2}}^{2}$, and we define a rejection region:

equation[equation omitted — 105 chars of source]

Then, in view of Corollary (ref),

equation[equation omitted — 162 chars of source]

To obtain $q_{1-\alpha}^{GF}$, we approximate the distribution of the random variable $\lVert \mathcal{C} \rVert_{L_{2}}^{2}$ by a weighted Riemann sum on a, possibly unequally spaced, grid $\{(u_{i},v_{j}) \}_{i,j=1}^{m} \subset [0,1]^{2}$, i.e. $\sum_{1 \leq i,j \leq m} w_{ij} \mathcal{C}(u_{i},v_{j})^{2}$, where $w_{ij}$ are weights induced by the grid and $\sum_{i,j} w_{ij} = 1$. Note that $\{ \mathcal{C}(u_{i}, v_{j}) \}_{i,j=1}^{m}$ is a sequence of $m^{2}$ normal random variates with explicit covariance matrix given by (ref) following Theorem (ref), which we denote by $\Gamma$. Therefore,

equation[equation omitted — 146 chars of source]

where $\chi_{k}^{2}$ are i.i.d. $\chi^{2}(1)$ random variables and $\pi_{k}$ are the eigenvalues of $W^{1/2} \Gamma W^{1/2}$ with $W = \mathrm{diag} \big( \{w_{ij} \}_{i,j=1}^{m} \big)$, for $k = 1, \dots, m^{2}$.

We draw repeated observations of the right-hand side of (ref) based on the eigenvalues from the estimator $W^{1/2} \widehat{ \Gamma} W^{1/2}$ to simulate the distribution of $\lVert G_{n,T} \rVert_{L_{2}}^{2}$ under $\mathcal{H}_{0}$, from which an appropriate quantile can be extracted.\footnote{Since $\widehat{ \Gamma}$ is an estimator of a covariance matrix, it can possess negative eigenvalues. We follow andersen-su-todorov-zhang:24a and retain only those terms in the sum associated with positive eigenvalues.}

To estimate the covariance matrix $\Gamma$, we rearrange the pairs $\{(u_{i}, v_{j})\}_{i,j=1}^{m}$ as $\{(u_{k}', v_{k}')\}_{k=1, \dots, m^{2}}$ and compute $\widehat{ \Gamma}_{k \ell} = \hat{g}_{k}^{ \top} \hat{A}_{k \ell} \hat{g}_{ \ell}$, where

equation[equation omitted — 274 chars of source]

and

equation[equation omitted — 227 chars of source]

It is also possible to test a composite null hypothesis:

eqnarray*[eqnarray* omitted — 116 chars of source]

where $\mathcal{S} = \left\{C_{ \theta} : \theta \in \Theta \subseteq \mathbb{R}^{p} \right\}$ is a statistical model for the copula function, and $\theta$ is a $p$-dimensional parameter vector.

In practice, the true value of the copula parameter, $\theta_{0}$, is unknown. A complication occurs if the parameter is estimated from the full sample, in which case the estimator, $\hat{ \theta}$, is typically $\sqrt{T}$-consistent matching the convergence rate of $\widehat{C}_{n,T}(u,v)$. Then, $G_{n,T}(u,v) = \sqrt{T}( \widehat{C}_{n,T}(u,v) - C_{ \hat{ \theta}}(u,v)) = \sqrt{T}(\widehat{C}_{n,T}(u,v) - C_{ \theta_{0}}(u,v)) - \sqrt{T}( \hat{ \theta} - \theta_{0})^{ \top} \partial_{ \theta} C_{ \theta_{0}}(u,v) + o_{p}(1)$, where $\partial_{ \theta} C_{ \theta_{0}}(u,v) = \partial_{ \theta} C_{ \theta}|_{ \theta = \theta_{0}}$, so $G_{n,T}$ converges in law to a Gaussian process plus a bias term that is a linear functional of $\sqrt{T}( \hat{ \theta} - \theta_{0}) \overset{d}{ \longrightarrow} \mathcal{W} \sim N(0, \Sigma_{ \theta})$. Hence, the $L^{2}$-statistic $\lVert G_{n,T} \rVert_{L_{2}}^{2} \overset{d}{ \longrightarrow} \lVert \mathcal{C} - \langle \partial_{ \theta} C_{ \theta_{0}}, \mathcal{W} \rangle \rVert_{L_2}^{2}$ is a quadratic form of a Gaussian process, which depends on unknown quantities. In other words, the asymptotic distribution is non-pivotal, as it depends on the parameter vector $\theta_{0}$, the gradient $\partial_{ \theta} C_{ \theta_{0}}(u,v)$ and the law of $(\mathcal{W}, \mathcal{C})$. However, we can easily resample the empirical process and the plug-in estimator of the covariance function using a parametric or multiplier bootstrap under $C_{\hat{ \theta}}$ to approximate the distribution of $\lVert G_{n,T} \rVert_{L_{2}}^{2}$ and construct critical values without analytically having to remove the dependence on $\theta_{0}$.

On the other hand, to avoid the effect of $\hat{\theta}$, we can estimate the copula parameter such that $\hat{\theta} - \theta_{0}$ achieves a convergence rate faster than $\widehat{C}_{n,T}(u,v) - C_{ \theta_{0}}(u,v)$.\footnote{This can be accomplished, for example, by reserving data in a shorter time interval, say $[0, T^{\eta}]$ with $0 < \eta < 1$, to construct $\widehat{C}_{n,T}$, but still employing data over the whole time interval $[0,T]$ to get $\hat{\theta}$.} In this case, the test statistic is asymptotically pivotal and we can replace $C_{0}$ by $C_{\hat{\theta}}$ in the previous testing procedure with a single null hypothesis and attain a valid inference that does not involve nuisance parameters and is shape-only. Hence, it suffices to prepare a single table of critical values.

Simulation analysis

In this section, we conduct a Monte Carlo study to explore the finite sample properties of our realized estimator of the empirical copula of stochastic volatility. We restrict attention to the bivariate setting and assume that a pair of asset log-price processes follow the motion

equation[equation omitted — 113 chars of source]

where

equation[equation omitted — 152 chars of source]

for $i = 1$ and $2$.

In the above, $W$ and $B$ are Brownian motions with leverage correlation $\text{corr}[\mathrm{d}W_{i,t} \mathrm{d}B_{i,t}] = \rho = -\sqrt{0.5}$. We assume a common speed of mean reversion $\kappa = 10$. The other parameters are $\mu_{1} = 0.10$, $\eta_{1} = \log(0.05)$, $\beta_{1} = 3$, $\mu_{2} = 0.05$, $\eta_{2} = \log(0.01)$, and $\beta_{2} = 3$. This design is intended to represent a risky asset with an expected return and volatility of 10% and 20% and a near risk-free asset with an expected return and volatility of 5% and 10%. It emulates our empirical work in Section (ref), where we look at high-frequency data from stock index and medium-duration treasury bond futures contracts.

A routine calculation for the stationary Gaussian Ornstein-Uhlenbeck process shows that the unconditional distribution of $\log V_{i,t}$ has the form

equation[equation omitted — 93 chars of source]

We denote the distribution function as $F_{ \log V_{i}}(v)$, for $- \infty < v < \infty$.

The univariate log-variance processes are tied together with a Gumbel copula:

equation[equation omitted — 178 chars of source]

which is parameterized by $\theta \in [1, \infty)$. A larger value indicates a stronger relationship. The independence copula---$C(u,v) = uv$---corresponds to $\theta = 1$, while the comonotonicity copula---$C(u,v) = \min(u,v)$---induced by the upper Fr\'{e}chet-Hoeffding bound is the limit as $\theta \rightarrow \infty$. In our empirical application, we explore the properties of this particular copula, because our goodness-of-fit test suggests that it provides an excellent description of the dependence between the log-variance of the equity and treasury bond market. In line with our maximum likelihood estimate from that section, we set $\theta = 2$.

To construct a realization of this process, we begin by drawing from the Gumbel copula to generate a pair of dependent random variables $(U_{1,0}, U_{2,0})$, where $U_{i,0} \sim U(0,1)$. These serve to initialize the log-variance processes via $\log V_{i,0} = F_{\log V_i}^{-1}(U_{i,0})$. We then apply an Euler approximation with constant step size $\Delta$ to iteratively generate the sample path of (ref)--(ref) forward in time, using the Metropolis--Hastings algorithm to ensure that each proposal update is consistent with the target distribution.\footnote{The symbol $\Delta$ does double-duty in this context. In isolation it is the time increment $\Delta = t - s$, but when it is applied to a stochastic process $Y$, it represents the difference operator $\Delta Y_{i,t} = Y_{i,t} - Y_{i,s}$.}

In the “continuous-time” version of the model, we set $\Delta = 1/N$ with $N = 23{,}400$. We assume data is captured over a time interval $[0,T]$ with $T = 50, 100, \dots, 2000$, but the log-price process is only observed at a much coarser equidistant grid $t_{i} = i/n$ with $n = 39, 78, 390$, for $i = 0, 1, \dots, nT$. This design is intended to imitate a security that is traded in a 6.5 hour window each day and whose price is recorded at the ten-, five-, and one-minute frequency over at most an eight-year period. This resembles the dataset analyzed in our empirical application. To estimate the latent spot volatility, associated with each $n$ we employ a bandwidth of $h_{n} = 36, 48, 120$. With the above interpretation, this amounts to a six-, four-, and two-hour window. Hence, as $n$ increases we collect a larger amount of high-frequency data while reducing the time span of the estimation window, as stipulated by the rate conditions in the asymptotic theory. These choices follow li-todorov-tauchen:13a and christensen-thyrsgaard-veliyev:19a.\footnote{In unreported results that are available at request, we conducted a more extensive set of numerical experiments with other bandwidth choices. This showed that the average squared estimation error of the spot volatility estimator was basically flat over an extended region, so long as $h_{n}$ was not too small or too large.}

A sample path of each log-variance process, together with the log-realized variance, is illustrated in Panel A of Figure (ref). In Panel B, we report a scatter of 500 randomly selected pairs of the foregoing series, after transforming them through their empirical distribution functions. We contrast this with selected contours of the probability density function of the Gumbel copula. The latter are increasing as we move toward the 45$^{\circ}$-line and eventually disconnect at the lower-left and upper-right corners. The main takeaway is that while the log-variance observations are dispersed as expected, the estimation error in the log-realized variance shifts several observations into regions where the Gumbel copula has almost no concentration of mass.

figure[figure omitted — 972 chars of source]

The outcome of the full-blown simulation analysis is summarized in Figure (ref). We perform $M = 1{,}000$ repetitions in total. In Panel A, we construct a decile plot of the bivariate distribution function of the Gumbel copula, which we contrast against the empirical copula of volatility---based on $T = 2{,}000$---along with the estimated one recovered by the realized version. We observe that even with a conservative amount of high-frequency data, the empirical copula is close to its stationary counterpart with only some minor deviations. In contrast, the contours of the realized version are generally farther away and located to the northeast of the empirical copula. The explanation is that the sampling variation in the realized variances render them less correlated and shifts mass in the direction of the independence copula.\footnote{The level curves of the independence copula are given by $u = c/v$ for $c \in (0,1)$.} As $n$ increases, however, the realized copula of volatility steadily trends toward the empirical counterpart, which is a validation of the asymptotic theory developed in Section (ref).

figure[figure omitted — 869 chars of source]

In Panel B of Figure (ref), we report the root mean squared error (RMSE):

equation[equation omitted — 152 chars of source]

where $C(u,v)$ is the Gumbel copula given by the expression in (ref) and $C_{m}^{x}(u,v)$ represents either the empirical ($x = \mathrm{Empirical}$) or realized ($x = \mathrm{Realized}$) estimate during the $m$th simulation. The figure reveals that, as expected, the RMSE of the empirical copula of volatility approaches zero, as $T \rightarrow \infty$, which is consistent with our theoretical exploration. It takes about a quadrupling of the long-span dimension to reduce the RMSE in half, which reflects its convergence rate. Turning to the realized copula of volatility, we see that the RMSE is generally higher, and while it too drops as $T$ increases, it tends to flatten out strictly above zero. This is because the estimation error does not vanish with a fixed $n$. Boosting $n$, however, lessens the discretization error, and by the time that $n = 390$ the realized version is within proximity of the empirical copula (its target in the in-fill limit as $n \rightarrow \infty$).

To gauge the functional CLT for the realized copula of volatility as an approximation to its finite sample behavior, we next explore the distribution of the asymptotic pivot:

equation[equation omitted — 181 chars of source]

for $u = v = 0.10, 0.25, 0.50, 0.75, 0.90$.

figure[figure omitted — 1,164 chars of source]

Figure (ref) presents the outcome. We base the analysis on $n = 390$ and $T = 2{,}000$, which adheres to our empirical work in the next section, and the estimator of the asymptotic variance introduced in Section (ref). In Panel A, we show kernel density estimates, while Panel B reports the associated quantile-quantile (Q--Q) plot. We include the sample average and sample standard deviation of $Z_{ \mathcal{C}}$ for each setting. We observe that the kernel density estimates are well-aligned with the standard normal curve, and the Q--Q plots follow the $45^{o}$-line reasonably. The dispersion is almost as predicted, but the test statistic does exhibit a minuscule downward bias of about 0.2-0.5 standard deviation. Importantly, these numbers are more or less identical to those we compute for the empirical copula of volatility.\footnote{The kernel density estimates and Q--Q plots for the empirical copula of volatility are excluded from the presentation for brevity, but they can be shared at request.}

At last, we turn our attention toward the goodness-of-fit test. We calculate rejection rates of the test statistic in (ref) from Corollary (ref) with $n = 390$ and $T = 250,500,1000,2000$. To curb the computational cost, we implement the bootstrapping algorithm for producing critical values with 5,000 draws on a $5 \times 5$ grid $(u,v) \in \{0.10, 0.25, 0.50, 0.75, 0.90 \}^{2}$, yielding 25 evaluation markers.\footnote{It is arguably preferable to use a more refined grid, or possibly even alter how the points are scattered over the domain. The intuition is that if the test statistic gets a better resolution of the realized copula of volatility, it has more scope to detect departures from the null, and this may improve rejection rates under the alternative. However, this also increases the dimensionality of the required covariance matrix, impacting memory usage and runtime speed. In line with this argument, we surmise that the power of the test reported here is conservative.}

The null hypothesis (size) is based on the Gumbel copula with $\theta = 2$, which has Kendall's tau given by $\tau = 1 - 1/ \theta$, so that our choice implies $\tau = 0.5$. The alternative hypothesis (power), is given by three distinct parametric models. The first benchmark is the independence copula, $C(u,v) = uv$. This corresponds to the absence of volatility codependency with $\tau = 0$ and represents a stark violation of $\mathcal{H}_{0}$. It provides a natural baseline for assessing the ability of the test to detect tail concentration. The second choice is the Clayton copula, $C(u,v) = \left(u^{- \alpha} + v^{- \alpha} - 1 \right)^{-1/ \alpha}$ with parameter $\alpha \geq 0$. This is another example of an Archimedean copula used to capture positive tail dependence, but in contrast to the Gumbel copula, the dependence is located in the other extreme of the distribution. The Clayton copula has $\tau = \alpha/( \alpha+2)$, and we set $\alpha = 2$---or $\tau = 0.5$---to be consistent with the Gumbel. This design intends to isolate departures from the null in the form of dependence in the opposite direction, while fixing the overall rank correlation summarized by Kendall's tau. At last, we complement the analysis with a power experiment based on an alternative member of the Gumbel class. Specifically, we examine a Gumbel copula with $\theta = 1.5$---or $\tau = 1/3$---that also has upper tail concentration, albeit slightly weaker compared to the null.

table[table omitted — 1,333 chars of source]

The results are reported in Table (ref). We inspect the nominal significance levels 10%, 5%, and 1%. The size analysis is shown in Panel A. We observe that the test is to a small extent oversized for the lowest value of $T$, but it rapidly settles around the expected level. Moving to Panel B, we learn that the power versus the independence copula is almost unity across scenarios, even at the smallest level of significance. Moreover, the test has a good ability to distinguish between the Gumbel null and the Clayton alternative, especially with larger sample sizes. Lastly, while the rejection rates against the Gumbel(1.5) alternative are initially higher than those of the Clayton copula, but they increase less fast with $T$. The power function is possibly flatter in the Gumbel direction, where the functional form is identical and only the strength of the upper tail dependence is different.

Overall, our simulations suggest that the realized copula of volatility can reliably infer the latent volatility copula at a sampling frequency that is relevant for practical analysis. Moreover, the small sample properties of the standardized statistic are close to those implied by Theorem (ref), which is assuring, because it means that confidence intervals and hypothesis tests based on the normal distribution are accurate. We also uncover that the goodness-of-fit test has appropriate size control and excellent power against many relevant alternatives. A word of caution, though. Our goodness-of-fit test merely provides evidence against the null hypothesis as a whole. Thus, a rejection of a simple null in favor of the alternative should not be interpreted as identifying the exact source of misspecification. Indeed, it can arise either because the assumed functional form of the copula is wrong, or because the family is correctly specified but the concrete parameter value is not.

Empirical application

The forensic analysis presented here aims to determine the codependency between the volatility processes of a pair of leading financial indicators. We delve into high-frequency transaction data from futures contracts that track the aggregate U.S. equity index and treasury bond market, namely the E-mini S&P 500 (ES) and 10-year treasury note (TY). The chosen instruments are listed on the Globex platform at the Chicago Mercantile Exchange (CME). They trade around the clock five days per week, but we restrict attention to the most active hours from 9:30am to 4:00pm EST that overlap with the NYSE trading session. The data at our disposal were purchased from Tickdata (\url{https://www.tickdata.com/}). The sample covers the period from March 18, 2010 to October 14, 2021. We remove short days associated with a reduced trading schedule, which leaves $T = 2{,}891$ days for our empirical investigation.\footnote{To construct a consecutive price series, we employ a built-in procedure in the extraction software delivered by Tickdata. It rolls the front contract over into the back contract about a week before expiration.}

table[table omitted — 886 chars of source]

In Table (ref), we report a few descriptive statistics of the retained data. The “$N$” column reports the average daily number of transactions (in 1000s). These futures contracts are very liquid. However, sampling the tick-by-tick data at the maximal resolution can be detrimental for volatility estimation, because of microstructure noise, such as price discreteness and bid-ask spread hansen-lunde:06b. To alleviate this concern, we downsample to a one-minute equidistant frequency using previous tick imputation, so there are $n = 390$ log-returns available per day to excavate the volatility process.\footnote{We follow the design of the simulation section in terms of choosing the various tuning parameters that are required to calculate the spot realized variance.}

figure[figure omitted — 610 chars of source]

In what follows, we study the log-spot realized variance, which is less influenced by extreme outliers. We remark that by the invariance principle the copula of a random vector does not change under a strictly increasing transformation, so the realized copula of volatility is unaffected by studying the log-transform. In Figure (ref), we show a relative frequency histogram for the log-realized variance of ES in Panel A, whereas Panel B is for TY. In agreement with the sample average and quantiles reported in Table (ref), we observe that volatility in the stock market is, on average, higher than in the bond market. Furthermore, and consistent with andersen-bollerslev-diebold-ebens:01a, andersen-bollerslev-diebold-labys:03a, the shapes are close to the Gaussian bell curve, although both exhibit a mild positive skewness and excess kurtosis compared to the normal distribution. As evident from the right-hand side of Table (ref), this effect is more pronounced in TY than ES.

To provide an initial assessment on the degree of volatility codependency in these time series, Figure (ref) reports the average daily log-spot realized variance for ES in Panel A and TY in Panel B. The smoothing attenuates the measurement error in the spot volatility estimator a bit to better capture any comovement. There is forceful evidence of a positive relationship, most notably during periods of elevated market distress. To support this claim, in Panel A of Figure (ref) we draw a scatter plot of the log-realized variance of ES (on the $x$-axis) against the log-realized variance of TY (on the $y$-axis). The graph reinforces the previous impression, which is further backed by Pearson's and Kendall's sample correlation coefficients.\footnote{Kendall's tau has been transformed as $\hat{ \rho}_{ \text{K}}^{a} = \sin( \pi/2 \hat{ \rho}_{ \text{K}}^{u})$, where $a$ ($u$) denotes the adjusted (unadjusted = raw) estimator. The correction---motivated by Greiner's equality---makes it unbiased for the correlation coefficient of a random sample from a bivariate normal distribution and, hence, its scale is comparable to the Pearson's measure.}

figure[figure omitted — 602 chars of source]

The realized copula of volatility is illustrated in Panel B of Figure (ref). We plot 500 randomly drawn pairwise observations of $\left( \widehat{F}_{n,T} \big( \log( \hat{V}_{t}^{ \text{ES}}) \big), \widehat{G}_{n,T} \big( \log( \hat{V}_{t}^{ \text{TY}}) \big) \right)$ and add decile contours based on the whole sample. We contrast the latter with those implied by a parametric Gumbel copula, where the single parameter was estimated by maximum likelihood to be $\hat{ \theta}_{ \text{MLE}} = 1.4786$. Overall, there is a close alignment between them.

To provide an alternative graphical representation of the realized copula of volatility, we construct a couple of descriptive statistics that are widely used in the literature. The first is an exceedance measure called tail concentration. It studies the amount of probability mass in the lower-left and upper-right quadrant of a bivariate distribution function, here in the form of a copula joe:93a,sibuya:60a. In particular, the lower (L) and upper (U) tail concentration functions are defined as

equation[equation omitted — 196 chars of source]

for any $z \in (0,1)$.\footnote{We note that the tail concentration function is symmetric, since its value does not change if we swap the position of $U$ and $V$, because the latter are uniformly distributed.}

figure[figure omitted — 1,107 chars of source]

At one end of the domain, $L$ is convergent toward $\lim_{z \rightarrow 1^{-}}L(z) = 1$, whereas $R$ has the limit $\lim_{z \rightarrow 0^{+}} R(z) = 1$. This holds for every copula, and so it is not very informative. In the other end, however, $L(z)$ and $R(z)$ can converge to any number in the interval $[0,1]$. Provided the limits exist, the lower and upper tail concentration parameters are defined as $\lambda_{L} = \lim_{z \rightarrow 0^{+}} L(z)$ and $\lambda_{U} = \lim_{z \rightarrow 1^{-}} R(z)$. If either limit is nonzero, we say there is asymptotic dependence in that tail. Furthermore, noting that $L(0.5) = R(0.5)$, it is common to plot $T(z) = \min(L(z),R(z))$ as the aggregate tail concentration function.

The second measure is Kendall's distribution function of a copula, which is defined as

equation[equation omitted — 74 chars of source]

Here, $C(U,V)$ is a univariate random variable, and $(U,V)$ are uniform with bivariate distribution function $C(u,v)$. This can be interpreted as the multivariate version of the probability integral transform genest-rivest:93a.\footnote{The connection between this probability measure and Kendall's tau is $\rho_{ \text{K}} = 3-4 \int_{0}^{1}K(z) \mathrm{d}z$ schweizer-wolff:81a.}

figure[figure omitted — 1,146 chars of source]

We plot the estimated tail concentration function of the realized copula between the ES and TY (log-)variance processes in Panel A of Figure (ref). In the graph, we also report the tail dependence of the Gumbel copula, again with $\hat{ \theta}_{ \text{MLE}} = 1.4786$. We recall that for the Gumbel copula, $\lambda_{L} = 0$ and $\lambda_{U} = 2 - 2^{1/ \theta}$, which are superimposed as horizontal dotted lines. Overall, there is compelling evidence in favor of asymmetric tail dependence. While the lower tail concentration fizzles out and converges to zero, there is a nontrivial upper tail parameter, which suggests that elevated levels of volatility in the equity and bond markets are likely to occur in parallel. This feature is embedded in the Gumbel copula. Indeed, the tail concentration function of the realized copula of volatility and the associated one of the Gumbel copula are extremely close with a minimal distance, except for a small deviation in the far upper-right tail. We observe that the Gumbel copula has a Kendall’s rank correlation of $\rho_{ \text{K}} = 1 - 1 / \theta$. It thus describes variables that exhibit positive comovement, which is stronger for extremely large values of $\theta$. We estimate an upper tail concentration of about $\hat{ \lambda}_{U} = 0.40$ and a rank correlation of $\hat{ \rho}_{ \text{K}}^{a} = 0.50$.

We corroborate this in Panel B of the figure by showing the estimated Kendall's distribution function. We overlay the population counterpart of the Gumbel copula. The latter is expressible in closed-form, i.e. $K(z) = z - \log(z^{1/ \theta})$. We observe a near-perfect evolution with this measure over the entire support.

To assess goodness-of-fit, we examine a composite null hypothesis with the Gumbel, Clayton, and Independence copulas as contenders for describing the data. The Clayton and Independence copulas are both heavily rejected with a $P\text{-value}$ close to zero, however. Conversely, the Gumbel copula is more appropriate with a $P\text{-value} = 0.1136$.

Taken together, the goodness-of-fit test coupled with the descriptive analysis from Figures (ref) -- (ref) provide overwhelming evidence that the Gumbel copula offers a very good approximation of the codependency between the variance processes of the ES and TY futures contracts.

Conclusion

We propose the realized copula of volatility as a nonparametric tool for studying the codependency in the stochastic volatility component of a multivariate continuous-time asset price process. In analogy with the classical copula framework, the aim is to decouple the marginal behavior of the univariate volatility processes from their dependence structure. The distinctive challenge in this setting is that volatility is latent and has to be recovered with a realized measure constructed from discretely observed high-frequency data of the price process. Thus, inference on the realized copula of volatility requires control of measurement error.

We start by proposing a multivariate extension of the realized distribution function of volatility by li-todorov-tauchen:13a, from which the realized copula of volatility is derived. Under in-fill asymptotics with a fixed time span, we then show that it affords a consistent estimator of the empirical copula of the latent stochastic volatility (measuring the codependency of volatility on the part of the sample path observed so far). In a double-asymptotic framework with the long-span dimension also going to infinity, our estimator converges to the stationary marginal copula of volatility, provided it exists. We derive a functional central limit theorem for the realized-copula process, which supports uniform inference on the volatility copula. To render the limit theory operational, we design a feasible estimator of the asymptotic covariance function. Based on this, we develop pointwise confidence intervals, uniform confidence bands, and a goodness-of-fit test to evaluate a hypothesis about the marginal copula of volatility.

A simulation study sheds light on the finite sample properties of the realized copula of volatility and the associated inference procedures under realistic sampling schemes. We show that it accurately recovers the empirical copula of volatility for a fixed time span and converges at the expected rate toward the marginal copula of volatility as the long-span dimension increases. The goodness-of-fit test is found to exhibit good size control and demonstrate good power.

In our empirical application, we construct the realized copula of volatility from transaction price data on futures contracts tracking the aggregate U.S. equity index and treasury bond market. The evidence presented points toward a Gumbel copula with asymptotic upper-tail dependence as being highly representative of the behavior of the realized variance processes in these data, in line with the notion that volatility across major asset classes tends to spike in tandem during financial market stress.

The asymptotic theory developed here for the latent volatility of semimartingale processes and the associated inference based on local realized proxies of spot variance can be extended in several directions in future research. First, the baseline analysis can be pursued in the presence of market microstructure noise by building on the univariate treatment of the realized distribution function of volatility in christensen-thyrsgaard-veliyev:19a. A second avenue is to allow the realized copula of volatility to vary with the state of the economy or with observable macro-financial indicators, so that dependence in volatility can differ across regimes or along the business cycle. A third extension is to embed the volatility copula in a higher-dimensional framework together with a model for time-varying correlations, thereby delivering a unified copula-based representation of the dynamics of the full covariance matrix. These developments require additional structure but can be built directly on the econometrics and inference procedures developed in this paper.