An interval-valued GARCH model for range-measured return processes

TL;DR

This paper introduces an interval-valued GARCH model that captures the information in range-measured returns, improving volatility estimation and prediction over traditional GARCH models.

Contribution

It proposes a novel interval-valued GARCH model for range data, with theoretical properties and empirical validation on financial data.

Findings

Int-GARCH outperforms traditional GARCH in volatility prediction

Model effectively captures information in range-measured returns

Parameters estimated with maximum likelihood, with established asymptotic properties

Abstract

Range-measured return contains more information than the traditional scalar-valued return. In this paper, we propose to model the [low, high] price range as a random interval and suggest an interval-valued GARCH (Int-GARCH) model for the corresponding range-measured return process. Under the general framework of random sets, the model properties are investigated. Parameters are estimated by the maximum likelihood method, and the asymptotic properties are established. Empirical application to stocks and financial indices data sets suggests that our Int-GARCH model overall outperforms the traditional GARCH for both in-sample estimation and out-of-sample prediction of volatility.

Tables

Table 1. Table 1: Average result of 100 repeated simulations from two Int-GARCH (1,1,1) models (I, II) and two Int-ARCH (1,1) models (III, IV).

			Mean	Mean	Empiriccal	Asymptotic
Model	Parameters		Estimate	Absolute Error	Standard Error	Standard Error
I	$k$	1.8147	1.8081	0.077	0.1061
	$\mu$	0.0906	0.0887	0.0072	0.0086	0.0082
	${\alpha }_1$	0.0318	0.0284	0.0184	0.0235	0.0211
	${\beta }_1$	0.374	0.3586	0.0171	0.0148	0.0149
	${\gamma }_1$	0.1265	0.1269	0.0314	0.0381	0.0314
II	$k$	1.2134	1.2251	0.0412	0.05
	$\mu$	0.071	0.0712	0.0068	0.0086	0.0094
	${\alpha }_1$	0.1833	0.1784	0.025	0.0318	0.0306
	${\beta }_1$	0.2334	0.2345	0.0152	0.0192	0.0207
	${\gamma }_1$	0.1732	0.174	0.0467	0.0579	0.0602
III	$k$	1.5139	1.5045	0.04	0.0506
	$\mu$	0.074	0.0738	0.0026	0.0034	0.0036
	${\alpha }_1$	0.037	0.0334	0.0185	0.0228	0.0208
	${\beta }_1$	0.3436	0.3426	0.0139	0.017	0.0174
	${\gamma }_1$	0
IV	$k$	1.3632	1.3621	0.038	0.0486
	$\mu$	0.0584	0.0593	0.0029	0.0037	0.0033
	${\alpha }_1$	0.1927	0.1962	0.0208	0.026	0.0269
	${\beta }_1$	0.322	0.3253	0.0161	0.0197	0.0195
	${\gamma }_1$	0

Table 2. Table 2: Fitted Int-(G)ARCH models for the period 2006-2011. For the variance parameters, the number in the parenthesis to the right of the estimate is the associated standard error based on the asymptotic normality.

			Parameter Estimates
			k	$\mu$	${\alpha }_1$	${\beta }_1$	${\gamma }_1$
Stocks	AAPL		1.5694	0.0043 (0.0004)	0.024 (0.0242)	0.4948 (0.0233)	0
	AXP		1.7676	0.0012 (0.0003)	0	0.5357 (0.0147)	0
	BA		1.6245	0.0031 (0.0005)	0.0616 (0.0239)	0.4820 (0.0258)	0
	BAC		1.511	0.00003 (0.0003)	0.0374 (0.0266)	0.6507 (0.0214)	0
	DD		1.7335	0.0019 (0.0004)	0.0034 (0.0240)	0.5078 (0.0226)	0
	JPM		1.7228	0.0006 (0.0004)	0.0139 (0.0203)	0.3399 (0.0096)	0.3894 (0.0354)
	KO		1.8602	0.0014 (0.0002)	0.0336 (0.0227)	0.4480 (0.0212)	0
	MSFT		1.7135	0.0029 (0.0004)	0.0708 (0.0257)	0.4485 (0.0238)	0
	T		1.8338	0.0022 (0.0003)	0	0.4510 (0.0181)	0
	WMT		1.8777	0.0017 (0.0003)	0	0.4540 (0.0200)	0
Indices	DJI		1.7755	0.0009 (0.0002)	0	0.5051 (0.0172)	0
	SPX		1.6472	0.0008 (0.0002)	0	0.5383 (0.0348)	0.0252 (0.0586)
	FTSE		1.7529	0.0008 (0.0002)	0	0.5100 (0.0183)	0
	CAC		1.4469	0.0014 (0.0003)	0	0.5800 (0.0382)	0.0265 (0.0608)

Table 3. Table 3: In-sample comparisons of Int-GARCH (1,1,1) and GARCH (1,1). Larger values of R 2 superscript 𝑅 2 R^{2} and smaller values of QLIKE and HMSE indicate higher ranking of the volatility model.

		${\text{R}}^2$		QLIKE		HMSE
		GARCH	Int-GARCH	GARCH	Int-GARCH	GARCH	Int-GARCH
Stocks	AAPL	0.4707	0.3822	-7.0084	-7.0484	0.802	0.5178
	AXP	0.5099	0.5692	-7.0235	-7.0851	0.36	0.3431
	BA	0.4949	0.5446	-7.2523	-7.2592	0.5623	0.4291
	BAC	0.5179	0.4889	-6.6973	-6.7652	1.104	0.665
	DD	0.4637	0.5778	-7.29	-7.3093	0.5374	0.4178
	JPM	0.3728	0.4608	-6.8266	-6.9093	0.6181	0.4411
	KO	0.4512	0.5287	-8.2397	-8.242	0.5816	0.3866
	MSFT	0.3831	0.4323	-7.4539	-7.4733	0.8664	1.1465
	T	0.4099	0.4522	-7.5155	-7.5569	1.5514	0.7272
	WMT	0.3644	0.4009	-7.8462	-7.8861	0.7694	0.7242
Indices	DJI	0.5006	0.4988	-8.5019	-8.5877	0.5695	0.5428
	SPX	0.5209	0.5286	-8.401	-8.4873	0.481	0.4839
	FTSE	0.4769	0.538	-8.4324	-8.5158	0.4122	0.359
	CAC	0.4171	0.4988	-8.1337	-8.1737	0.4158	0.3854

Table 4. Table 4: Comparisons of 1-step-ahead predictions of Int-GARCH (1,1,1) and GARCH (1,1).

1-step-ahead Prediction
		${\text{R}}^2$		QLIKE		HMSE
		GARCH	Int-GARCH	GARCH	Int-GARCH	GARCH	Int-GARCH
Stocks	AAPL	0.0485	0.298	-7.5694	-7.6838	0.5426	0.4262
	AXP	0.2591	0.5507	-7.3522	-7.4218	0.3678	0.248
	BA	0.2909	0.4847	-7.4518	-7.5054	0.3293	0.2911
	BAC	0.2461	0.4203	-6.6467	-6.7399	0.4693	0.443
	DD	0.3053	0.5058	-7.3962	-7.4448	0.3275	0.2614
	JPM	0.3587	0.4318	-7.0733	-7.1648	0.392	0.3583
	KO	0.2844	0.4279	-8.3893	-8.4008	0.3503	0.445
	MSFT	0.1032	0.1894	-7.6106	-7.6581	0.7244	0.7286
	T	0.153	0.2121	-8.1725	-8.1884	0.4431	0.6729
	WMT	0.2538	0.4033	-6.8641	-8.3991	0.7947	0.3363
