EconBase
← Back to paper

Adaptive Random Bandwidth for Inference in CAViaR Models

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.

84,574 characters · 10 sections · 35 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.

Adaptive Random Bandwidth for Inference in CAViaR Models

abstractThis paper investigates the size performance of Wald tests for CAViaR models engle2004caviar. We find that the usual estimation strategy on test statistics yields inaccuracies. Indeed, we show that existing density estimation methods cannot adapt to the time-variation in the conditional probability densities of CAViaR models. Consequently, we develop a method called adaptive random bandwidth which can approximate time-varying conditional probability densities robustly for inference testing on CAViaR models based on the asymptotic normality of the model parameter estimator. This proposed method also avoids the problem of choosing an optimal bandwidth in estimating probability densities, and can be extended to multivariate quantile regressions straightforward. JEL Codes: C22 Keywords: covariance matrix estimation in quantile regressions, CAViaR models, bandwidth choice, stability conditions for CAViaR DGPs.

Introduction

Financial risk management is at the heart of banks' and financial institutions' activities to guide them in their investment plans, supervisory decisions, risk capital allocations and for external regulations. The use of quantitative risk measures has become essential in financial risk management. One of the most popular risk measures associated with financial portfolios is the value at risk (VaR hereafter). The VaR at probability $\tau \in (0,1)$ of a portfolio is defined as the minimum potential loss that the portfolio may suffer in the worst $\tau$ portion of all possible outcomes over a given time horizon. VaR is very intuitive duffie1997overview and has for instance been incorporated into the 1996 Amendment to the Capital Accord for measuring the market risk in financial positions of each financial institution. Therefore, VaR is still a widely used risk measure even though many approaches to measuring market and credit risks have been proposed in the literature.

Generally, there are three ways to estimate VaR: (i) historical simulations, (ii) semi-parametric approaches and (iii) fully parametric frameworks. Within the class of semi-parametric approaches, it typically includes extreme value theory analyses and quantile regression techniques. In this paper, we focus on quantile regressions for the VaR estimation as quantile regressions are straightforward in studying one quantile of interest and numerically efficient without imposing parametric distributional assumptions.

Despite that the VaR is just a particular quantile of future portfolio losses conditional on present information, it is essentially a part of the underlying conditional distribution. VaR models are supposed to embrace features of the empirical conditional distributions of returns, such as time-variation and conditional heteroskedasticity. Drawing on (G)ARCH specifications which capture the presence of time-varying conditional heteroskedasticity in time series, engle2004caviar have proposed to estimate conditional autoregressive value at risk by regression quantiles (CAViaR). It is appealing to consider CAViaR models for estimating VaR as CAViaR models associate the conditional quantile of interest with observable variables as well as the implicit information on lagged conditional quantiles.

This paper carefully investigates the size performance of Wald tests for CAViaR models. Having an accurate test statistic is important to obtain reliable models in financial applications. Several specifications are nested within a CAViaR specification, such as static quantile regressive models and quantile autoregressive models koenker2006quantile,hecq2020selecting. Moreover, there exists several models nested within the general CAViaR specification that have been proposed in the literature. For instance, asymmetric slope CAViaR models engle2004caviar that split the effect of positive and negative yesterday's news shocks. Wald tests are used to test the null of a symmetric news impact. However, we find that the usual estimation strategy yields inaccuracies. Indeed, we show that existing density estimation methods cannot adapt to the time-variation in the conditional probability densities of CAViaR models. The method that we develop in this paper is able to adapt to time-varying conditional probability densities and produces much more reliable results than the existing ones for inference testing on CAViaR models based on the asymptotic normality of the model parameter estimator. This proposed method also avoids the haunting problem of choosing an optimal bandwidth in estimating probability densities, and can be extended to multivariate quantile regressions straightforward in theory.

The remainder of this paper is structured as follows. In Section (ref), stability conditions for CAViaR data generating processes (DGPs) to be non-explosive are derived. In Section (ref), we investigate the size performance of Wald tests for CAViaR models and find large size distortions by the usual estimation strategy. So we introduce a method called adaptive random bandwidth. An empirical study on stock returns is performed in Section (ref). Finally Section (ref) concludes this paper.

The CAViaR model

Let us consider a stationary time series process $ \left\{ y_{t}\right\} _{t=1}^{T}$ for instance the return of an asset or a portfolio, and denote $\boldsymbol{x}_{t}$ a vector of observable variables at time $t$ and $\mathcal{F}_t$ the information set up to time $t$ which is the $\sigma$-algebra generated by $\left\{\boldsymbol{x}_{t}, y_t, \boldsymbol{x}_{t-1}, y_{t-1}, \ldots \right\}$. The $\tau$-th quantile ($\tau \in (0,1)$) or the opposite $\text{VaR}_{\tau}$ of $y_{t}$ conditional on $\mathcal{F}_{t-1}$ is denoted as $ f_{t}(\boldsymbol{\beta }_{\tau },\boldsymbol{x}_{t-1})$ (or simply $f_{t}(\boldsymbol{\beta _{\tau }})$ when $\boldsymbol{x}_{t-1}$ is taken in obviously). A generic CAViaR specification proposed by engle2004caviar is

equation[equation omitted — 207 chars of source]

where $\boldsymbol{\beta _{\tau }}^{\prime }:=\left[ \beta _{0},\beta _{1},\ldots ,\beta _{p}\right] $ collects the $p=q+r$ slope parameters, and $l$ is a function of a finite number of lagged observable variables, for instance the lagged returns entering potentially with different weights for positive and negative past lagged returns. As described in engle2004caviar the autoregressive terms $\beta _{i}f_{t-i}(\boldsymbol{\beta }_{\tau })$ can ensure that the quantile changes smoothly over time. The quantile autoregressive model (QAR) of Koenker and Xiao (2006) is nested in the CAViaR specification by restricting $\beta _{1}=...=\beta _{q}=0$ in CAViaR. The role of $l(\boldsymbol{x}_{t-j})$ is to account for the association of $f_{t}(\boldsymbol{\beta } _{\tau })$ with observable variables in $\mathcal{F}_{t-1}$. CAViaR models as a generalization of QAR models are able to capture the time-variation in the conditional quantile in a way similar to GARCH models in explaining time-varying volatility and volatility clustering in financial time series in addition to ARCH models.

The CAViaR model (ref) is nonlinear in parameters as long as there exists a nonzero $\beta _{i},i\in \left\{ 1,\ldots ,q\right\} $ which leads to $\frac{\partial f_{t}(\boldsymbol{\beta }_{\tau })}{\partial \beta _{i}}=f_{t-i}(\boldsymbol{ \beta }_{\tau })+\beta _{i}\frac{\partial f_{t-i}(\boldsymbol{\beta }_{\tau })}{\partial \beta _{i}}$ not independent of $\beta _{i}$.\footnote{In Appendix (ref), the gradient and the Hessian matrix of CAViaR models are illustrated to emphasize that the nonlinearity of model parameters makes CAViaR models different from other linear quantile regression models.} The algorithm to estimate CAViaR models is given in Section (ref).

For illustration, we simulate samples from the following three CAViaR DGPs in (ref) and plot Figure (ref) (a). \footnote{All the simulations of CAViaR DGPs in this paper follow the procedure given in Appendix (ref).} In Figure (ref) (a), we see a decreasing trend in CAViaR DGP 1.a mainly due to the negative term $-0.5|y_{t-1}|$ in $f_{t}(\boldsymbol{\beta }_{\tau})$ compared with CAViaR DGP 1.b. Comparing CAViaR DGP 1.b with 1.c, we find that CAViaR DGP 1.b has a larger spread due to a higher slope of $ f_{t-1}(\tau)$ in $f_{t}(\boldsymbol{\beta }_{\tau})$. A similar finding further applies on Figure (ref) (b) which plots simulated samples of CAViaR DGP 2.a, 2.b and 2.c in (ref) respectively.

equation[equation omitted — 502 chars of source]

where $\{u_{t}\} $ is i.i.d. in the standard uniform distribution (denoted as $\mathcal{U}(0,1)$) and $F_{t(3)}^{-1}(\cdot)$ is the inverse function of Student's t-distribution with $3$ degrees of freedom ($t(3)$ hereafter).

equation[equation omitted — 517 chars of source]

where $u_t \overset{i.i.d.}{\sim} \mathcal{U}(0,1)$, $t=1,2,\ldots, T$.

figure[figure omitted — 218 chars of source]

The stability conditions for CAViaR models

The stationarity of CAViaR time series is required for the model estimation consistency engle2004caviar. After simulating a CAViaR DGP, we can view its behaviour such as explosiveness in the long run. We know that a time series is explosive if and only if at least one conditional quantile of the time series with nonzero probability density to occur is explosive. So we derive stability conditions for the conditional $\tau$-th ($\tau\in(0,1)$) quantile of a CAViaR DGP $\left\{ y_{t}\right\} $ specified as follows:

equation[equation omitted — 224 chars of source]

where $u_t \overset{i.i.d.}{\sim} \mathcal{U}(0,1)$, and $\boldsymbol{\beta _{u_{t}}}^{\prime }:=\left[ \beta _{0}(u_{t}),\beta _{1}(u_{t}),\ldots ,\beta _{p}(u_{t})\right] $ with $p=q+r$. There is a monotonicity requirement on this model which is that $f_{t}( \boldsymbol{\beta }_{u_{t}})$ is monotonically increasing in $u_{t}$ so that the $\tau $-th quantile ($\tau \in (0,1)$) of $ y_{t}$ conditional on $\mathcal{F}_{t-1}$ can be expressed as $f_{t}( \boldsymbol{\beta _{\tau }})$.

Assume the conditional $\tau $-th quantile of $\left\{ y_{t}\right\} $ follows the model (ref) with nonzero probability density to occur at each time. Without loss of generality, there is a time $t\in \left\{ 1,\ldots ,T\right\} $ such that $$ y_{t}=f_{t}(\boldsymbol{\beta }_{\tau })=\beta _{0}+\sum\limits_{i=1}^{q}\beta _{i}f_{t-i}(\boldsymbol{\beta } _{\tau })+\sum\limits_{j=1}^{r}\beta _{q+j}\,y_{t-j}. $$ Now let us derive the value of $y_{t}$. First we have the following equation from (ref).

equation*[equation* omitted — 820 chars of source]

where the second line is obtained by substituting the specification (ref) of $f_{t-i-i_1}(\boldsymbol{\beta}_{\tau})$ into the first line, and $L$ is the lag operator. Further rewrite the above equation, and we have

equation*[equation* omitted — 514 chars of source]

We continue to rewrite the lagged terms of $f_{t}(\boldsymbol{\beta }_{\tau })$ on the right-hand side of the above equation, and then organize the equation such that only the left-hand side contains terms of $y_t$. Therefore, we obtain that

equation[equation omitted — 786 chars of source]

Now we can get the first necessary condition for $\left\{ y_{t}\right\} $ to be nonexplosive, which is

equation[equation omitted — 87 chars of source]

Under the condition (ref), we can simplify the equation (ref) when letting $n\rightarrow \infty $ as follows:

equation*[equation* omitted — 217 chars of source]

Now we obtain the autoregressive polynomial $g(x)$ of $y_{t}$ which is

equation[equation omitted — 147 chars of source]

So the second necessary condition for $\left\{ y_{t}\right\} $ to be nonexplosive is that the roots of $g(x)$ are outside the unit circle. When there exists at least one $\beta_i \neq 0, i\in\{1,\ldots,q\}$, this second condition is equivalent to require that the roots of $g_{1}(x):=1-\sum\limits_{i=1}^{q}\beta _{i}x^{i}-\sum\limits_{j=1}^{r}\beta_{q+j}\,x^{j}$ and the common roots of $g_{2}(x):=1 - \sum\limits_{i=1}^{q}\beta _{i}x^{i}$ and $g_{3}(x):= 1 - \sum\limits_{j=1}^{r}\beta_{q+j}\,x^{j}$ all are outside the unit circle. More examples of CAViaR DGPs are illustrated in Appendix (ref), among which we can find explosive DGPs which break the condition on the roots of $g_1(x)$ but meet the condition on the common roots of $g_{2}(x)$ and $g_{3}(x)$.

We can review Figure (ref) with the above stability conditions. CAViaR DGP 1.b, 1.c, 2.b and 2.c meet the above conditions and we also see their nonexplosive behaviours in the plots. The nonexplosiveness of CAViaR DGP 2.a can also be ensured since it has a narrower spread in theory in comparison with CAViaR DGP 2.b. On the other hand, we know that CAViaR DGP 1.a has a downward trend due to the negative term $-0.5|y_{t-1}|$ and hence is explosive.

Estimation algorithm

The estimation for CAViaR models can be achieved by the differential evolutionary genetic algorithm storn1997differential used by engle2004caviar. Suppose the model $f_{t}(\boldsymbol{\beta })$ is specified as (ref) for data $\left\{ y_{t}\right\} _{t=1}^{T}$. We want to obtain the parameter estimator $\widehat{\boldsymbol{\beta }}$ by the following optimization:

equation[equation omitted — 318 chars of source]

where $S_T(\bm{\beta})$ is the objective function in quantile regressions, and $\rho _{\tau }(x):=x\left( \tau -\bm{1}\{x<0\}\right) $ is called check function koenker2005quantile with the indicator function $\bm{1}\{\cdot \}$.

Following the steps below, we can obtain $\widehat{\boldsymbol{\beta }}$ in (ref).

enumerate• Generate $n$ (say $10^4$) trial vectors independently from a uniform distribution $\mathcal{U}(\bm{b}_L,\bm{b}_p)$ as $n$ parameter initial trials, where $\bm{b}_L$ and $\bm{b}_p$ are $(p+1)\times 1$ vectors roughly covering the lower and upper bounds of the true parameter vector $\boldsymbol{\beta }_{\tau}^o$ of the underlying process in our belief. It is worth mentioning that the values of $\left\{ f_{1-i}(\boldsymbol{\beta}_{\tau}), i = 1,\ldots,q \right\} $ and $\left\{ y_{1-j}, j = 1,\ldots,r \right\} $ acting as initial conditions are also input-demanded in order to calculate $\{f_{t}(\boldsymbol{\beta })\}_{t=1}^T$ for any $\bm{\beta}\in \mathbb{R}^{p+1}$. For instance, as used by engle2004caviar $f_{0}(\boldsymbol{\beta}_{\tau})$ is given as the estimated $\tau$-th quantle of $\left\{ y_t\right\}_{t=1}^{\lfloor 0.1\,T\rfloor}$ and is fixed in the optimization.\footnote{ $\lfloor\cdot\rfloor$ is known as the floor function (or the greatest integer function) and $\lfloor\cdot\rfloor: \mathbb{R} \to \mathbb{Z}$ of a real number $x$ denotes the greatest integer less than or equal to $x$.} • Each parameter initial is used to kick off a minimization routine\footnote{ The Nelder Mead simplex algorithm is used in our minimization routine.} on the objective function $S_T(\bm{\beta})$, and the returned value of $\widehat{\boldsymbol{\beta }}$ from the routine and its objective function value are stored. • Select $m$ (say 10) returned vectors of $\widehat{\boldsymbol{\beta }}$ which result in the lowest $m$ values among the $n$ stored objective function values. • Denote the $m$ selected vectors as $\widehat{\boldsymbol{\beta}}^{(1)}, \ldots, \widehat{\boldsymbol{\beta}}^{(m)}$ and use them as initials to restart the minimization routine individually, and update $\widehat{\boldsymbol{ \beta}}^{(1)}, \ldots, \widehat{\boldsymbol{\beta}}^{(m)}$ with the newly returned vectors respectively. • Repeat Step 4 $a$ (say $5$) times. • Calculate $S_T(\widehat{\boldsymbol{\beta}}^{(i)}), i=1,\ldots,m$. And set the solution to be $\widehat{\boldsymbol{\beta }} = \operatorname*{arg\,min}\limits_{i=1,\ldots,m} S_T(\widehat{\boldsymbol{\beta}}^{(i)})$.