Indices	DJI	0.3204	0.5334	-8.5295	-8.6221	0.3875	0.3421
	SPX	0.3413	0.5882	-8.3776	-8.461	0.3984	0.3436
	FTSE	0.3421	0.5555	-8.3617	-8.3958	0.4161	0.3795
	CAC	0.3673	0.5285	-7.8501	-7.8694	0.4095	0.3711

Table 5. Table 5: Comparisons of 2-step-ahead predictions of Int-GARCH (1,1,1) and GARCH (1,1).

2-step-ahead Prediction
		${\text{R}}^2$		QLIKE		HMSE
		GARCH	Int-GARCH	GARCH	Int-GARCH	GARCH	Int-GARCH
Stocks	AAPL	0.0132	0.1574	-7.5383	-7.5879	0.6343	0.4554
	AXP	0.1275	0.2865	-7.3106	-7.3747	0.635	0.3224
	BA	0.2132	0.3291	-7.4326	-7.4655	0.3797	0.3143
	BAC	0.1367	0.1879	-6.5818	-6.6704	0.8272	0.6055
	DD	0.2097	0.3287	-7.3671	-7.4149	0.4624	0.2928
	JPM	0.2601	0.2432	-7.0329	-7.1009	0.6407	0.6817
	KO	0.1908	0.2407	-8.3579	-8.3372	0.4782	0.6853
	MSFT	0.0596	0.1203	-7.5891	-7.6448	0.7986	0.5681
	T	0.0888	0.1128	-8.1397	-8.1643	0.6592	0.4805
	WMT	0.1925	0.2499	-6.7841	-8.3713	0.8143	0.347
Indices	DJI	0.2046	0.3392	-8.4939	-8.5731	0.525	0.3758
	SPX	0.2061	0.3736	-8.3409	-8.4095	0.5378	0.3789
	FTSE	0.2581	0.385	-8.3375	-8.3469	0.4742	0.4309
	CAC	0.2666	0.3509	-7.8209	-7.8204	0.4894	0.4406

Table 6. Table 6: Comparisons of 5-step-ahead (one week) predictions of Int-GARCH (1,1,1) and GARCH (1,1).

5-step-ahead Prediction
		${\text{R}}^2$		QLIKE		HMSE
		GARCH	Int-GARCH	GARCH	Int-GARCH	GARCH	Int-GARCH
Stocks	AAPL	0.002	0.0141	-7.4968	-7.4147	0.6861	0.5453
	AXP	0.0278	0.0488	-7.2142	-7.2724	1.6016	0.5616
	BA	0.1281	0.1872	-7.3887	-7.4081	0.5966	0.3523
	BAC	0.0625	0.0189	-6.401	-6.4962	3.9111	1.4395
	DD	0.1072	0.1175	-7.2938	-7.3466	0.946	0.4013
	JPM	0.1679	0.0986	-6.9515	-7.0213	1.2665	0.7303
	KO	0.0875	0.0743	-8.2977	-8.2983	1.2164	1.0582
	MSFT	0.0336	0.0416	-7.551	-7.5808	1.0182	0.7925
	T	0.052	0.0717	-8.0929	-8.0853	0.95	0.4313
	WMT	0.1068	0.0463	-6.5603	-8.3041	0.8564	0.4626
Indices	DJI	0.1092	0.1266	-8.4081	-8.4504	1.1346	0.6644
	SPX	0.1083	0.1302	-8.2571	-8.2834	1.1027	0.6043
	FTSE	0.1334	0.138	-8.2455	-8.2624	1.154	0.5915
	CAC	0.1696	0.1709	-7.7553	-7.7696	0.7602	0.5749

Equations

$$r_t$$

$$A+B=\left\{a+b:a\in A,b\in B\right\},\lambda A=\left\{\lambda a:a\in A\right\},$$

$$s_A\left(u\right)=\underset{a\in A}{\sup }⟨u,a⟩,\forall u\in S^{d-1},$$

$${\rho }_2(A,B)={\left[d{\int }_{S^{d-1}}|s_A(u)-s_B(u){|}^2\mu (du)\right]}^{\frac{1}{2}}.$$

$$\text{E}(X)=\left\{\text{E}\xi :\xi \in Xa.s.,\text{E}\Vert \xi \Vert <\infty \right\}.$$

$$\text{E}{\rho }^2(X,{\text{E}}_F(X))=\underset{A\in {\mathcal{K}}_{\mathcal{C}}}{\inf }\text{E}{\rho }^2(X,A).$$

$$\text{Var}(X)=\text{E}{\rho }_2^2(X,\text{E}(X)).$$

$$\text{E}(X)=\left[\text{E}(X^C)-\text{E}(X^R),\text{E}(X^C)+\text{E}(X^R)\right],$$

$${\rho }_2(x,y)={\left[(x^C-y^C{)}^2+(x^R-y^R{)}^2\right]}^{\frac{1}{2}}.$$

$$\text{Var}(X)=\text{E}\left[(X^C-\text{E}(X^C){)}^2+(X^R-\text{E}(X^R){)}^2\right]=\text{Var}(X^C)+\text{Var}(X^R).$$

$$\text{E}(X|\mathcal{F})=\left[\text{E}(X^C|\mathcal{F})-\text{E}(X^R|\mathcal{F}),\text{E}(X^C|\mathcal{F})+\text{E}(X^R|\mathcal{F})\right].$$

$$\text{Var}(X|\mathcal{F})$$

$$r_t=[{\lambda }_t-{\delta }_t,{\lambda }_t+{\delta }_t],t=1,2,\cdots ,T.$$

$$\text{E}(r_t)$$

$$\text{Var}(r_t)$$

$$H_t^2=\text{Var}(r_t|{\mathcal{F}}_{t-1})=\text{Var}({\lambda }_t|{\mathcal{F}}_{t-1})+\text{Var}({\delta }_t|{\mathcal{F}}_{t-1}).$$

$$H_t^2=g\left(H_s^2,{\hat{H}}_s^2(r_s):s\le t-1\right),$$

$$H_t^2=\mu +\sum _{i=1}^p{\alpha }_i\left[{\lambda }_{t-i}^2+{\delta }_{t-i}^2-c^2\right]+\sum _{i=1}^q{\beta }_i{H}_{t-i}^2,$$

$$H_t^2=\mu +\sum _{i=1}^p{\alpha }_i\left[{\lambda }_{t-i}^2+{\delta }_{t-i}^2\right]+\sum _{i=1}^q{\beta }_i{H}_{t-i}^2.$$

$$r_t$$

$$v_t$$

$${ϵ}_t$$

$${\eta }_t$$

$$h_t$$

$$H_t^2=\text{Var}(h_t{ϵ}_t|{\mathcal{F}}_{t-1})+\text{Var}(h_t{\eta }_t|{\mathcal{F}}_{t-1})=h_t^2(1+k).$$

$$\text{Var}\left({ϵ}_t|{\mathcal{F}}_{t-1}^{ϵ,\eta }\right)=1,$$

$$\text{Var}\left({\eta }_t|{\mathcal{F}}_{t-1}^{ϵ,\eta }\right)=k.$$

$${\rho }_W^2\left(x,y\right)={\left(x^C-y^C\right)}^2+{\left(x^R-y^R\right)}^2{\int }_{[0,1]}{\left(2\lambda -1\right)}^{2}dW(\lambda ),x,y\in {\mathcal{K}}_{\mathcal{C}}(\mathbb{R}),$$

$$r_t(\omega )=h_t\left({ϵ}_t+\omega {\eta }_t\right),\omega \in [-1,1].$$

$${\sigma }_t^2\left(\omega \right)=\text{Var}\left(r_t(\omega )|{\mathcal{F}}_{t-1}\right)=\left(1+{\omega }^{2}k\right)h_t^2.$$

$${\sigma }_t^2=\frac{{\int }_{-1}^1{\sigma }_t^2(\omega )d(\omega )}{{\int }_{-1}^{1}d(\omega )}=\frac{1}{2}{\int }_{-1}^1\left(1+{\omega }^{2}k\right)h_t^{2}d(\omega )=\left(1+\frac{1}{3}k\right)h_t^2.$$

$${\hat{h}}_t\left(1\right)=\mu +\sum _{i=1}^p{\alpha }_i|{\lambda }_{t+1-i}|+\sum _{i=1}^q{\beta }_i{\delta }_{t+1-i}+\sum _{i=1}^w{\gamma }_i{h}_{t+1-i}.$$

$$h_{t+2}$$

Full text

Keywords: GARCH; volatility; price range; random sets; random interval; interval-valued time series.

1 Introduction

The issue of assets’ volatility plays an essential role in modern finance. It provides a measured variability for the asset price over a certain period of time, and is a key parameter in many financial applications such as financial derivatives pricing, risk assessment, and portfolio management. The squared return of log prices, being the classical estimator of variance, has often been used as the “ideal” proxy of volatility. As a result, many volatility models built on returns were proposed and, for a long time, have been very popular and successful. The celebrated ARCH ([9]) and GARCH ([4]) models are examples of this type, which use information of traditional return series. Recently, as the high-frequency transaction data become widely available, price changes can be practically monitored in a continuous way, and the traditional low-frequency (e.g., daily) return is no longer quite representative of the volatility. For example, a small low-frequency return does not necessarily imply low volatility, as the price may fluctuate a lot and close at a similar level to the opening. On the other hand, a big return could only be the result of a very different opening price from the previous day’s closing. Therefore, the traditional return-based models, using lesser information, are likely to produce inefficient or even incorrect estimates of the volatility.

In fact, since closing price is only a “snapshot” among numerous prices during a day, there is nothing “special” about it and return need not be calculated solely based on it. With the wide availability of high-frequency data, more information can obviously be utilized. For example, the difference between any two observed (log) prices in two consecutive days can be viewed as a proxy of volatility. This idea leads us to propose a naturally generalized concept of daily return, which is an interval that includes all the “snapshot” returns. Let $y_{t}\left(s\right)$ be the log price of an asset at time $s$ on day $t$. Ideally $s$ should be a continuous time index. But since price can only be observed at discrete times even for high-frequency data, $s$ is assumed to be a discrete index here. We define the interval-valued daily return as

[TABLE]

Namely, $r_{t}$ is the range of “snapshot” returns during one day, which is expected to contain richer information about the daily volatility than the single closing-to-closing return. Our goal is to build a volatility model that reveals the dynamics of the entire return range $r_{t}$ as a whole. To this end, we propose to extend the GARCH model to allow for interval-valued return input, and the resulting interval-valued model is called the Int-GARCH model.

The essential idea of our Int-GARCH model is to enrich the return based volatility model (i.e., GARCH) with range information. Such an idea is not new. There has been a great deal of effort in the literature on volatility modeling using price range of either low-frequency or high-frequency data ([12], [27], [29], [20], [1], [6], [10], [5], [7]). In these methods, price range is essentially viewed as an exogenous variable and is included into the modeling via functions of certain forms. In a more systematic way, our approach integrates the level and range information by modeling the [low, high] return range $r_{t}$ as a whole in the framework of random interval. As is shown in Section 3, based on the fitted model, our final prediction of volatility is made such that all the information contained in $r_{t}$ is systematically accounted for.

In addition to volatility prediction, our Int-GARCH model in a broader sense makes a pioneer contribution to the study of conditional heteroscedasticity for interval-valued time series. It has been a while since interval-valued time series was investigated and there are a handful of results mainly on the mean models, i.e., models that aim at making prediction of the interval-valued mean, and their applications (e.g.,[22], [14]). The focus of our model, by comparison, is the conditional variance. Variance model for interval-valued time series, as far as the authors are aware, largely remains an unexploited area, and our study is among the very first investigations.

Our theoretical results are two-fold. We first establish that under certain conditions our Int-GARCH model achieves weak stationarity that is characterized by a time invariant mean and variance. Then, under the assumption of weak stationarity, we define and obtain the explicit formula of the autocorrelation function (ACF) of the Int-GARCH process. We propose to estimate the model parameters by the method of maximum likelihood and provide the associated asymptotic properties. Simulation shows that the results are consistent with our theoretical findings. For empirical study, we analyze several stocks and indices data that are representative of the market. Our Int-GARCH model is compared to GARCH using the RV as the market proxy. Based on both in-sample and out-of-sample comparisons, our Int-GARCH model overall outperforms GARCH with higher correlations and reduced errors to RV.

The rest of the paper is organized as follows. Random sets preliminaries are provided in Section 2. We formally introduce our Int-GARCH model and its general properties in Section 3. Stationarity of Int-GARCH (1,1,1) is presented in Section 4. Section 5 discusses the maximum likelihood estimator for the model parameters and carefully investigates its performances by a simulation study. Empirical study with the stocks and indices data, as well as a detailed discussion, are reported in Section 6. We finish with concluding remarks in Section 7. Proofs and useful lemmas are provided in the Appendix.

2 Preliminaries of random sets

Throughout the paper, we will view $r_{t}$ as a random interval and model its dynamics under the framework of random sets. To facilitate our presentation, we briefly introduce the basic notations and definitions in the random set theory. For more details we refer the readers to [16], [23], [2], [25], [32], among others.

Let $(\Omega,\mathcal{L},P)$ be a probability space. Denote by $\mathcal{K}\left(\mathbb{R}^{d}\right)$ or $\mathcal{K}$ the collection of all non-empty compact subsets of $\mathbb{R}^{d}$. In the space $\mathcal{K}$, a linear structure is defined by Minkowski addition and scalar multiplication, i.e.,

[TABLE]

$\forall A,B\in\mathcal{K}$ and $\lambda\in\mathbb{R}$. A random compact set is a Borel measurable function $A:\Omega\rightarrow\mathcal{K}$, $\mathcal{K}$ being equipped with the Borel $\sigma$-algebra induced by the Hausdorff metric. For each $A\in\mathcal{K}\left(\mathbb{R}^{d}\right)$, the function defined on the unit sphere $S^{d-1}$:

[TABLE]

is called the support function of A. If $A(\omega)$ is convex almost surely, then $A$ is called a random compact convex set. Especially, a one-dimensional random compact convex set is called a random interval. Let $\mathcal{K}_{\mathcal{C}}(\mathbb{R}^{d})$ denote the space of non-empty compact convex subsets of $\mathbb{R}^{d}$. An $L_{2}$ metric in $\mathcal{K}_{\mathcal{C}}(\mathbb{R}^{d})$ is given via the support function by

[TABLE]

Much of the random sets theory has focused on compact convex sets via their support functions. For a random compact convex set $X$, the well accepted expectation is defined by the Aumann integral of set-valued function ([3]) as

[TABLE]

Alternatively, [11] gave a general definition for the expectation of a random element $X$ in the metric space $(\mathcal{K}_{\mathcal{C}}(\mathbb{R}^{d}),\rho)$ as the solution of

[TABLE]

If we choose the metric $\rho$ to be $\rho_{2}$, then the Fréchet expectation $\text{E}_{F}(X)$ coincides with the Aumann expectation $\text{E}(X)$, and the variance is further defined as

[TABLE]

See [21], [17], and [18]. Restricting to $\mathcal{K}_{\mathcal{C}}(\mathbb{R})$, the Aumann expectation of a random interval $X$ is simply