We implement the above estimation algorithm throughout this paper for CAViaR model parameter estimations. There might be a concern if the artificial input of the initial values $\left\{ f_{1-i}(\boldsymbol{\beta}_{\tau}), i = 1,\ldots,q \right\} $ and $\left\{ y_{1-j}, j = 1,\ldots,r \right\} $ affects the parameter estimator. In fact, the effect usually is small and can be neglected when the sample size is large enough because the fitted conditional quantiles $\{f_{t}(\widehat{\bm{\beta }})\}_{t=1}^T$ are kept close to the true ones $\{f_{t}(\bm{\beta })\}_{t=1}^T$ such that it can minimize the objective function despite some burn-in period.

Adaptive random bandwidth method for CAViaR covariance matrix estimation

Consistency and asymptotic normality of CAViaR model parameters have been proved by engle2004caviar. After regressing data onto a CAViaR model, we would like to implement an inference testing on whether the model is correctly specified. In this section we first investigate how we result in the asymptotic normality of CAViaR model parameter estimators. We focus on the elements of the asymptotic covariance matrix to highlight their roles in connecting sample elements with the corresponding limit behaviours. Next, we check whether existing estimation strategies can perform robustly and satisfactorily for Wald tests on CAViaR models. Finally, we propose a new method called adaptive random bandwidth for CAViaR models.

Asymptotics of CAViaR

Consider a time series $\{y_t\}$ of random variables $y_t$ on a complete prbability space $(\Omega, \mathcal{F}, P)$ \footnote{See the assumption C0 of engle2004caviar. We also apply this assumption throughout this paper. That is to say, all the random variables considered in this paper are assumed on a complete prbability space $(\Omega, \mathcal{F}, P)$. }. For applying a generic CAViaR model (ref) on $\{y_t\}$, the consistency and asymptotic normality of the estimator $\widehat{\boldsymbol{\beta}}:= \operatorname*{arg\,min}\limits_{\bm{\beta}}\sum\limits_{t=1}^{T} \rho_{\tau}\left( y_t - f_t(\boldsymbol{\beta}) \right)$ has been derived out by engle2004caviar:

theorem[Asymptotics given by engle2004caviar] \ \\ For a data generating process $\{ y_t \}$ with its time $t$ conditional $\tau$-th quantile following a generic CAViaR model as (ref) parametrized by $\bm{\beta}^{o}$, it satisfies the regularity conditions (C0,\ldots, C7, AN1,\ldots, AN7) in the proof of engle2004caviar. Then \begin{equation} \sqrt{T}A_T^{-1/2}D_T \left(\widehat{\bm{\beta}} - \bm{\beta}^{o} \right) \overset{\mathcal{D}}{\sim} N(\bm{0},\bm{I}_{(p+1)\times(p+1)}), \end{equation} where \begin{equation} \begin{aligned} \widehat{\boldsymbol{\beta}} & := \operatorname*{arg\,min}\limits_{\bm{\beta}\in\mathbb{R}^{p+1}}\sum\limits_{t=1}^{T} \rho_{\tau}\left( y_t - f_t(\boldsymbol{\beta}) \right), \\ A_T & := \mathbb{E}\left[ T^{-1} \tau(1-\tau)\sum^{T}_{t=1} \nabla'f_t(\bm{\beta^{o}} ) \nabla f_t(\bm{\beta^{o}} )\right], \\ D_T & := \mathbb{E}\left[ T^{-1} \sum^{T}_{t=1} h_t\left(0 \middle| \mathcal{F}_{t-1} \right) \nabla'f_t(\bm{\beta^{o}} ) \nabla f_t(\bm{\beta^{o}} )\right], \\ \epsilon_{\tau\,t} & := y_t - f_t(\bm{\beta}^{o}), \end{aligned} \end{equation} and $h_t(0| \mathcal{F}_{t-1})$ is denoted as the probability density of $\epsilon_{\tau\,t}$ evaluated at $0$ conditional on the information set $\mathcal{F}_{t-1}$. $\bm{I}_{(p+1)\times(p+1)}$ is the $(p+1)\times(p+1)$ identity matrix.

The above theorem is useful for quantile model (mis)specification tests. For instance, Wald tests can be used to check whether the current model is correctly specified by testing the validity of a more parsimonious nested model. To perform such a quantile model specification test, it often requires to estimate $A_T$, $h_t(0| \mathcal{F}_{t-1})$ and $D_T$. When using traditional estimates $\widehat{A}_T$, $\{\widehat{h}_t(0| \mathcal{F}_{t-1})\}$, $\widehat{D}_T$ of $A_T$, $\{h_t(0| \mathcal{F}_{t-1})\}$ and $D_T$ respectively, we found considerable size distortions in inference tests on CAViaR models in general. We will show that the reason lies in the inaccuracy of $\{\widehat{h}_t(0| \mathcal{F}_{t-1})\}$ in the next subsection. In order to spot the discrepancy in approximating $\{h_t(0| \mathcal{F}_{t-1})\}$, we need a clear picture on how $\left\{ h_{t}(0|\mathcal{F}_{t-1})\right\} $ comes up into the asymptotic normality of the model parameter estimator. Doing so, we can see the role of $\left\{ h_{t}(0|\mathcal{F}_{t-1})\right\} $ and whether a sequence $\left\{ \widehat{h_{t}}(0|\mathcal{F} _{t-1})\right\} $ is capable to achieve the same role in practice. Let us review the proof of engle2004caviar for Theorem (ref) below.

The proof of engle2004caviar is obtained by applying Theorem 3 of huber1967behavior onto $T^{-1/2}\displaystyle\sum^{T}_{t=1} \left( \bm{1}\left\{ y_t \leq f_t(\widehat{\boldsymbol{\beta}}) \right\} - \tau\right)\nabla' f_t(\widehat{\boldsymbol{\beta}})$ and the central limit theorem onto $T^{-1/2}\displaystyle\sum^{T}_{t=1} \left( \bm{1}\left\{ y_t \leq f_t(\bm{\beta}^{o}) \right\} - \tau\right)\nabla' f_t(\bm{\beta}^{o})$. Huber's conditions are verified in the proof before applying Huber's theorem. Denote

equation[equation omitted — 213 chars of source]