[TABLE]

where $(\cdot)^{C}$ and $(\cdot)^{R}$ denote the corresponding center and radius, respectively. The $\rho_{2}$ distance between two intervals $x$ and $y$ is ([17])

[TABLE]

Thus, the variance of a random interval is consequently calculated as

[TABLE]

Let $\mathcal{F}$ be any $\sigma$-filed on $\Omega$. The conditional expectation of a random set given $\mathcal{F}$ is defined by the integral and conditional expectation of multivalued functions ([15]). Specifically for a random interval $X$,

[TABLE]

Based on the conditional expectation, the conditional variance of a random set is defined as the conditional mean squared distance from its conditional expectation ([26]). According to the definition, the conditional variance of $X$ with respect to the distance $\rho_{2}$ is given by

[TABLE]

3 The Int-GARCH Model

3.1 Model specification

We assume observing a range-measured return series $\left\{r_{t}\right\}_{t=1}^{T}$ of the form

[TABLE]

That is, $\left\{\lambda_{t}\right\}_{t=1}^{T}$ and $\left\{\delta_{t}\right\}_{t=1}^{T}$ are the associated center and radius processes, both of which are observable. According to (5) and (7), the mean and variance of $r_{t}$ as a random interval are

[TABLE]

Let $\mathcal{F}_{t}$ denote the information set up to time t, i.e. $\mathcal{F}_{t}=\sigma\left\{r_{s}:s\leq t\right\}.$ We are concerned with the conditional variance $H_{t}^{2}$ of $r_{t}$ given $\mathcal{F}_{t-1}$. According to (9), $H_{t}^{2}$ with respect to $\rho_{2}$ is computed as

[TABLE]

The GARCH model depicts the conditional variance of a point-valued return process as a linear function of the past squared returns and variances. This was inspired by the fact that assets returns usually exhibit volatility clustering: large variations in prices tend to cluster together, resulting in separate dynamic and tranquil periods of the market. Extending this spirit to the interval-valued process $\left\{r_{t}\right\}$, one would expect, conceptually, a model like

[TABLE]

where $\hat{H}_{s}^{2}(r_{s})$ denotes a return range based proxy for $H_{s}^{2}$, and $g$ is linear in $H_{s}^{2}$ and $\hat{H}_{s}^{2}(r_{s})$. Assuming $\text{E}(\lambda_{t}|\mathcal{F}_{t-1})=0$ and $\text{E}(\delta_{t}|\mathcal{F}_{t-1})=c>0$, a reasonable function $g$ in (12) seems to imply

[TABLE]

where $p>0,q\geq 0$. The constant $c$ can be absorbed into the parameter $\mu$ and the above equation is simplified to

[TABLE]

To give more flexibility to our model, we allow for different degrees of dependence of $H_{t}$ on the past centers and radii. In addition, adopting the idea of more robust modeling of volatility by [33] and [30], we propose to model the conditional standard deviation $H_{t}$ directly, instead of via the conditional variance $H_{t}^{2}$ ([8] also considered such a specification as a special case of their A-PARCH model).

Given the above discussion, our Int-GARCH $(p,q,w)$ model for the return range process is specified as

[TABLE]

where $p>0,q>0,w\geq 0$, and $\left\{\alpha_{i}:i=1,\cdots,p\right\},\left\{\beta_{i}:i=1,\cdots q\right\},\left\{\gamma_{i}:i=1,\cdots,w\right\}$ are positive constants. In addition, the error terms $\epsilon_{t}$ and $\eta_{t}$ are assumed to be independent. In (13), “$\cdot$” denotes the scalar multiplication. Under this specification, the conditional variance of $r_{t}$ is seen to be

[TABLE]

Although we impose parametric assumptions on the random errors $\epsilon_{t}$ and $\eta_{t}$ to simplify our presentation here, they are not really necessary. In practice, it is best to use the true data generating distributions, which vary from data to data. So, in replacement of (15)-(16), a relaxed yet sufficient specification for $\epsilon_{t}$ and $\eta_{t}>0$ is

[TABLE]

Remark 1**.**

A more general metric for $\mathcal{K}_{\mathcal{C}}(\mathbb{R})$ was proposed by [13], which essentially takes the form

[TABLE]

where $W$ is any non-degenerate symmetric measure on $[0,1]$. Compared to the $\rho_{2}$ metric in (6), the flexibility of $\rho_{W}$ lies in its choice of a weight between the center and radius. Since we use different parameters $\alpha$’s and $\beta$’s for the center and radius, respectively, this flexibility is in fact accounted for in our specification of $h_{t}$ in (17).

3.2 Volatility forecasting

The $h_{t}$ in the model (13)-(17) relates to the conditional standard deviation $H_{t}$ of the return range $r_{t}$ according to equation (18). However, $H_{t}$ is not exactly the daily volatility as in the literature. To see the relation between $h_{t}$ and the daily volatility $\sigma^{2}_{t}$, we notice that $r_{t}$ contains all possible returns of day $t$, using different prices during days $t$ and $t-1$. Consider an arbitrary return at position $\omega$ in $r_{t}$, denoted by $r_{t}(\omega)$, as the volatility proxy, i.e.

[TABLE]

Volatility based on $r_{t}(\omega)$ is calculated as

[TABLE]

Our Int-GARCH volatility $\sigma^{2}_{t}$ is defined as the average of $\left\{\sigma^{2}_{t}(\omega):\omega\in[-1,1]\right\}$. Assuming equal weight for each point return, $\sigma^{2}_{t}$ can be calculated as

[TABLE]

Therefore, volatility forecast is essentially made by predicting $h_{t}$. The 1-step-ahead prediction of $h_{t}$ is immediately defined by equation (17) as

[TABLE]

To make a 2-step-ahead prediction, notice that

[TABLE]

Replacing $|\epsilon_{t+1}|$ and $\eta_{t+1}$ with their expectations, and $h_{t+1}$ with $\hat{h}_{t}\left(1\right)$, we get

[TABLE]

The general $l$-step-ahead prediction can be calculated recursively. In particular, the formula for Int-GARCH (1,1,1) model is given by

[TABLE]

3.3 Mean stationarity

We provide the necessary and sufficient conditions of mean stationarity for the Int-GARCH (p,q,w) model in the following theorem.

Theorem 1**.**

Consider the general Int-GARCH model (13)-(17). Define

[TABLE]

and

[TABLE]

where $i=1,2,\cdots,m=\max\left\{p,q,w\right\}$. Assume $\left\{r_{t}\right\}$ starts from its infinite past with a finite mean. Then, $Eh_{t}<\infty$ if and only if $\sum_{i=1}^{k}\mu_{i}<1$. When this condition is satisfied,

[TABLE]

and

[TABLE]

4 Stationarity of Int-GARCH (1,1,1)

Similar to the traditional GARCH model, the Int-GARCH(1,1,1) process is a simple but effective model for analyzing interval-valued time series with conditional heteroskedasticity. In this section, we derive several important distributional properties of Int-GARCH (1,1,1). Before we present our theoretical results, we first notice that for the Int-GARCH (1,1,1) process,

[TABLE]

Defining the i.i.d. random variables $x_{t}=\alpha_{1}|\epsilon_{t-1}|+\beta\eta_{t-1}+\gamma_{1}$, $t\in\mathbb{N}$, $h_{t}$ can be re-written as

[TABLE]

We will use (28) throughout this section to simplify notations.

4.1 Weak stationarity

It is derived in [17] that the covariance between two random intervals with respect to the $\rho_{2}$ metric is the sum of the covariances between the two centers and two radii. This implies

[TABLE]

We are ready to extend the notion of weak stationarity to interval-valued time series in the obvious way.

Definition 1**.**

An interval-valued time series $\left\{r_{t}\right\}$ is said to be weakly stationary, or second-moment stationary, if its unconditional mean $\text{E}\left(r_{t}\right)$ and covariance $\text{Cov}\left(r_{t},r_{t+s}\right)$ exist and are independent of time $t$ for all integers $s$, where $\text{E}\left(r_{t}\right)$ and $\text{Cov}\left(r_{t},r_{t+s}\right)$ are given in (10) and (29), respectively.

The existences of $\text{E}\left(r_{t}\right)$ and $\text{Var}\left(r_{t}\right)$ are closely related to those of $\text{E}\left(h_{t}\right)$ and $\text{E}\left(h_{t}^{2}\right)$, respectively. In fact, $\text{E}(h_{t}^{2})<\infty$ implies the existences of the first two moments of $r_{t}$. We give precise conditions in the following two theorems.

Theorem 2**.**

Consider the Int-GARCH model (13)-(17) with $p=q=w=1$. Assume $\left\{r_{t}\right\}$ starts from its infinite past with a finite mean. Then, $\text{E}(h_{t})<\infty$ if and only if $\text{E}(x_{t})<1$, i.e.

[TABLE]

When this condition is satisfied,

[TABLE]

and

[TABLE]

Theorem 3**.**

Consider the Int-GARCH(1,1,1) model $\left\{r_{t}\right\}$ as in Theorem 2. Assume $\left\{r_{t}\right\}$ starts from its infinite past with a finite variance. Then $E\left(h_{t}^{2}\right)<\infty$ if and only if $E\left(x_{t}^{2}\right)<1$, i.e.

[TABLE]

When this condition is satisfied,

[TABLE]

and

[TABLE]

where $E\left(h_{t}\right)$ is given in (30), and $C_{1}=E\left(x_{t}\right)$, $C_{2}=E\left(x_{t}^{2}\right)$.

So far we have found equivalent conditions for the existence of the unconditional mean and variance. In order to achieve weak stationarity, we still need to find conditions under which the (unconditional) covariances are finite and time-invariant. According to the following Theorem 4, these conditions turn out to be the same as those for the existence of variance. This is not surprising as

[TABLE]

We summarize this conclusion in Corollary 1 following Theorem 4.

Theorem 4**.**

Consider the Int-GARCH(1,1,1) process $\left\{r_{t}\right\}$. Under the assumptions of Theorem 3, the covariance of any two random intervals $r_{t}$ and $r_{t+s}$ is given by

[TABLE]

where $E\left(h_{t}\right)$ and $E\left(h_{t}^{2}\right)$ are given in (30) and (32), respectively, and $E\left(h_{t}h_{t+s}\eta_{t}\right)$ is calculated explicitly in Lemma 1 (see Appendix).

Corollary 1**.**

The Int-GARCH(1,1,1) process is weakly stationary, or second-moment stationary, if and only if $E\left(x_{t}^{2}\right)<1$.

4.2 Auto-correlation function (ACF)

The notion of the variance and covariance for compact convex random sets were naturally extended to the correlation coefficient of two random sets $A$ and $B$, which is defined as

[TABLE]

Based on this definition, we calculate the auto-correlation function (ACF) of the Int-GARCH(1,1,1) process and give the result in the corollary below.

Corollary 2**.**

Under the assumptions of Theorem 3, the auto-correlation function of the Int-GARCH(1,1,1) process $\left\{r_{t}\right\}$ is

[TABLE]

We plot the ACF for a specific Int-GARCH(1,1,1) model (Model I in the simulation) in Figure 1. We see that the centers are uncorrelated. This has been verified by (42) in the proof of Theorem 4. The radii, or the lengths of the intervals, have a relatively persistent auto-correlation, which coincides with the phenomenon of “volatility clustering”. This long-term dependence of radii carries over to the intervals as a whole, and results in a slow-dying ACF of the interval-valued process.

5 Parameter Estimation

In this section, we develop parameter estimation of our Int-GARCH model and the asymptotic properties of the estimators. There are two groups of parameters: the error distribution parameter $k>0$ and the variance parameters

[TABLE]

We first provide a sample estimate for $k$, which is shown to be strongly consistent. Then the variance parameters $\boldsymbol{\theta}$ are estimated by the method of maximum likelihood, and the asymptotic normality for the MLE is established. A simulation study is presented that shows the empirical performance of the estimators.

5.1 Two-stage estimation of $k$ and $\boldsymbol{\theta}$

The error distribution parameter $k$ can be conveniently estimated by the method of moments. Notice that

[TABLE]

and consequently,

[TABLE]

Replacing $E\left(\delta_{t}\right)$ and $E\left(|\lambda_{t}|\right)$ in (35) by their sample estimates, we obtain the method of moments estimator of $k$ as

[TABLE]

Under the condition that $\left\{h_{t}\right\}$ is strictly stationary and ergodic, and additionally that $E\left(h_{t}\right)<\infty$, the sample means $\overline{\delta_{t}}$ and $\overline{|\lambda_{t}|}$ are both strongly consistent, and so is $\hat{k}$.

Given $k$, the conditional likelihood function of $\left\{r_{t},t=1,\cdots,T\right\}$ is

[TABLE]

Thus, the conditional log-likelihood function up to a constant is

[TABLE]

where $h_{t}=\mu+\sum\alpha_{i}|\lambda_{t-i}|+\sum\beta_{i}\delta_{t-i}+\sum\gamma_{i}h_{t-i}$. Then, the maximum likelihood estimate $\hat{\boldsymbol{\theta}}_{T}$ is defined as

[TABLE]

and it can be computed by the scoring algorithm

[TABLE]

Under the condition of strict stationarity and ergodicity, and some moment requirement for $\left\{h_{t}\right\}$, $\hat{\boldsymbol{\theta}}_{T}$ is consistent and asymptotically normal. We state the result in the following theorem. Details of the proof are deferred to the Appendix.

Theorem 5**.**

*Assume the process $\left\{h_{t}\right\}$ is strictly stationary and ergodic, and $E\left(h_{t}^{2}\right)<\infty$. Assume in addition that the parameter space $\Theta$ is compact. Then, the maximum likelihood estimator defined in (36) satisfies:

(i) $\hat{\boldsymbol{\theta}}_{T}\stackrel{{\scriptstyle\mathcal{P}}}{{\to}}\boldsymbol{\theta}$ as $T\to\infty$;

(ii) $T^{\frac{1}{2}}\left(\hat{\boldsymbol{\theta}}_{T}-\boldsymbol{\theta}_{0}\right)\stackrel{{\scriptstyle\mathscr{D}}}{{\to}}N\left(\boldsymbol{0},\boldsymbol{I}^{-1}(\boldsymbol{\theta}_{0})\right)$, where $\boldsymbol{I}(\boldsymbol{\theta}_{0})=-E\left[\nabla^{2}l_{t}(\boldsymbol{\theta}_{0})\right]$ is the Fisher information matrix evaluated at $\boldsymbol{\theta}_{0}$.*

As an immediate consequence of the theorem, the asymptotic covariance matrix of $\hat{\boldsymbol{\theta}}_{T}$ is $\frac{1}{T}\boldsymbol{I}^{-1}(\boldsymbol{\theta}_{0})$. It can be consistently estimated by the inverse Hessian $-\left[\nabla^{2}l\left(\hat{\boldsymbol{\theta}}_{T}\right)\right]^{-1}$, which is easily obtained from the scoring algorithm. Especially for Int-GARCH(1,1,1), the condition for strict stationarity and ergodicity is the same as that for mean stationarity, that is, $E(x_{t})<\infty$, and consequently $E(h_{t})<\infty$ (see Lemma 3 in the Appendix). So, assuming mean stationarity, the method of moments estimator $\hat{k}$ is strongly consistent. If we further assume covariance stationarity, i.e., $E(x_{t}^{2})<\infty$ in view of Corollary 1, the maximum likelihood estimator $\hat{\boldsymbol{\theta}}$ is consistent and asymptotically normal. We summarize the conclusion in the following corollary.