$\text{Hit}_t(\bm{\beta}) $ gives value $-\tau$ every time $y_t$ exceeds $f_t(\boldsymbol{\beta})$ and $1-\tau$ otherwise. With the true underlying parameter $\bm{\beta}^{o}$, $\{ \text{Hit}_t(\bm{\beta}^{o}) \}$ is a martingale difference sequence with respect to $\{\mathcal{F}_{t-1}\}$. It is easy to get that $T^{-1/2}\displaystyle\sum^{T}_{t=1} \text{Hit}_t(\bm{\beta}^{o}) g_t(\bm{\beta}^{o})$ follows the central limit theorem because $\left\{ \text{Hit}_t(\bm{\beta}^{o}) g_t(\bm{\beta}^{o}) \right\}$ is a martingale difference sequence with the assumption AN1 of engle2004caviar on its uniformly bounded second moment. So we get that

equation[equation omitted — 159 chars of source]

It has also been proved by engle2004caviar that

equation[equation omitted — 153 chars of source]

Next, we are going to manifest $\left\{h_t(0|\mathcal{F}_{t-1})\right\}$ in the proof in a way which makes the appearance of $\left\{h_t(0|\mathcal{F}_{t-1})\right\}$ more intuitive. We rewrite $\text{Hit}_t(\widehat{\bm{\beta}})\, g_t(\widehat{\bm{\beta}})$ as follows:

equation[equation omitted — 478 chars of source]

Take expectation on the both sides of Equation (ref) and get

equation[equation omitted — 1,868 chars of source]

where $\norm{\cdot}_{\infty}$ is the supremum norm of vectors. And

equation[equation omitted — 1,261 chars of source]

where $F_t\left( \cdot\middle| \mathcal{F}_{t-1} \right) $ is the probability density function of $y_t$ conditional on $\mathcal{F}_{t-1}$, and $ h_t(0|\mathcal{F}_{t-1})= F_t'\left( f_t(\bm{\beta}^{o} ) \middle| \mathcal{F}_{t-1} \right) $. Substituting (ref) into (ref) gives

equation[equation omitted — 543 chars of source]

\\ Success in applying Huber's theorem gives

equation[equation omitted — 350 chars of source]

Therefore, the asymptotic normality of $ T^{1/2}\left(\widehat{\bm{\beta}} - \bm{\beta}^{o} \right)$ is obtained by substituting (ref) and (ref) into (ref).

From the above derivation, it is clear that the role of $ h_t(0|\mathcal{F}_{t-1})$ is actually an approximation to $F_t'\left( f_t(\bar{\bm{\beta}} ) \middle| \mathcal{F}_{t-1} \right)$ in which $\bar{\bm{\beta}}$ is between $\bm{\beta}^{o}$ and $\widehat{\bm{\beta}}$. This role comes to the surface of (ref) using the fact that

equation[equation omitted — 340 chars of source]

by the Mean Value Theorem. This approximating role of $h_t(0|\mathcal{F}_{t-1})$ sets a clear mission of any $\widehat{h_t}(0|\mathcal{F}_{t-1})$ supposed to achieve, which can be used to examine an estimator for $h_t(0|\mathcal{F}_{t-1})$ as well as to propose an improved estimation method. In next subsection, we are going to examine the performances of some existing methods for estimating $h_t(0|\mathcal{F}_{t-1})$ and the role of $h_t(0|\mathcal{F}_{t-1})$ will help to find out the intrinsic defects of those methods.

Existing methods for CAViaR covariance matrix estimation

Based on the literature on quantile regressions, in general there are two ways to estimate $\left\{ h_t\left(0 \middle| \mathcal{F}_{t-1} \right) \right\}$ in $D_T$ with $\{ \epsilon_{\tau\,t}\}$ being potentially non-i.i.d.. One is referred to as the Hendricks Koenker Sandwich Approach hendricks1992hierarchical, koenker2005quantile analogous to the finite difference idea resulting in the estimator $\widehat{h_t}^{fd}(0 | \mathcal{F}_{t-1} ) $ for $h_t(0 | \mathcal{F}_{t-1} )$ as follows:

equation[equation omitted — 242 chars of source]

where $\Delta\tau_T$ is subject to $ 0<\tau \pm \Delta\tau_T <1$ with $\Delta\tau_T \rightarrow 0$ as $T \rightarrow \infty$. The other one is referred to as the Powell Sandwich powell1991estimation, koenker2005quantile based on the kernel density estimation idea resulting in the estimator $\widehat{h_t}^{kernel}(0 | \mathcal{F}_{t-1} ) $ for $h_t(0 | \mathcal{F}_{t-1} )$ as follows:

equation[equation omitted — 442 chars of source]

where $K(\cdot)$ is a suitable kernel function with bandwidth $2\,c_T$ and $c_T \rightarrow 0$ as $T \rightarrow \infty$. As we can see in (ref), one kernel function is applied throughout $\{y_t\}$ with $ y_t - f_{t}(\boldsymbol{\beta}_{\tau})$ being the only distinguishable information for $\widehat{h_t}^{kernel}(0 | \mathcal{F}_{t-1} ) $. Therefore, this kernel method does not capture sufficient information to distinguish time-varying conditional distributions of $\{y_t\}$, and consequently cannot fully adapt to the time-variations. Additionally, the choice of the kernel function $K(\cdot)$ and the bandwidth parameter $c_T$ are still in a lot of nettlesome questions in practice. A similar issue in the Hendricks Koenker Sandwich Approach is on choosing $\Delta\tau_T $ and extra error resulted from estimating $ f_{t}(\boldsymbol{\beta}_{\tau + \Delta\tau_T}) $ and $ f_{t}(\boldsymbol{\beta}_{\tau - \Delta\tau_T}) $.

The estimation method adopted by engle2004caviar is a form of the Powell Sandwich as follows:

equation[equation omitted — 196 chars of source]

As suggested by koenker2005quantile and machado2013quantile, the bandwidth $\widehat{c}_T$ generally adopted is defined as follows:

equation[equation omitted — 121 chars of source]

where $m_T$ is defined as

equation[equation omitted — 203 chars of source]

with $\Phi(\cdot)$ and $\phi(\cdot)$ being the cumulative distribution and probability density functions of $N(0,1)$ respectively. And $\widehat{k}_T $ is defined as the median absolute deviation of the conditional $\tau$-th quantile regression residuals.

Wald tests are applied in this subsection to check the performances of the above estimation methods for CAViaR models.