Corollary 3**.**

*Consider the Int-GARCH(1,1,1) model.

(i) If $E(x_{t})<\infty$, then $\hat{k}$ is strongly consistent;

(ii) If in addition $E(x_{t}^{2})<\infty$, then $\hat{\boldsymbol{\theta}}$ satisfies both consistency and asymptotic normality.*

Computation of the maximum likelihood estimate requires initial values of $\boldsymbol{\theta}$ and starting values $\left\{h_{0},\cdots,h_{-(m-1)}\right\}$ and $\left\{r_{0},\cdots,r_{-(m-1)}\right\}$. Recall from Theorem 1 that

[TABLE]

So $\mu$ can be initialized by replacing $E\left(h_{t}\right)$ in (37) by its sample mean $\overline{h_{t}}$ and a rough guessing of $1-\sqrt{\pi/2}\sum_{i=1}^{p}\alpha_{i}-k\sum_{i=1}^{q}\beta_{i}-\sum_{i=1}^{w}\gamma_{i}$, for example, 0.4. Additionally, initial values of $\alpha_{i}$’s, $\beta_{i}$’s, and $\gamma_{i}$’s can be obtained by setting each of $\sqrt{\pi/2}\sum_{i=1}^{p}\alpha_{i}$, $k\sum_{i=1}^{q}\beta_{i}$, and $\sum_{i=1}^{w}\gamma_{i}$ to be a small value, e.g., 0.2, and by letting $\alpha_{1}=\cdots=\alpha_{p}$, $\beta_{1}=\cdots=\beta_{q}$, $\gamma_{1}=\cdots=\gamma_{w}$. Finally, we assume that the process $\left\{r_{t}\right\}$ starts from its infinite past with a finite mean and variance, so it is reasonable to let $h_{t}=E\left(h_{t}\right)$, $t=0,\cdots,-(m-1)$. An alternative is to let $h_{t}=0$, $t=0,\cdots,-(m-1)$, assuming $\left\{r_{t}\right\}$ starts from constant intervals $\left\{r_{0},\cdots,r_{-(m-1)}\right\}$. In either case, we let $r_{t}=E\left(r_{t}\right)$, $t=0,\cdots,-(m-1)$. Based on our experience, the initial values lead to very fast convergence, with only a couple of iterations.

Alternative to the two-stage estimation, one can also compute the maximum likelihood estimates of $k$ and $\boldsymbol{\theta}$ simultaneously, which jointly maximize the conditional likelihood

[TABLE]

or equivalently the conditional log-likelihood up to a constant

[TABLE]

For the iterative algorithm, the two-stage estimates can be naturally input as initial values. Due to the involvement of the gamma function in the conditional likelihood, computation of this joint maximum likelihood estimates will be more complicated. Additionally, the required conditions to establish the asymptotic normality are also expected to be much more restrictive.

5.2 Simulation and finite sample performances

Our simulation study considers four distinct Int-GARCH models with specific sets of parameter values. From the empirical analysis detailed in Section 6, there seem to be certain ranges for the parameters: 1) $k$ is between $1$ and $2$; 2) $\mu$ and $\alpha_{1}$ are very small; 3) $\beta_{1}$ is generally above $0.3$; 4) $\gamma_{1}$ is either very small or equal to [math]. We design our models such that the parameter values are inside these ranges, with a little added flexibility to account for potential possibilities not included in our real data. We are particularly interested in the Int-ARCH model, since very often in the real data analysis $\gamma_{1}$ is estimated to be [math]. So, after examining two full Int-GARCH models (Models I and II), we also consider two Int-ARCH models with $\gamma_{1}$ set to [math] (Models III and IV). The parameter values for the two Int-GARCH models are generated as follows.

[TABLE]

For Int-ARCH models, $\beta_{1}$ is generally larger due to the removal of $\gamma_{1}$, so we generate the parameters for the two Int-ARCH models with $k,\mu,\alpha_{1}$ being the same as the Int-GARCH models but $\beta_{1}$ slightly larger.

[TABLE]

All of the generations are subject to the constraint that each combination will result in a weakly stationary Int-(G)ARCH process that achieves consistency and asymptotic normality for the maximum likelihood estimators. The exact parameter values are listed in Table 1. Realizations of the designed models are simulated using the initial values $h_{0}=0$ and $r_{0}=E\left(r_{t}\right)$. Plots of two simulated data sets are shown in Figure 2. Denote

[TABLE]

Recall that the theoretical ACF of $\left\{r_{t}\right\}$ is

[TABLE]

We consequently define the sample ACF of $\left\{r_{t}\right\}$ as

[TABLE]

where $\hat{\gamma_{\lambda}}(s)$ and $\hat{\gamma_{\delta}}(s)$ are the sample auto-covariance functions of $\left\{\lambda_{t}\right\}$ and $\left\{\delta_{t}\right\}$, respectively. Figures 3 and 4 show the sample ACF’s for simulated data sets with 1000 observations from an Int-GARCH(1,1,1) (Model I) and an Int-ARCH(1,1) (Model III), respectively.

For each of the four models, we simulate a data set with 1000 observations and estimate the parameters using the proposed two-stage procedure. The process of simulation and estimation is repeated for 100 times independently, and the average results are reported in Table 1. We see that, with a very moderate sample of size $1000$, both the moment estimate for $k$ and the maximum likelihood estimate for $\boldsymbol{\theta}$ are very close to the true values with small mean absolute errors. In addition, the asymptotic standard errors for $\boldsymbol{\theta}$ based on the asymptotic normality are compared to the empirical ones, and they match very well.

6 Empirical Application

In this section, we apply our Int-GARCH model to analyze 10 typical stocks and 4 indices data. The symbols of these stocks and indices are listed, for example, in Table 2. We obtained from Thomson Reuters both the 5-minute intraday data and the daily-closing data of totally 1511 trading days, from January 3, 2006 to December 30, 2011. As an example for demonstration purpose, the sample ACF of the BA (Bank of America) return range is displayed in Figure 5. It is very similar to the theoretical ACF of model I in Figure 1 from Section 5.2, which indicates the feasibility of our Int-GARCH model.

6.1 Data cleaning and ranking metrics

In order to compare the empirical (in-sample and out-of-sample) performances of our Int-GARCH model and the GARCH model, we use realized variance (RV) as the proxy for the “true” market variance. The RV using our 5-minute intraday data is computed as:

[TABLE]

As before, $y_{t}(s)$ denotes the log price of an asset at time $s$ on day $t$, and $y_{t}(s-1)$ denotes the log price of the previous “snapshot” which was taken 5 minutes ago.

Following [31], we applied 4 filtration rules to clean the data.

1. When multiple quotes have the same timestamp, we replace all these with a single entry with the median bid and median ask price.
2. Delete entries for which the spread is negative.
3. Delete entries for which the spread is more than 50 times the median spread on that day.
4. Delete entries for which the mid-quote deviated by more than 10 mean absolute deviations from a rolling median centred but excluding the observation under consideration of 50 observations (i.e. 25 observations before and 25 after).

Denote the RV by $v_{t}$ and estimated volatility by $\hat{\sigma}^{2}_{t}$. It is shown in [24] that the ranking of volatility forecasts based on the $R^{2}$ from the Mincer-Zarnowitz regression

[TABLE]

is robust to noises in $v_{t}$. It is shown in [28] that the ranking using the negative quasi-likelihood (QLIKE) as the loss function is also robust to noises in the volatility proxy. Additionally, the heteroskedasticity-adjusted MSE (HMSE), although not robust to noises, is a popularly used loss function to compare volatility models. The definitions for the QLIKE and HMSE loss functions are given as:

[TABLE]

where $N$ is the length of sampling period.

6.2 In-sample and out-of-sample results

First, we compare the in-sample volatility estimation with that of GARCH (1,1). For this purpose, we fit an Int-GARCH model to the entire data from January 3, 2006 to December 30, 2011. The parameters are estimated by the two-stage estimation method we proposed in Section 5 and results are listed in Table 2. Occasionally, the maximum likelihood estimate of the full model contains negative values. This means that the maximum of the likelihood function subject to the nonnegative constraints happens on the boundary, i.e., at least one of the estimates is [math]. This is equivalent to removing the corresponding parameter(s) and running the unconstrained estimation for the sub-model (see, e.g. [34]). For such cases, the removed parameters are displayed as [math] exactly without a standard error. In general, the fitted models show that $\alpha$ is much smaller than $\beta$ both in magnitude and in statistical significance, indicating that, consistent to our expectation, the return range is of much more importance than a single “snapshot” return, in terms of their contributions to the volatility. Additionally, $\gamma$ also tends to be small, so most likely an Int-ARCH model is sufficient. The comparisons of Int-GARCH (1,1,1) and GARCH (1,1) using all the three criteria ($R^{2}$, QLIKE and HMSE) are shown in Table 3, which generally favor Int-GARCH, except for a few exceptions (i.e., $R^{2}$ of BA, HMSE of T).

Next, to examine the out-of-sample performance, we use the first 5 years (2006-2010) as training period and leave the last year (2011) for out-of-sample predictions. We consider 3 forecasting horizons: 1-step, 2-step, and 5-step, corresponding to one day, two day, and one week ahead predictions, respectively. The model parameters are updated daily and the h-step-ahead predictions (h=1,2,5) are calculated according to equations (20)-(21). The predicts are compared to those from GARCH (1,1) model for the entire year of 2011, and the average results are shown in Tables 4, 5, and 6, for h=1,2,5, respectively. Similar to the in-sample comparison, the out-of-sample comparisons in Tables 4, 5, and 6 indicate the Int-GARCH(1,1,1) perform better than GARCH (1,1) for 1-step, 2-step, and 5-step ahead predictions, for most of the stocks (indices) and ranking metrics.

Finally, we also notice that occasionally the GARCH performs better than the Int-GARCH model in our experiments. For example, in Table 3, $\text{R}^{2}$ of GARCH is higher for BAC, and HMSE of GARCH is smaller for MSFT. In table 4, the Int-GARCH is dominated by GARCH for KO, MSFT, and T. Theoretically, the comparison between the GARCH and the Int-GARCH are complicated. Intuitively understanding, the advantage of our Int-GARCH model is that it makes full use of all possible values in the return range as volatility proxies and summarizes the information by integrating the GARCH mechanism with the random set theory. By comparison, the traditional GARCH model only considers use of the information of a single (closing-to-closing) return as the volatility proxy. Generally, the “summary” that accounts for more information from the data is expected to reflect the volatility more accurately. However, noting that different measures of model performance characterize different features of the model fitting the data, it is therefore possible occasionally that a simpler GARCH model with the point-valued return may outperform the more involved “summary” of Int-GARCH in terms of one performance measure. With regard to this, more extensive experiments with much more and larger data sets are needed to thoroughly investigate and optimize the performance of the Int-GARCH model in the near future.

7 Conclusion

In this paper, we have developed an interval-valued GARCH model for analyzing range-measured return processes. It can be viewed as an extension of the point-valued GARCH model that allows for interval-valued returns to produce “information-richer” estimation of the volatility. Inferences of our Int-GARCH model can be made based on the maximum likelihood method, which has been shown to have both consistency and asymptotic normality. Our empirical study of stocks and indices data has demonstrated the advantages of Int-GARCH model over the GARCH for both in-sample and out-of-sample performances.

An equally important contribution of our Int-GARCH model lies in that it forms a pioneer study in the area of conditional heteroscedasticity for interval-valued time series. Further improvements to the model can be made by incorporating interactions between the level and the range. Körner and Näther ([19]) proposed for a more general space a generalized $L_{2}$ metric, which when restricted to $\mathcal{K}_{\mathcal{C}}(\mathbb{R})$ is

[TABLE]

where $K$ is a symmetric positive definite kernel. It can be represented by the center-radius form as

[TABLE]

where $A$ as a function of $K$ is a symmetric positive definite matrix. So $\rho_{K}$ is a further generalization of $\rho_{W}$ in (19) that takes into account the correlation between the center and the radius, and can be utilized to model the level-range interaction.

8 Appendix

8.1 Proof of Theorem 1

Proof.

By the definition of (24), (17) can be rewritten as

[TABLE]

Expanding $h_{t}$ recursively, we obtain

[TABLE]

Notice that $x_{i,t}$ and $x_{j,s}$ are independent, $\forall i,j\in\mathbb{N}$ and $t\neq s$. Taking expectations on both sizes of (40), we get

[TABLE]

The first term

[TABLE]

The second term

[TABLE]

Therefore, $E\left(h_{t}\right)<\infty$ if and only if $\sum_{i=1}^{k}\mu_{i}<1$. When it is satisfied,

[TABLE]

by the finiteness of $E\left(h_{-\infty}\right)$. The formula for $E\left(r_{t}\right)$ follows immediately from the Aumann expectation. ∎

8.2 Proof of Theorem 2

Proof.

By recursive calculations,

[TABLE]

Taking expectations on both sides, we get

[TABLE]

for all $N\in\mathbb{N}$. Letting $N\to\infty$,

[TABLE]

since $|Ex_{t}|<1$. On the other hand, if $Ex_{t}\geq 1$, then $Eh_{t}>\mu\sum_{i=0}^{\infty}(Ex_{t})^{i}=\infty$. Therefore, $Eh_{t}<\infty$ if and only if

[TABLE]

and when this is satisfied, the Aumann expectation of $r_{t}$ is found to be

[TABLE]

∎

8.3 Proof of Theorem 3

Proof.

Recall, from the proof of Theorem 2, that

[TABLE]

Consequently,

[TABLE]

Note that

[TABLE]

and

[TABLE]

Therefore,

[TABLE]

Let $C_{1}=E\left(x_{t}\right)$ and $C_{2}=E\left(x_{t}^{2}\right)$. By the independence of $h_{t-(N+1)}$ and $x_{t-j}$, $\forall j\leq N+1$ and the fact that $\left\{x_{t}\right\}$ are iid, we obtain,

[TABLE]

The geometric series $\sum_{i=1}^{N}C_{2}^{i}$ and $\sum_{i=1}^{N}C_{1}^{i}$ converge if and only if $\left|C_{1}\right|<1$ and $\left|C_{2}\right|<1$. Since

[TABLE]

only $|C_{2}|<1$ is sufficient here. Under this assumption, the fourth term in the above expression,

[TABLE]

And similarly,

[TABLE]

Therefore, $E\left(h_{t}^{2}\right)<\infty$ if and only if $|C_{2}|<1$, or $C_{2}<1$, since it is positive. Under this assumption and in addition that $Eh^{2}_{-\infty}<\infty$, letting $N\to\infty$, we find the second moment of $h_{t}$ to be

[TABLE]

Consequently, the unconditional variance of $r_{t}$ is

[TABLE]

This completes the proof. ∎