First, we consider the following candidate model specifications for the conditional $\tau$-th ($\tau\in (0,1)$) quantile of a time series $\left\{ y_t \right\}$ with $f_t(\boldsymbol{\beta}_{\tau})$ denoted as the $\tau$-th quantile of $y_t$ conditional on the information set $\mathcal{F}_{t-1}$. \[ \left\{

tabular[tabular omitted — 1,357 chars of source]

\right. \]

The models (ref) and (ref) are nested within model (ref). Now let us consider the Wald test on models (ref) and (ref) first. Simulate a time series $\left\{ y_t \right\}$ with its DGP specified as the model (ref) with the underlying parameter vector $\bm{\beta}_{u_t}^{R1} = [F^{-1}_{N(0,1)}(u_t), 0.2, 0.3]'$, where $\left\{ u_t \right\}\overset{i.i.d.}{\sim} \mathcal{U}(0,1)$ and $F^{-1}_{N(0,1)}(\cdot)$ is the inverse standard normal probability distribution function. The sample size of each simulated sample is $4000$. Conditional $50\%$-th quantiles are estimated for each of total $1000$ simulated samples in this DGP by regressing the sample onto the full model (ref). The Wald test implemented here consists of the null hypothesis of the form $\mathit{H_0}: R\bm{\beta}_{\tau}^{FM}= \gamma$, where $R = [ 0 , 0, 1, -1] $, $\gamma = 0$, and $\widehat{\bm{\beta}}_{\tau}$ is the estimator of the full model parameter vector in (ref). The Wald test statistic denoted by $W_T$ is formulated weiss1991estimating as follows:

equation[equation omitted — 229 chars of source]

where $\widehat{A}_T$ and $\widehat{D}_T$ are estimates for $A_T$ and $D_T$ in (ref) respectively. It is straightforward to obtain $\widehat{A}_T$ and $\widehat{D}_T$ by plugging in $\widehat{\bm{\beta}}_{\tau}$ and $\left\{ \widehat{h_t}\left(0 \middle| \mathcal{F}_{t-1} \right) \right\}$, i.e., $$ \left\{

aligned\widehat{A}_T & = T^{-1} \tau(1-\tau)\sum^{T}_{t=1} \nabla'f_t( \widehat{\bm{\beta}}_{\tau} ) \nabla f_t( \widehat{\bm{\beta}}_{\tau} ), \\ \widehat{D}_T & = T^{-1} \sum^{T}_{t=1} \widehat{h_t}\left(0 \middle| \mathcal{F}_{t-1} \right) \nabla'f_t(\widehat{\bm{\beta}}_{\tau} ) \nabla f_t(\widehat{\bm{\beta}}_{\tau} ).

\right. $$

Notations on $\widehat{D}_T$ to distinguish different estimators used for $\left\{ h_t\left(0 \middle| \mathcal{F}_{t-1} \right) \right\}$ are given by

equation[equation omitted — 295 chars of source]
equation[equation omitted — 337 chars of source]

where $\widehat{c}_T$ is determined as (ref).

We are going to examine each element in the estimation of $D_T$. The analytic solution to $ h_t\left(0 \middle| \mathcal{F}_{t-1} \right)$ can be obtained as follows:

equation[equation omitted — 552 chars of source]

where $\beta_0'(\tau):=\frac{\partial \beta_0(\tau)}{\partial \tau}$. The last line is obtained by knowing $\abs{\beta_1} < 1$. The analytic solution to $ h_t\left(0 \middle| \mathcal{F}_{t-1} \right)$ is used to help identify inaccurate elements in $\widehat{D}_T$ by comparing the test performances of using $\widehat{D}_T^{ker}$, $\widehat{D}_T^{fd}$ and the following

equation[equation omitted — 254 chars of source]

The test performances of using $\widehat{D}_T^{ker}$, $\widehat{D}_T^{fd}$ and $\widehat{D}_T^{h_0}$ are shown in Table (ref) and (ref), which are compared together with the Wald test result using the true underlying parameter vector $\bm{\beta}^{FM}_{\tau} = [F^{-1}_{N(0,1)}(\tau), 0.2, 0.3]', \tau\in(0,1)$ into

equation[equation omitted — 449 chars of source]

where $\phi(\cdot)$ is the probability density function of $N(0,1)$.

The size performances of the Wald tests on the models (ref) and (ref) using different $D_T$ estimators are listed in Table (ref) in which each estimated size is obtained by the percentage rejection rate among the 1000 samples of $T=4000$ in the DGP (ref). Analogously, we implement the Wald test on models (ref) and (ref) with the underlying DGP $\left\{ y_t \right\}$ specified as the model (ref) with the underlying parameter vector $\bm{\beta}_{u_t}^{R2} = [F^{-1}_{N(0,1)}(u_t), 0.2, 0.3]'$, where $\left\{ u_t \right\}\overset{i.i.d.}{\sim} \mathcal{U}(0,1)$. The number of observations in each stimulated sample from this DGP is $4000$. Conditional $50\%$-th quantiles are estimated for each of $1000$ simulated samples by regressing the sample onto the full model (ref). The Wald test implemented in this case consists of the null hypothesis of the form $\mathit{H_0}: R\bm{\beta}_{\tau}^{FM}= \gamma$, where $R = [ 0 , 0, 1, 1] $, $\gamma = 0$, and $\widehat{\bm{\beta}}_{\tau}$ is the estimator of the full model regression (ref). In result, the size performances of the Wald tests on (ref) and (ref) are listed in Table (ref).

From Table (ref) and (ref), we can see large size distortions with $\widehat{D}_T^{fd}$, unlike $\widehat{D}_T^{ker}$, $\widehat{D}_T^{h_0}$ or $\widehat{D}_T^{0}$ that are performing in line with the nominal size. This comparison points out the crucial element estimation to the accuracy of $\widehat{D}_T$ which is $\left\{ \widehat{h_t}\left(0 \middle| \mathcal{F}_{t-1} \right) \right\}$. To check whether $\{ \widehat{h_t}^{ker}\left(0 \middle| \mathcal{F}_{t-1} \right)\}$ is capable to achieve the role of $\{ h_t\left(0 \middle| \mathcal{F}_{t-1} \right) \}$ robustly for time-varying conditional probability densities, we consider the following DGP:

align[align omitted — 529 chars of source]

where $\left\{ u_t \right\}\overset{i.i.d.}{\sim} \mathcal{U}(0,1)$ and the underlying parameters are given as $\bm{\beta}_{u_t}^{R3} = [F^{-1}_{N(0,1)}(u_t), 0.2, 0.3]'$. The analytic form of the corresponding conditional probability density $ h_t\left(0 \middle| \mathcal{F}_{t-1} \right)$ of $y_t$ at its $\tau$-th quantile $ f_t(\bm{\beta}_{\tau})$ given $\mathcal{F}_{t-1}$ can be derived out as follows:

equation[equation omitted — 552 chars of source]

where the first equation is obtained by iteratively rewriting $ \frac{\partial f_{t-i}(\boldsymbol{\beta}_{\tau})}{\partial \tau} $ at each $i$ and knowing $|\beta_1|<1$. This analytic form of $h_t\left(0 \middle| \mathcal{F}_{t-1} \right) $ in (ref) shows that $\{h_t\left(0 \middle| \mathcal{F}_{t-1} \right) \}$ indeed is time-varying and nonzero with probability one.

We simulate $1000$ samples from the DGP (ref) with $T=5000$, and estimate the conditional $50\%$-th quantiles of each sample by regressing the sample onto the full model specification (ref). The Wald test described as (ref) with $R = [ 0 , 0, 1, -1] $ is performed on these $1000$ samples and the size performance is presented in Table (ref). We see a large size distortion with the kernel method $\widehat{D}_T^{ker}$ in Table (ref). More tests are conducted for different DGPs and together with the results are presented in Appendix (ref). Based on our test results, we see that the kernel method for estimating $\left\{ h_t\left(0 \middle| \mathcal{F}_{t-1} \right) \right\}$ is not robust and cannot fully adapt to time-varying conditional probability densities.

Estimating $\left\{ h_t\left(0 \middle| \mathcal{F}_{t-1} \right) \right\}$ robustly has to be achieved in order to ensure the reliability of CAViaR analysis based on the asymptotic properties of CAViaR model parameter estimators. In seeking for improving the accuracy of $\left\{ \widehat{h_t}\left(0 \middle| \mathcal{F}_{t-1} \right) \right\}$, we bear in mind two guidances. One is the role of $\left\{ h_t\left(0 \middle| \mathcal{F}_{t-1} \right) \right\}$ on how it links sample elements with the corresponding limit behaviours, see Section (ref). The other guidance is the fundamental flaws of $\{ \widehat{h_t}^{ker}\left(0 \middle| \mathcal{F}_{t-1} \right)\}$ and $\{ \widehat{h_t}^{fd}\left(0 \middle| \mathcal{F}_{t-1} \right)\}$ in their accuracy. In terms of $\{ \widehat{h_t}^{fd}\left(0 \middle| \mathcal{F}_{t-1} \right)\}$, $\Delta\tau_T$ needs to be determined properly and two more quantile regressions need to be preformed in order to obtain $\widehat{\bm{\beta}}_{\tau + \Delta \tau_T} $ and $\widehat{\bm{\beta}}_{\tau - \Delta \tau_T}$. The effect of this extra estimation error is crucial to the performance of $\{ \widehat{h_t}^{fd}\left(0 \middle| \mathcal{F}_{t-1} \right)\}$. Although $\{ \widehat{h_t}^{ker}\left(0 \middle| \mathcal{F}_{t-1} \right)\}$ does not need extra quantile regressions, it still requires a proper choice on the kernel function $K(\cdot)$ and the bandwidth $\widehat{c}_T$. Remarkably, $\{ \widehat{h_t}^{ker}\left(0 \middle| \mathcal{F}_{t-1} \right)\}$ does not differentiate the observations within the bandwidth regardless of the number of the observations in the bandwidth while using the kernel function $\bm{1}\{ | y_t - f_t(\widehat{\bm{\beta}}_{\tau} ) | < \widehat{c}_T \} $. Therefore, it is desirable to get rid of choosing bandwidth $\Delta\tau_T$ or $c_T$ and the kernel function $K(\cdot)$ in the estimation. In the next subsection, a robust estimation method for $\left\{ h_t\left(0 \middle| \mathcal{F}_{t-1} \right) \right\}$ is developed up without the need in choosing a bandwidth or a kernel function.

Adaptive random bandwidth method

We have noticed that the accuracy of the $\left\{ h_t(0|\mathcal{F}_{t-1}) \right\}$ estimation is crucial to the performance of inference tests based on the asymptotic normality of CAViaR model parameter estimators. It is also well known that $\{ \widehat{h_t}^{fd}\left(0 \middle| \mathcal{F}_{t-1} \right)\}$ suffers both from the error in estimating $f_t(\bm{\beta}_{\tau + \Delta\tau})$ and $f_t(\bm{\beta}_{\tau - \Delta\tau})$ and from choosing a proper $\Delta\tau_T$. On the other hand, $\{ \widehat{h_t}^{ker}\left(0 \middle| \mathcal{F}_{t-1} \right)\}$ has some fundamental problems. First of all, $\{ \widehat{h_t}^{ker}\left(0 \middle| \mathcal{F}_{t-1} \right)\}$ cannot fully adapt to time-varying conditional distributions of time series due to the fact that the same kernel function $K(\cdot)$ and only timely information $(y_t - f_t)$ are used in estimating $h_t(0|\mathcal{F}_{t-1})$ for all t. Second, finding a proper kernel function $K(\cdot)$ with a proper bandwidth $c_T$ still faces a lot nettlesome problems in practice. Neither of these two methods is practically robust. The goal in this subsection is to develop an estimation method for $\left\{ h_t(0|\mathcal{F}_{t-1}) \right\}$ which can adapt to time-variation characteristics of CAViaR DGPs and is robust in practice without the need to determine a proper bandwidth. We name this estimation method as the adaptive random bandwidth (ARB) method which can reliably bridge asymptotic properties of CAViaR models in theory with CAViaR applications.

The idea of this method is inspired by viewing the role of $\left\{ h_t \left(0 \middle| \mathcal{F}_{t-1} \right) \right\}$ on how it links sample elements with the corresponding limit behaviours, see Section (ref). Reviewing equation (ref), we can explicitly formulate $\left\{ h_t \left(0 \middle| \mathcal{F}_{t-1} \right) \right\}$ as follows:

equation[equation omitted — 339 chars of source]

which actually is a conditional expectation taken with respect to random variables $ y_t $ and $\widehat{\bm{\beta}} $. We use the subscript in $\mathbb{E}$ to clarify the expectation is taken with respect to specific random variable(s) hereafter. Considering this role of $\left\{ h_t \left(0 \middle| \mathcal{F}_{t-1} \right) \right\}$ as well as equation (ref), we are enlightened to use random bandwidth $\nabla'f_t(\widehat{\bm{\beta}}) \left( \bm{b}_i - \widehat{\bm{\beta}} \right)$ with $ \sqrt{T}\left( \bm{b}_i - \widehat{\bm{\beta}} \right) \overset{\mathcal{D}}{\sim} N(\bm{0},\bm{V_d})$ and $i = 1,2,\ldots, n$. We can set $\bm{V_d} = I_{(p+1)\times (p+1)}$ to start with. After sufficient $n$ times Monte Carlo simulating $\bm{b}_i - \widehat{\bm{\beta}} $ from $ N(\bm{0},\bm{V_d}) $, an estimator of $h_t \left(0 \middle| \mathcal{F}_{t-1} \right)$ can be achieved as follows:

equation[equation omitted — 420 chars of source]

After achieving the above $\widehat{h_t}(0|\mathcal{F}_{t-1})$, we can estimate $\widehat{D}_T$ so as to update $\bm{V_d} = \widehat{D}_T^{-1}\widehat{A}_T\widehat{D}_T^{-1}$. Redo the simulation of $\{\bm{b}_i - \widehat{\bm{\beta}} \}_{i=1}^{n}$ with the updated $\bm{V_d}$. We can estimate $\widehat{h_t}(0|\mathcal{F}_{t-1})$ and $\widehat{D}_T$ again. This estimation repetition can mitigate the influence of an arbitrary chosen $\widehat{h_t}(0|\mathcal{F}_{t-1})$ in ARB.

Compared to the Powell Sandwich estimation (ref) with $c_T$, our proposed method uses random bandwidth $\nabla'f_t(\widehat{\bm{\beta}}) \left( \bm{b}_i - \widehat{\bm{\beta}} \right)$ and Monte Carlo simulations such that it can adapt to time-varying conditional distributions of CAViaR DGPs by approaching to the role of $\left\{ h_t \left(0 \middle| \mathcal{F}_{t-1} \right) \right\}$ as in (ref) and in (ref). The adaptive random bandwidth method can remarkably outperform the Powell Sandwich method in the applications on DGPs of time-varying conditional distributions, as shown in Table (ref). In theory, the adaptive random bandwidth method is valid as long as $\bm{b}_i - \widehat{\bm{\beta}}$ and $\bm{\beta}^{o} - \widehat{\bm{\beta}}$ have the same order of magnitude. We formally establish this adaptive random bandwidth method in Theorem (ref).

theorem[Adaptive Random Bandwidth Method] {\ \\} Assume the conditions and the asymptotic normality result in Theorem (ref). Choose an arbitrary positive definite symmetric matrix $\bm{V_d}$. Under the condition that \begin{equation} \sqrt{T}\left( \bm{b}_i - \widehat{\bm{\beta}} \right) \overset{i.i.d.}{\sim} N(\bm{0},\bm{V_d}), \qquad i = 1,\ldots,n, \end{equation} and \begin{equation*} \biggl|\nabla'f_t(\widehat{\bm{\beta}}) \left( \bm{b}_i - \widehat{\bm{\beta}} \right) \biggr| \neq 0, \end{equation*} the adaptive random bandwidth estimator for $h_t(0|\mathcal{F}_{t-1})$ is formulated as follows: \begin{equation} \resizebox{1.08\textwidth}{!}{ $ \widehat{h_t}(0|\mathcal{F}_{t-1}) = \left\{ \begin{aligned} & n^{-1}\displaystyle\sum^{n}_{i=1} \frac{ \bm{1}\left\{ y_t \leq f_t(\widehat{\bm{\beta}}) + \nabla'f_t(\widehat{\bm{\beta}}) \left( \bm{b}_i - \widehat{\bm{\beta}} \right) \right\} - \bm{1}\left\{ y_t \leq f_t(\widehat{\bm{\beta}}) \right\} }{\nabla'f_t(\widehat{\bm{\beta}}) \left( \bm{b}_i - \widehat{\bm{\beta}} \right) }, \quad & \text{when } y_t \neq f_t(\widehat{\bm{\beta}}), \\ & 0, \quad & \text{when } y_t = f_t(\widehat{\bm{\beta}}), \end{aligned} \right. $} \end{equation} such that $$ \mathbb{E}_{y_t,\widehat{\bm{\beta}}}\biggl[ \widehat{h_t}(0|\mathcal{F}_{t-1}) \biggm|\mathcal{F}_{t-1} \biggr] \overset{p}{\longrightarrow} h_t(0|\mathcal{F}_{t-1}) $$ as $n \rightarrow \infty$. \footnote{ We regard the least $(p+1)$ absolute residuals in $\{ | y_t - f_t(\widehat{\bm{\beta}})| \}_{t=1}^T$ as zeros. In fact, iterations of a simplex-based direct search method like the Nelder–Mead method for optimizing $(p+1)$ parameters terminates at the vertices of a simplex in the parameter space lagarias1998convergence. That is to say, the iterations in optimizing the $\tau$-th quantile regression objective function terminate with $(p+1)$ elements of $\{ ( \tau -\bm{1}\{ y_t - f_t(\bm{\beta} ) < 0 \} ) (y_t - f_t(\bm{\beta} ) ) \}$ solved to be zeros. Therefore, we set $\widehat{h_t}(0|\mathcal{F}_{t-1}) = 0$ at the least $(p+1)$ absolute residuals in $\{ | y_t - f_t(\widehat{\bm{\beta}})| \}_{t=1}^T$ in all the tests throughout this paper.}
proofSee Appendix (ref).

\\

We separate the case of $y_t = f_t(\widehat{\bm{\beta}})$ from others to maintain the convergence of the ARB estimator due to $\lim_{x\rightarrow 0} \frac{1}{x} = \infty $. Zero given to $\widehat{h_t}(0|\mathcal{F}_{t-1})$ at $y_t = f_t(\widehat{\bm{\beta}})$ also enables the ARB estimator to approximate $h_t(0|\mathcal{F}_{t-1})$ from the left and from the right in half weights respectively in expectation, see the proof of Theorem (ref). The convergence property of the partial sum in the sequence $\{\widehat{h_t}(0|\mathcal{F}_{t-1})\}$ by ARB is given in Corollary (ref).

corollaryUnder the conditions of Theorem (ref), the adaptive random bandwidth estimator $\{\widehat{h_t}(0|\mathcal{F}_{t-1})\}$ has the following property: \begin{equation} \frac{1}{T}\sum^{T}_{t=1} \widehat{h_t}(0|\mathcal{F}_{t-1}) {\overset{m.s.}{\longrightarrow}}\, \frac{1}{T}\sum^{T}_{t=1} h_t(0|\mathcal{F}_{t-1}), \end{equation} as $T, n \rightarrow \infty$.
proofSee Appendix (ref).

\\

It is clear that both $\widehat{\epsilon}_t := y_t - f_t(\widehat{\bm{\beta}})$ and $\nabla'f_t(\widehat{\bm{\beta}})$ are taken into account by ARB to approximate $h_t(0|\mathcal{F}_{t-1})$. In order to identify how $\widehat{\epsilon}_t$ and $\nabla'f_t(\widehat{\bm{\beta}})$ jointly shape $\widehat{h}_t(0|\mathcal{F}_{t-1})$, we would like to formulate $\widehat{h}_t(0|\mathcal{F}_{t-1})$ in Theorem (ref) into an analytic expression in terms of $\widehat{\epsilon}_t$ and $\nabla'f_t(\widehat{\bm{\beta}})$ so as to manifest the relationship. The analytic form of $\widehat{h_t}(0|\mathcal{F}_{t-1})$ by ARB described in Theorem (ref) is presented in Corollary (ref).

corollaryUnder the conditions of Theorem (ref), we can get the analytic form of $\widehat{h_t}(0|\mathcal{F}_{t-1})$ as follows: \begin{equation} \widehat{h_t}(0|\mathcal{F}_{t-1}) = \left\{ \begin{aligned} & \frac{1}{2\delta_{\nabla_t}\,\sqrt{2\pi}} E_1\left( \frac{ \widehat{\epsilon}_t^2 }{2\delta_{\nabla_t}^2} \right), & \quad when \widehat{\epsilon}_t \neq 0, \\ & 0, & \quad when \widehat{\epsilon}_t = 0, \end{aligned} \right. \end{equation} where $\widehat{\epsilon}_t := y_t - f_t(\widehat{\bm{\beta}})$, $\delta_{\nabla}:=\sqrt{\frac{\nabla'f_t(\widehat{\bm{\beta}}) \bm{V_d}\nabla f_t(\widehat{\bm{\beta}})}{T}} = T^{-\frac{1}{2}}\norm{\nabla f_t}_2$, and $ E_1(s):= \int_{s}^{\infty} x^{-1}e^{-x} d\,x$ is a special integral known as the exponential integral or the incomplete gamma function $\Gamma(0,s)$.
proofSee Appendix (ref).

\\

For visually checking the roles of $\widehat{\epsilon}_t $ and $ \delta_{\nabla_t} $ in the analytic $\widehat{h_t}(0|\mathcal{F}_{t-1})$ in Corollary (ref), we present a level plot of the analytic $\widehat{h_t}(0|\mathcal{F}_{t-1}) $ over $\widehat{\epsilon}_t $ and $ \delta_{\nabla_t} $ in Figure (ref) which uses colors to differentiate different ranges of $\widehat{h_t}(0|\mathcal{F}_{t-1}) $. It is straightforward to get that the analytic $\widehat{h_t}(0|\mathcal{F}_{t-1})$ is decreasing in $|\widehat{\epsilon}_t |$ as also shown in Figure (ref). However, $\delta_{\nabla_t} $, or say $T^{-\frac{1}{2}}\norm{\nabla f_t}_2$, can shift $\widehat{h_t}(0|\mathcal{F}_{t-1})$ by reflecting on how rare an $\widehat{\epsilon}_t $ is observed given the information set $\mathcal{F}_{t-1}$ and the model specification. That is how the information of $\delta_{\nabla_t} $ in ARB shapes $\widehat{h_t}(0|\mathcal{F}_{t-1})$ adaptively to time-varying conditional probability densities.

figure[figure omitted — 293 chars of source]

The ARB estimator $\widehat{h_t}(0|\mathcal{F}_{t-1})$ via simulations in Theorem (ref) performs as robustly as the analytic ARB estimator $\widehat{h_t}(0|\mathcal{F}_{t-1})$ in Corollary (ref), as shown in Table (ref), (ref) and (ref). The analytic way is faster than the simulation one. However, the ARB estimator via simulations is more intuitive and more flexible to adapt to a very different distribution for simulating $\{\bm{b}_i - \widehat{\bm{\beta}} \}_{i=1}^{n}$.

$D_T$ need to be estimated consistently for inference tests on CAViaR models based on the asymptotic normality of the model parameter estimator. $\{\widehat{h_t}(0|\mathcal{F}_{t-1})\}$ by ARB facilitates our estimation on $D_T$ by just plugging in $ \widehat{\bm{\beta}}$ and $\{\widehat{h_t}(0|\mathcal{F}_{t-1})\}$. The resulted estimator $\widehat{D}_T^{arb}$ has the consistency property presented in Theorem (ref).

theoremUnder the conditions of Theorem (ref), we can get that \begin{equation} \widehat{D}_T^{arb} {\overset{p}{\longrightarrow}}\, D_T, \end{equation} as $T \rightarrow \infty$ and $n \rightarrow \infty$, where $\widehat{D}_T^{arb} := T^{-1} \displaystyle\sum^{T}_{t=1} \widehat{h}_t\left(0 \middle| \mathcal{F}_{t-1} \right) \nabla'f_t( \widehat{\bm{\beta}}_{\tau} ) \nabla f_t(\widehat{\bm{\beta}}_{\tau} ) $ and $\widehat{h}_t\left(0 \middle| \mathcal{F}_{t-1} \right)$ is the adaptive random bandwidth estimator shown in (ref).
proofSee Appendix (ref).

\\

The adaptive random bandwidth (ARB) method is intuitive, robust and simple in practice, which can adapt to time-varying conditional distributions without a specific bandwidth or kernel function. A comparison of size performances of Wald tests using ARB with other competing methods are presented in Tables (ref), (ref) and (ref). We also find that updating $\bm{V_d}$ improves the size performance with use of $\alpha$ levels in the interquartile range around but not much for $\alpha$ levels like $1\%, 5\%$. More test results are presented in Appendix (ref) with changing sample size, quantile index and varying DGPs. The performance of ARB is robust. ARB can also be easily generalized to apply on multivariate quantile regressions, which is beyond the scope of this paper but in the interest of multivariate quantile regressions for future research. ARB also has the potential to achieve the second-order accuracy to Wald tests of nonlinear restrictions phillips1988formulation,de1993corrections in quantile regressions, which we would like to leave for future research.

table[table omitted — 1,322 chars of source]
table[table omitted — 1,314 chars of source]
table[table omitted — 1,139 chars of source]

Empirical Results

We study four US stock prices which are the Dow Jones Composite Average (DJCA), the NASDAQ 100 Index (NASDAQ100), the S&P 500, and the Wilshire 5000 Total Market Index (Will5000ind). We implement inference tests using the adaptive random bandwidth method with $n=1000$ and $\bm{V_d} = \bm{I}_{(p+1)\times (p+1)}$ which is not updated in simulations in this section. Each stock price time series has 2448 daily prices, ranging from 8th April 2010 to 30th December 2019. The price data were converted to return rates by multiplying 100 with the difference of the natural logarithm of the daily prices. The obtained return time series of each stock contains 2447 observations which of the last 400 observations are used for the out-of-sample testing after the first 2047 observations are used to estimate the model.

The $5\%$ 1-day VaRs of a return time series are the opposite conditional $5\%$ 1-day quantiles of this time series. There are four different CAViaR models considered in this section to model the conditional quantiles of the stock return time series. The $5\%$ 1-day VaRs are estimated via the four different CAViaR specifications and the estimation results are shown in Table (ref), (ref), (ref) and (ref) respectively. Each table contains the estimated parameters in a specified model, the corresponding standard errors obtained by the adaptive random bandwidth method with $n=1000$ and $\bm{V_d} = \bm{I}_{(p+1)\times (p+1)}$, the resulted two-sided p-values on parameter significance, the optimized value of the quantile regression objective function (RQ), the percentage of times the VaR is exceeded, and the p-values of dynamic quantile (DQ) tests, both in-sample and out-of-sample. The model estimations, the in-sample DQ tests as well as the out-of-sample DQ tests in this empirical study are set up in the same way of Section 6 of engle2004caviar. \[ \left\{

tabular[tabular omitted — 958 chars of source]

\right. \] The above four CAViaR specifications have been defined as the adaptive CAViaR, the symmetric absolute value CAViaR, the asymmetric slope CAViaR, and the indirect GARCH$(1,1)$ respectively in the Section 3 of engle2004caviar. In the implementation of the adaptive model in this emprical study, we follow engle2004caviar and set $G = 10$.

table[table omitted — 1,585 chars of source]
table[table omitted — 1,390 chars of source]
table[table omitted — 1,371 chars of source]
table[table omitted — 959 chars of source]

Comparing with the results in Section 6 of engle2004caviar, we can see the standard errors obtained by the adaptive random bandwidth method is much smaller relatively to the size of estimated parameters. We use significance level 5% to reject a parameter equal to zero as well as DQ tests. \lq\lq\, * \rq\rq\, denotes the rejections in Table (ref), (ref), (ref) and (ref). Each of the four models shows almost the same rejection results for the stock return time series. Remarkably, it is observed that the coefficient $\beta_1$ of the VaR autoregressive term is highly significant from zero in all the four models for each stock return time series. This further supports the standpoint of CAViaR specifications, confirming that the phenomenon of volatility clustering can be associated with the autoregressive VaR behaviour. The VaR exceedance in percentage indicates the realized risk level in applications. Dynamic quantile (DQ) tests based on the independence information regarding $\{ \text{Hit}_t\} $ are used to test model misspecification. We see a rejection in the in-sample DQ test on the symmetric absolute value model for the S&P500 but the realized VaR exceedances (in-sample and out-of-sample) are much close to 5% in Table (ref). So it can be complementary to judge CAViaR model specifications by looking at both VaR exceedances and inference tests like DQ tests.

In contract to the significance of $\beta_1$, the coefficient $\beta_2$ of $(y_{t-1})^{2}$ is insignificant in the indirect GARCH(1,1) model for all the stock return time series, see Table (ref). And the coefficient $\beta_2$ of $(y_{t-1})^{+}$ is insignificant in the asymmetric slope model, see Table (ref). Although the coefficient of $y_{t-1}$ is significant in the symmetric absolute model for all the stock return time series (see Table (ref)), it is mainly due to the significant explanatory role of $(y_{t-1})^{-}$ based on the results of the asymmetric slope model which the symmetric absolute model is nested in. The significance results of $\beta_1$ in the adaptive model for each stock return time series suggest that the 5% 1-day VaR can be associated with its 1-day lagged VaR violation which equals one if $y_{t-1}\leq f_{t-1}$ and zero otherwise. The significance results together implies that negative movements of a stock is significantly influential on its 5% 1-day VaR in the next day.

In terms of the model goodness of fit, we look at the RQ results. The asymmetric slope model presents the lowest RQ result for each stock return time series among the four models despite that it has the most coefficients.

Overall, all the four stock return time series present the same strong associations with the lagged 5% 1-day VaR in interpreting the present 5% 1-day VaR. The asymmetric slope model and the adaptive CAViaR are satisfying for all the four stock returns in terms of data interpretation and model performance concerns.

Conclusions

We found that the inference test performance in CAViaR models is not robust and unsatisfying due to the estimation of the conditional probability densities of time series. We found that the existing density estimation methods cannot fully adapt to time-varying conditional probability densities of CAViaR time series. So in this paper we have developed a method called adaptive random bandwidth which can robustly approximate the time-varying conditional probability densities of CAViaR time series by Monte Carlo simulations. This method not only avoids the haunting problem of choosing an optimal bandwidth but also ensures the reliability of CAViaR analysis based on the asymptotic normality of the model parameter estimator. In theory, our proposed method can be extended to general quantile regressions including multivariate cases easily and robustly. This method also has the potential to achieve the second-order accuracy to Wald tests of nonlinear restrictions phillips1988formulation,de1993corrections in quantile regressions.