8.4 Proof of Theorem 4

Proof.

First we notice, $\forall t,s\in\mathbb{N}$,

[TABLE]

and therefore,

[TABLE]

The first term

[TABLE]

since $\left\{\varepsilon_{t}\right\}$ are i.i.d.

Similarly, the second term becomes

[TABLE]

Plugging (42) and (43) into (41), we obtain

[TABLE]

where $E\left(h_{t}\right)=\dfrac{\mu}{1-C_{1}}$, $E\left(h_{t}^{2}\right)=\mu^{2}\dfrac{C_{1}+1}{\left(C_{2}-1\right)\left(C_{1}-1\right)}$, and

[TABLE]

(see Lemma 1, Section 8.8). ∎

8.5 Proof of Theorem 5

Proof.

Let $l_{t}(\boldsymbol{\theta})=-(k+1)\log(h_{t})-\frac{\lambda_{t}^{2}}{2h_{t}^{2}}-\frac{\delta_{t}}{h_{t}}$, that is, $l(\boldsymbol{\theta})=\sum\limits_{t=1}^{T}l_{t}(\boldsymbol{\theta})$. The elements in $\nabla l_{t}$ are calculated to be

[TABLE]

In addition, the elements in $\nabla^{2}l_{t}$ are

[TABLE]

(i) To prove consistency, we follow Weiss (1986) to verify the three conditions in Basawa et al. (1976) that guarantees the existence of a consistent root of the equation $\partial L(\boldsymbol{\theta})/\partial\boldsymbol{\theta}=0$.

(1) $\frac{1}{T}\sum_{t}\nabla l_{t}\left(\boldsymbol{\theta}_{0}\right)\stackrel{{\scriptstyle\mathcal{P}}}{{\to}}0$ as $n\to\infty$;

(2) There exists a nonrandom positive definite matrix $M\left(\boldsymbol{\theta}_{0}\right)0$ such that $\forall\epsilon>0$:

[TABLE]

(3) There exists a constant $M<\infty$ such that

[TABLE]

For (1), notice that evaluated at $\boldsymbol{\theta}_{0}$,

[TABLE]

Thus,

[TABLE]

Namely, $E\left[\nabla l_{t}(\boldsymbol{\theta}_{0})\right]=0$. By the ergodic theorem, $\frac{1}{T}\sum_{t=1}^{T}\nabla l_{t}(\boldsymbol{\theta}_{0})\stackrel{{\scriptstyle\mathcal{P}}}{{\to}}0$, $T\to\infty$.

For (2), in view of the ergodic theorem, it is sufficient to verify $E\left[-\nabla^{2}l_{t}(\boldsymbol{\theta}_{0})\right]<\infty$. Notice

[TABLE]

Therefore,

[TABLE]

Since $h_{t}^{-1}\frac{\partial h_{t}}{\partial\boldsymbol{\theta}_{i}}(\boldsymbol{\theta}_{0})$ is bounded by $\max\left\{\frac{1}{\alpha_{i}},\frac{1}{\beta_{j}},\frac{1}{\gamma_{k}},\frac{1}{\mu}\right\}$, we have that $E\left[-\nabla^{2}l_{t}(\boldsymbol{\theta}_{0})\right]<\infty$, and hence condition (2) is satisfied.

For (3), we note that

[TABLE]

By the compactness of $\Theta$, $h_{t}^{-1}\frac{\partial h_{t}}{\partial\boldsymbol{\theta}_{i}}$ is uniformly bounded, and so is $h_{t}^{-3}\frac{\partial h_{t}}{\partial\boldsymbol{\theta}_{i}}\frac{\partial h_{t}}{\partial\boldsymbol{\theta}_{j}}\frac{\partial h_{t}}{\partial\boldsymbol{\theta}_{k}}$. Similarly, $h_{t}^{-1},h_{t}^{-2}$ are also uniformly bounded since $h_{t}>\mu$. So it suffices to have $E\left(\lambda_{t}^{2}\right)<\infty$, which is equivalent to $E\left(h_{t}^{2}\right)<\infty$ and guaranteed by the condition of weak stationarity.

(ii) The requirements in Basawa et al. (1976) for the asymptotic normality of the maximum likelihood estimate are:

(1) $\frac{1}{\sqrt{T}}\sum_{t}\nabla l_{t}(\boldsymbol{\theta}_{0})\stackrel{{\scriptstyle\mathscr{D}}}{{\to}}N(\boldsymbol{0},B)$, $T\to\infty$, for a nonrandom positive definite matrix $B$;

(2) -$\frac{1}{T}\sum_{t}\nabla^{2}l_{t}(\boldsymbol{\theta}_{0})\stackrel{{\scriptstyle\mathcal{P}}}{{\to}}A$, $T\to\infty$, for a nonrandom positive definite matrix $A$;

(3) Condition (3) for consistency.

For (1), since $E\left(\nabla l_{t}(\boldsymbol{\theta}_{0})\right|\mathcal{F}_{t-1})=0$, by the martingale central limit theorem,

[TABLE]

where $\boldsymbol{I}(\boldsymbol{\theta}_{0})=E\left[\nabla\nabla^{T}l_{t}(\boldsymbol{\theta}_{0})\right]$ is the Fisher information matrix. Condition (2) is already shown by the ergodic theorem in the preceding argument. Hence, the proof is completed.

∎

8.6 Proof of Corollary 1

Proof.

It is immediate from Theorem 2, 3, and 4. ∎

8.7 Proof of Corollary 2

Proof.

The Auto-correlation Function (ACF) of $\left\{r_{t}\right\}$ is defined to be $\rho(s)=\mbox{Corr}\left(r_{t},r_{t+s}\right)=\dfrac{\gamma(s)}{\gamma(0)}$. Then the ACF of $\left\{r_{t}\right\}$ is

[TABLE]

∎

8.8 Lemmas

Lemma 1**.**

[TABLE]

Proof.

(i)

[TABLE]

(ii) First, we expand $h_{t+s}$ recursively:

[TABLE]

Consequently,

[TABLE]

Then the expected value is found to be

[TABLE]

where $C_{1}=E\left(x_{t}\right)$.

Finally, remembering that $E\left(h_{t}^{2}\right)=\mu^{2}\dfrac{C_{1}+1}{\left(C_{2}-1\right)\left(C_{1}-1\right)}$, the above calculation is simplified to

[TABLE]

∎

Lemma 2**.**

[TABLE]

Proof.

[TABLE]

∎

Lemma 3**.**

Consider the Int-GARCH(1,1,1) model and assume the condition for mean stationarity in Theorem 2. Then, the process $\left\{h_{t}\right\}$ is strictly stationary and ergodic. Assume in addition that $E\left(x_{t}^{p}\right)<1$, $p>0$, then $E\left(h_{t}^{p}\right)<\infty$.

Proof.

The process $\left\{h_{t}\right\}$ is defined by the stochastic recurrence equation (28), where $x_{t}=\alpha_{1}|\epsilon_{t-1}|+\beta\eta_{t-1}+\gamma_{1}$, $t\in\mathbb{N}$ is an i.i.d. process whose moments are all finite. Therefore, $\left\{h_{t}\right\}$ is strictly stationary and ergodic if

[TABLE]

See, e.g., Nelson (1990). Since $E\left[\log\left(x_{t}\right)\right]\leq\log\left[E\left(x_{t}\right)\right]$, (44) is automatically true under the assumption of mean stationarity, i.e., $E\left(x_{t}\right)<1$. Furthermore, $E\left(h_{t}^{p}\right)<\infty$ iff $E\left(x_{t}^{p}\right)<1$, $p>0$ (Nelson 1990, Theorem 3). This finishes the proof. ∎