A Modified Multiplicative Thinning-Based INARCH Model: Properties, Saddlepoint Maximum Likelihood Estimation, and Application

Abstract

In this article, we propose a modified multiplicative thinning-based integer-valued autoregressive conditional heteroscedasticity model and use the saddlepoint maximum likelihood estimation (SPMLE) method to estimate parameters. A simulation study is given to show a better performance of the SPMLE. The application of the real data, which is concerned with the number of tick changes by the minute of the euro to the British pound exchange rate, shows the superiority of our modified model and the SPMLE.

1. Introduction

In practice, we can often observe a series of integer-valued data that have their own distinguishing characteristics, and many models were proposed for modeling integer-valued time series, such as the integer-valued autoregressive (INAR) process introduced by McKenzie (1985) [1], and Al-Osh and Alzaid (1987) [2]; the integer-valued moving average process proposed by Al-Osh and Alzaid (1988) [3]; the integer-valued autoregressive moving-average model defined by McKenize (1988) [4]; and the integer-valued generalized autoregressive conditional heteroscedasticity (INGARCH) model proposed by Ferland et al. (2006) [5], among others. Here we focus on two kinds of the models above: one is the INAR process, which was introduced as a convenient way to transfer the usual autoregressive structure to a discrete-valued time series, and a p-order model, which is defined as follows:

$$X_t=\sum _{i=1}^p{\alpha }_i\circ X_{t-i}+{\epsilon }_t,$$

where ${\alpha }_i\in [0,1)$ for $i=1,\dots ,p,$ and $\left\{{\epsilon }_t\right\}$ is a sequence of independent and identically distributed (i.i.d.) non-negative integer-valued random variables with $E({\epsilon }_t)=\mu$ and $\mathrm{Var}({\epsilon }_t)={\sigma }_{\epsilon }^2$. The binomial thinning operator ∘ is defined by Steutel and Van Harn (1979) [6] as:

$$\alpha \circ X=\sum _{i=1}^X{Y}_i,\mathrm{if}X>0\mathrm{and}0\mathrm{otherwise},$$

where $Y_i$ are i.i.d. Bernoulli random variables, independent of X, with a success probability are defined by $\alpha$. This model has been generalized by Qian and Zhu (2022) [7], and Huang et al. (2023) [8], among others.

The other is the INGARCH model which was proposed by Ferland et al. (2006) [5] to model the observations of integer-valued time series which exist heteroscedasticity; this INGARCH$(p,q)$ model with a Poisson deviate is defined as:

$$X_t|{\mathcal{F}}_{t-1}:P({\lambda }_t),{\lambda }_t={\alpha }_0+\sum _{i=1}^p{\alpha }_i{X}_{t-i}+\sum _{j-1}^q{\beta }_j{\lambda }_{t-j},$$

where ${\alpha }_0>0,{\alpha }_i\ge 0,{\beta }_j\ge 0,i=1,\dots ,p,j=1,\dots ,q,p\ge 1,q\ge 0$, and ${\mathcal{F}}_{t-1}$ is the $\sigma$-field generated by $\{X_{t-1},X_{t-2},\dots \}$. This model has been generalized by Hu (2016) [9], Liu et al. (2022) [10], and Weiß et al. (2022) [11], among others. Weiß (2018) [12] and Davis et al. (2021) [13] gave recent reviews. According to definitions of INAR and INGARCH models, we noticed that the INAR model is thinning-based, while the INGARCH model is specified by a conditional distribution with a time-varying mean depending on past observations. Combining the thinning-based stochastic equations and the INGARCH model, Aknouche and Scotto (2022) [14] proposed a multiplicative thinning-based INGARCH (MthINGARCH) model to model the integer-valued time series with high overdispersion and persistence. Furthermore, it fits well with heavy-tailed data regardless of the choice of innovation distribution and does not require recourse to complex random coefficient equations. The MthINGARCH model is denoted by:

$$\left\{\begin{array}{c}{X}_t={\lambda }_t{\epsilon }_t,\hfill \\ {\lambda }_t=1+\omega \circ m+\sum _{i=1}^q{\alpha }_i\circ X_{t-i}+\sum _{j=1}^p{\beta }_j\circ {\lambda }_{t-j},\hfill \end{array}\right.$$

(1)

where the symbol ∘ stands for the binomial thinning operator, and $0\le \omega \le 1$, $0\le {\alpha }_i<1$ and $0\le {\beta }_j<1(i=1,\dots ,q,j=1,\dots ,p)$, m is a fixed positive integer number that was introduced for more flexibility. Since there is no explicit probability mass function for the series $\left\{X_t\right\}$, then the traditional maximum likelihood estimation (MLE) cannot be applied to estimate the parameters; therefore, Aknouche and Scotto (2022) [14] used a two-stage weighted least squares estimation instead.

Note that the probability mass function of the random variables cannot be given directly for the likelihood function in some cases; to solve this problem, saddlepoint approximation has been proposed. Daniel (1954) [15] introduced saddlepoint techniques into the statistical field, which have been extended by Field and Ronchetti (1990) [16], Jensen (1995) [17], and Butler (2007) [18]. Saddlepoint techniques have been used successfully in many applications because of the high accuracy with which they can approximate intractable densities and tail probabilities. Pedeli et al. (2015) [19] proposed an alternative approach based on the saddlepoint approximation to log-likelihood, and the saddlepoint maximum likelihood estimation (SPMLE) was used to estimate the parameters of the INAR model, which demonstrates the usefulness of this technique. Thus, through combining the MthINGARCH model of Aknouche and Scotto (2022) [14] and the saddlepoint approximation, we propose a modified multiplicative thinning-based INARCH model for modeling high overdispersion, before applying the saddlepoint method to the estimated parameters. Although the two-stage weighted least squares estimation could be used to estimate the parameters of our modified model, we still adopted the SPMLE as it was still expected to have a better performance than the two-stage weighted least squares estimation in practice. Here, we just consider the INARCH model instead of the INGARCH model because it is difficult and complex to give the conditional cumulant-generating function of random variables for the latter model when applying the saddlepoint approximation.

This article has the following structure. A modified multiplicative thinning-based INARCH model is given, alongside some related properties in Section 2. Moreover, we use the Poisson distribution and geometric distribution for innovations. Section 3 discusses the SPMLE and its asymptotic properties, then simulation studies for both models with SPMLE are also given. A real data example is analyzed with our modified models in Section 4, and comparisons with existing models are made. In-sample and out-of-sample forecasts are used to show the superiority of the SPMLE and our modified model. The conclusion is given in Section 5. Some details of SPMLE and proof of some theorems are presented in the Appendix A.

2. A Multiplicative Thinning-Based INARCH Model

Note that $\mathbb{N}=\{0,1,2,\dots \}$ and $\mathbb{Z}=\{\dots ,-1,0,1,\dots \}$ are the set of non-negative integers and integers, respectively. It can be supposed that $\{{\epsilon }_t,t\in \mathbb{Z}\}$ is a sequence of i.i.d. random variables with a mean of one and finite variance of ${\sigma }^2$. The modified multiplicative thinning-based INARCH (denoted by the MthINARCH$(q)$) model, which we deal with in this paper, is defined by

$$X_t={\lambda }_t{\epsilon }_t,{\lambda }_t=\omega \circ m+\sum _{i=1}^q{\alpha }_i\circ X_{t-i},$$

(2)

where $0<\omega \le 1$, $0\le {\alpha }_i<1,i=1,\dots ,q$, m is a fixed positive integer number. In real applications, we can set m as the upper integer part of the sample mean. It is assumed that the Bernoulli terms corresponding to the binomial variables $\omega \circ m$ and ${\alpha }_i\circ X_{t-i}$ are mutually independent and independent of the sequence $\{{\epsilon }_t,t\in \mathbb{Z}\}$. The reason that we defined the new model in this way can be explained as follows. The additive term 1 in ${\lambda }_t$ and in (1) is unnatural, and is posed to ensure ${\lambda }_t>0$, but we can achieve this by adjusting the range of $\omega$; therefore, we adopted a simple version of ${\lambda }_t$ in (2).

Now that we discuss the conditional mean and conditional variance of $X_t$. Note that ${\mathcal{F}}_{t-1}$ is the $\sigma$-field generated by $X_{t-1},X_{t-2},\dots$. For $E({\epsilon }_t)=1$, let ${\mu }_t:=E(X_t|{\mathcal{F}}_{t-1})=E({\lambda }_t{\epsilon }_t|{\mathcal{F}}_{t-1})=E({\epsilon }_t)E({\lambda }_t|{\mathcal{F}}_{t-1})=E({\lambda }_t|{\mathcal{F}}_{t-1})=\omega m+{\sum }_{i=1}^q{\alpha }_i{X}_{t-i}$. Then we can obtain the conditional variance; first, let ${\nu }_t:=\mathrm{Var}({\lambda }_t|{\mathcal{F}}_{t-1})$ and ${\sigma }_t^2:=\mathrm{Var}(X_t|{\mathcal{F}}_{t-1})$. For $E({\epsilon }_t)=1,\mathrm{Var}({\epsilon }_t)={\sigma }^2$, so $E({\epsilon }_t^2)={\sigma }^2+1$. Therefore,

$$\begin{array}{cc}\hfill {\nu }_t:& =\mathrm{Var}({\lambda }_t|{\mathcal{F}}_{t-1})=\omega (1-\omega )m+\sum _{i=1}^q{\alpha }_i(1-{\alpha }_i)X_{t-i},\hfill \\ \hfill {\sigma }_t^2:& =\mathrm{Var}(X_t|{\mathcal{F}}_{t-1})=E(X_t^2|{\mathcal{F}}_{t-1})-{\left[E(X_t|{\mathcal{F}}_{t-1})\right]}^2=E({\lambda }_t^2|{\mathcal{F}}_{t-1})E({\epsilon }_t^2)-{\mu }_t^2\hfill \\ \hfill & =[\mathrm{Var}({\lambda }_t|{\mathcal{F}}_{t-1})+{(E({\lambda }_t|{\mathcal{F}}_{t-1}))}^2]E({\epsilon }_t^2)-{\mu }_t^2\hfill \\ \hfill & =({\sigma }^2+1)({\nu }_t+{\mu }_t^2)-{\mu }_t^2=({\sigma }^2+1){\nu }_t+{\sigma }^2{\mu }_t^2.\hfill \end{array}$$

Proposition 1.

The necessary and sufficient condition for the first-order stationarity of $X_t$ defined in (2) is that all roots of $1-{\sum }_{i=1}^q{\alpha }_i{z}^i=0$ should lie outside the unit circle.

Proposition 2.

The necessary and sufficient condition for the second-order stationarity of $X_t$ defined in (2) is that $({\sigma }^2+1){\sum }_{i=1}^q{\alpha }_i^2<1.$

Proofs of Propositions 1 and 2 are similar to the proofs of Theorems 2.1 and 2.2 in Aknouche and Scotto (2022) [14], so we omit the details.

For convenience, we need to specify the distribution of $\left\{{\epsilon }_t\right\}$ in (2). First, we let ${\epsilon }_t\sim P(1)$, then $E({\epsilon }_t)=\mathrm{Var}({\epsilon }_t)=1$, and this model is denoted by PMthINARCH$(q)$. It is easy to obtain

$${\mu }_t=\omega m+\sum _{i=1}^q{\alpha }_i{X}_{t-i},{\sigma }_t^2=2{\nu }_t+{\mu }_t^2.$$

Second, let ${\epsilon }_t\sim Ge(p^{*})$. The mean of ${\epsilon }_t$ is $(1-p^{*})/p^{*}=1$, so we have $p^{*}=0.5$ and the variance is ${\displaystyle \mathrm{Var}({\epsilon }_t)=2}$. This model is denoted by GMthINARCH$(q)$, then we have

$${\mu }_t=\omega m+\sum _{i=1}^q{\alpha }_i{X}_{t-i},{\sigma }_t^2=3{\nu }_t+2{\mu }_t^2.$$

3. Parameter Estimation

In this section, we will consider the SPMLE and its asymptotic properties, and a simulation study will be conducted to assess the performance of this estimator.

3.1. Saddlepoint Maximum Likelihood Estimation

Let $\theta ={(\omega ,{\alpha }_1,\dots ,{\alpha }_q)}^{\mathrm{T}}$ be the unknown parameter vector. Note that according to the condition on ${\epsilon }_t$, ${\sigma }^2$ is no longer an unknown parameter. The maximum likelihood estimator of $\theta$ was obtained by maximizing the conditional log-likelihood function

$$l(\theta )=\sum _{t=1}^{n}logP(X_t=x_t|X_{t-1}=x_{t-1},\dots ,X_{t-q}=x_{t-q}),$$

(3)

giving $\widehat{\theta }=arg{max}_{\theta }l(\theta ).$ But the above procedure is challenging to implement because it is difficult to give the likelihood function due to the thinning operations.

Now we discuss the SPMLE. The conditional moment generating function of $X_t$ is

$$\begin{array}{cc}\hfill E({\mathrm{e}}^{u{X}_t}|X_{t-1}=x_{t-1},\dots ,X_{t-q}=x_{t-q})& =E({\mathrm{e}}^{u{\lambda }_t{\epsilon }_t}|X_{t-1}=x_{t-1},\dots ,X_{t-q}=x_{t-q})\hfill \\ \hfill & =E({\mathrm{e}}^{u(\omega \circ m+{\sum }_{i=1}^q{\alpha }_i\circ X_{t-i}){\epsilon }_t}|X_{t-1}=x_{t-1},\dots ,X_{t-q}=x_{t-q})\hfill \\ \hfill & =E({\mathrm{e}}^{u(\omega \circ m){\epsilon }_t})\prod _{i=1}^{q}E({\mathrm{e}}^{u({\alpha }_i\circ x_{t-i}){\epsilon }_t}).\hfill \end{array}$$

Remark 1.

Here we just consider the INARCH model instead of the INGARCH model because for the INGARCH model, the conditional cumulant-generating function of $X_t$ should be given by $E({\mathrm{e}}^{u{X}_t}|X_{t-1}=x_{t-1},\dots ,X_{t-q}=x_{t-q})=E({\mathrm{e}}^{u(\omega \circ m+{\sum }_{i=1}^q{\alpha }_i\circ X_{t-i}+{\sum }_{j=1}^p{\beta }_j\circ {\lambda }_{t-i}){\epsilon }_t}|X_{t-1}=x_{t-1},\dots ,X_{t-q}=x_{t-q}).$ Notice that $X_t$ and ${\lambda }_t$ are correlated, it is difficult and complex to show the conditional cumulant-generating function.

Using the binomial theorem ${(a+b)}^n={\sum }_{k=0}^n{C}_n^k{a}^{n-k}{b}^k$, we have

$$\begin{array}{cc}\hfill E({\mathrm{e}}^{u(\omega \circ m){\epsilon }_t})& =E\left[E({\mathrm{e}}^{u(\omega \circ m){\epsilon }_t}|{\epsilon }_t)\right]=E{(\omega {\mathrm{e}}^{u{\epsilon }_t}+(1-\omega ))}^m\hfill \\ \hfill & =E\left[\sum _{r=0}^m{C}_m^r{(1-\omega )}^r{\omega }^{m-r}{\mathrm{e}}^{u(m-r){\epsilon }_t}\right]=\sum _{r=0}^m{C}_m^r{(1-\omega )}^r{\omega }^{m-r}E({\mathrm{e}}^{u(m-r){\epsilon }_t}).\hfill \end{array}$$

Similarly, we also have

$$E({\mathrm{e}}^{u({\alpha }_i\circ x_{t-i}){\epsilon }_t})=\sum _{r=0}^{x_{t-i}}{C}_{x_{t-i}}^r{(1-{\alpha }_i)}^r{\alpha }_i^{x_{t-i}-r}E({\mathrm{e}}^{u(x_{t-i}-r){\epsilon }_t}).$$

Therefore, for the PMthINARCH$(q)$ model, we have

$$\begin{array}{cc}\hfill E({\mathrm{e}}^{u(\omega \circ m){\epsilon }_t})& =\sum _{r=0}^m{C}_m^r{(1-\omega )}^r{\omega }^{m-r}{\mathrm{e}}^{({\mathrm{e}}^{u(m-r)}-1)},\hfill \\ \hfill E({\mathrm{e}}^{u({\alpha }_i\circ x_{t-i}){\epsilon }_t})& =\sum _{r=0}^{x_{t-i}}{C}_{x_{t-i}}^r{(1-{\alpha }_i)}^r{\alpha }_i^{x_{t-i}-r}{\mathrm{e}}^{({\mathrm{e}}^{u(x_{t-i}-r)}-1)},\hfill \end{array}$$

while for the GMthINARCH$(q)$ model, we have

$$\begin{array}{cc}\hfill E({\mathrm{e}}^{u(\omega \circ m){\epsilon }_t})& =\sum _{r=0}^m{C}_m^r{(1-\omega )}^r{\omega }^{m-r}\frac{1}{2-{\mathrm{e}}^{u(m-r)}},\hfill \\ \hfill E({\mathrm{e}}^{u({\alpha }_i\circ x_{t-i}){\epsilon }_t})& =\sum _{r=0}^{x_{t-i}}{C}_{x_{t-i}}^r{(1-{\alpha }_i)}^r{\alpha }_i^{x_{t-i}-r}\frac{1}{2-{\mathrm{e}}^{u(x_{t-i}-r)}}.\hfill \end{array}$$

Thus the conditional cumulant-generating function of $X_t$ is:

$$\begin{array}{c}\hfill K_t(u)=log\left[E({\mathrm{e}}^{u{X}_t}|X_{t-1}=x_{t-1},\dots ,X_{t-q}=x_{t-q})\right]=logE({\mathrm{e}}^{u(\omega \circ m){\epsilon }_t})+\sum _{i=1}^{q}logE({\mathrm{e}}^{u({\alpha }_i\circ x_{t-i}){\epsilon }_t}).\end{array}$$

A highly accurate approximation to the conditional mass function of $X_t$ at $x_t$ is provided by the saddlepoint approximation:

$${\tilde{f}}_{X_t|X_{t-1}=x_{t-1},\dots ,X_{t-q}=x_{t-q}}(x_t)={\left(2\pi K_t^{″}({\tilde{u}}_t)\right)}^{-\frac{1}{2}}exp\{K_t({\tilde{u}}_t)-{\tilde{u}}_t{x}_t\},$$

(4)

where ${\tilde{u}}_t$ is the unique value of u which satisfies the saddlepoint equation $K_t^{\prime }(u)=x_t,$ with $K_t^{\prime }$ and $K_t^{″}$ represent the first and second order derivatives of $K_t$ with respect to u. Notice that it is difficult to solve the saddlepoint equation $K_t^{\prime }(u)=x_t$ analytically; similar to that mentioned in Pedeli et al. (2015) [19], we can use the Newton–Raphson method to solve this equation.

The log-likelihood function (3) can be approximated by summing the logarithms of the corresponding density approximations (4), yielding:

$${\tilde{L}}_n(\theta )=\sum _{t=1}^n{\tilde{l}}_t(\theta ):=\sum _{t=1}^{n}log{\tilde{f}}_{X_t|X_{t-1}=x_{t-1},\dots ,X_{t-q}=x_{t-q}}(x_t).$$

(5)

The value $\theta$ maximizing this expression is called the saddlepoint maximum likelihood estimator (SPMLE).

3.2. Asymptotic Properties of the SPMLE

Now we discuss the asymptotic properties of the SPMLE. First we give the first-order Taylor expansion of $K_t^{\prime }(u)$ at $u=0$ yields,

$$K_t^{\prime }(u)=K_t^{\prime }(0)+u{K}_t^{″}(0)+o(u)={\mu }_t(\theta )+u{\sigma }_t^2(\theta )+o(u),$$

(6)

where ${\mu }_t(\theta )$ and ${\sigma }_t^2(\theta )$ are the conditional mean and conditional variance of $X_t$. Notice that ${\tilde{u}}_t$ can be given by $K_t^{\prime }({\tilde{u}}_t)=x_t$, so with the Taylor series expansion of $K_t^{\prime }(u)$ in (6), we have:

$${\tilde{u}}_t=\frac{x_t-{\mu }_t(\theta )}{{\sigma }_t^2(\theta )}+o(1),t=q+1,\dots ,n.$$

(7)

Then, we can obtain the second-order Taylor expansion of $K_t(u)$ at $u=0$, which is:

$$K_t(u)\approx u{K}_t^{\prime }(0)+\frac{u^2}{2}{K}_t^{″}(0)=u{\mu }_t(\theta )+\frac{u^2}{2}{\sigma }_t^2(\theta ).$$

(8)

Focusing on the exponent of the saddlepoint approximation (4), Equation (8) gives

$$K_t(u)-u{x}_t\approx u({\mu }_t(\theta )-x_t)+\frac{u^2}{2}{\sigma }_t^2(\theta ).$$

Then using Equation (7), we have

$$K_t({\tilde{u}}_t)-{\tilde{u}}_t{x}_t\approx -\frac{{[x_t-{\mu }_t(\theta )]}^2}{2{\sigma }_t^2(\theta )}.$$

(9)

Hence, we can derive from (8) and (9) that the first-order saddlepoint approximation to the conditional probability mass function is approximately:

$$\begin{array}{cc}\hfill & {\tilde{f}}_{X_t|X_{t-1}=x_{t-1},\dots ,X_{t-q}=x_{t-q}}(x_t)={\left(2\pi K_t^{″}({\tilde{u}}_t)\right)}^{-\frac{1}{2}}\hfill \\ \hfill & \times exp\left[-\frac{{(x_t-\omega m-{\sum }_{i=1}^q{\alpha }_i{x}_{t-i})}^2}{2\left[({\sigma }^2+1)(\omega (1-\omega )m+{\sum }_{i=1}^q{\alpha }_i(1-{\alpha }_i)x_{t-i})+{\sigma }^2{(\omega m+{\sum }_{i=1}^q{\alpha }_i{x}_{t-i})}^2\right]}\right].\hfill \end{array}$$

Therefore, ${\tilde{L}}_n(\theta )={\sum }_{t=1}^n{\tilde{l}}_t(\theta )={\sum }_{t=1}^{n}log{\tilde{f}}_{X_t|X_{t-1}=x_{t-1},\dots ,X_{t-q}=x_{t-q}}(x_t)$ is the quasi-likelihood function for the estimation of $\theta$. To establish the large-sample properties, we have

$$L_n(\theta )=\sum _{t=1}^n{l}_t(\theta )=\sum _{t=1}^{n}log{f}_{X_t|X_{t-1}=x_{t-1},\dots ,X_{t-q}=x_{t-q}}(x_t),$$

which is the ergodic approximation of ${\tilde{L}}_n(\theta )$. The first and second derivatives of the quasi-likelihood function are given in the Appendix A. The strong convergence and asymptotic normality for the SPMLE ${\widehat{\theta }}_n$ are established in the following theorems.

First of all, the assumptions for Theorems 1 and 2 are listed as follows.

Assumption 1.

The solution of the MthINARCH process is strictly stationary and ergodic.

Assumption 2.

Θ is compact and ${\theta }_0\in \stackrel{°}{\Theta }$, where $\stackrel{°}{\Theta }$ denotes the interior of Θ. For technical reasons, we assumed the lower and upper values of each component of parameters as $0<{\omega }_L\le \omega \le {\omega }_U\le 1$ and $0\le {\alpha }_L\le {\alpha }_i\le {\alpha }_U<1$, $i=1,\dots ,q.$

Theorem 1.

Let ${\widehat{\theta }}_n$ be a sequence of SPMLEs satisfying ${\widehat{\theta }}_n=arg\underset{\theta \in \Theta }{max}{\tilde{L}}_n(\theta )$, then under Assumptions 1 and 2, ${\widehat{\theta }}_n$ converges to ${\theta }_0$ almost as surely, as $n\to \infty .$

Theorem 2.

Under Assumptions 1 and 2, there exists a sequence of maximizers ${\widehat{\theta }}_n$ of ${\tilde{L}}_n(\theta )$ such as that of $n\to \infty$,

$$\sqrt{n}({\widehat{\theta }}_n-{\theta }_0)\stackrel{d}{⟶}N(0,{\Sigma }^{-1}),$$

where

$$\Sigma =-E_{{\theta }_0}\left(\frac{{\partial }^2{l}_t({\theta }_0)}{\partial \theta \partial {\theta }^{\mathrm{T}}}\right),$$

and Σ is positively definite.

3.3. Simulation Study

In this section, simulation studies of PMthINARCH$(q)$ and GMthINARCH$(q)$ models for finite sample size are given, where $q=2.$ Here, we used several combinations to show the performance of SPMLE, and the mean absolute deviation error (MADE) ${\displaystyle \frac{1}{s}\sum _{j=1}^s|\widehat{{\theta }_j}-{\theta }_j|}$ was used as the evaluation criterion; here, s is the number of replications. The sample size is $n=100,200,500$, and the number of replications is $s=200$. We used the following combinations of ${(\omega ,{\alpha }_1,{\alpha }_2)}^{\mathrm{T}}$ as the true values to generate the random sample: A1 $={(0.65,0.4,0.4)}^{\mathrm{T}}$, A2 $={(0.9,0.5,0.3)}^{\mathrm{T}}$ for the PMthINARCH$(2)$ model, and B1 $={(0.8,0.4,0.4)}^{\mathrm{T}}$, B2 $={(0.65,0.3,0.5)}^{\mathrm{T}}$ for the GMthINARCH$(2)$ model.

Table 1. Mean and MADE of estimates for PMthINARCH$(2)$ model with SPMLE.

Model				$\mathit{\omega }$	${\mathit{\alpha }}_1$	${\mathit{\alpha }}_2$
A1	m = 3	n = 100	Mean	0.6069	0.5356	0.3569
			MADE	0.3681	0.2866	0.2510
		n = 200	Mean	0.5722	0.5026	0.3952
			MADE	0.3557	0.2434	0.2243
		n = 500	Mean	0.6436	0.4888	0.4140
			MADE	0.2724	0.1287	0.1005
A2	m = 8	n = 100	Mean	0.7782	0.5076	0.4750
			MADE	0.2533	0.2752	0.3007
		n = 200	Mean	0.7935	0.5161	0.4701
			MADE	0.2318	0.2527	0.2778
		n = 500	Mean	0.8703	0.5170	0.4677
			MADE	0.1752	0.2155	0.2390

and

Table 2. Mean and MADE of estimates for GMthINARCH$(2)$ model with SPMLE.

Model				$\mathit{\omega }$	${\mathit{\alpha }}_1$	${\mathit{\alpha }}_2$
B1	m = 4	n = 100	Mean	0.7821	0.2930	0.2870
			MADE	0.1195	0.1499	0.1766
		n = 200	Mean	0.8190	0.3611	0.3185
			MADE	0.1121	0.1425	0.1640
		n = 500	Mean	0.8456	0.3610	0.3298
			MADE	0.0601	0.1331	0.1414
B2	m = 6	n = 100	Mean	0.4718	0.2086	0.3811
			MADE	0.1965	0.1466	0.1463
		n = 200	Mean	0.5186	0.2632	0.5080
			MADE	0.1607	0.1198	0.1412
		n = 500	Mean	0.5468	0.2874	0.4896
			MADE	0.1415	0.1050	0.0770

show the results of these simulations. Notice that as the sample sizes become larger, the MADEs become smaller, and the estimates seem to be close to the true values. Therefore, the SPMLE performs well.

4. A Real Example

Here, we considered the number of tick changes by the minute of the euro to the British pound exchange rate (ExRate for short) on December 12th from 9.00 a.m. to 9.00 p.m. The dataset is available at the website http://www.histdata.com/ (accessed on 17 January 2023). The series comprises of 720 observations with a sample mean of 13.2153 and a sample variance of 224.2498. Obviously, the sample variance is much larger than the sample mean, which shows high overdispersion, and this high overdispersion can also be seen in Figure 1a. Figure 1b,c are the plots of the autocorrelation function (ACF), and the partial autocorrelation function (PACF) means that we know the tick changes are correlated.

We analyzed the data using the PMthINARCH$(3)$ model, GMthINARCH$(3)$ model, Poisson INAR$(3)$ (here denoted by PINAR$(3)$ for short) model, and the INARCH$(3)$ model. The Poisson INAR model is mentioned in Pedeli et al. (2015) [19], and the SPMLE was used to estimate the parameters. Here, the innovations in the PINAR model were assumed to be Poisson with a mean of one. The INARCH model with a Poisson deviate was proposed by Ferland et al. (2006) [5], and the MLE was used to estimate the parameters. According to Aknouche and Scotto (2022) [14], in real applications, we can set m as the upper integer part of the sample mean. Here the sample mean is 13.2153, so m is set to the value of 14.

Table 3. Estimation results: AIC and BIC values for PMthINARCH$(3)$, GMthINARCH$(3)$, PINAR$(3)$ and INARCH$(3)$ models.

PMthINARCH(3)	$\omega$	${\alpha }_1$	${\alpha }_2$	${\alpha }_3$	AIC	BIC
	0.3242	0.5214	0.1945	0.0842	1395.296	1413.613
GMthINARCH(3)	$\omega$	${\alpha }_1$	${\alpha }_2$	${\alpha }_3$	AIC	BIC
	0.4904	0.2532	0.2155	0.2392	1402.472	1420.789
PINAR(3)	${\alpha }_1$	${\alpha }_2$	${\alpha }_3$		AIC	BIC
	0.1335	0.4116	0.3901		1572.806	1586.544
INARCH(3)	$\omega$	${\alpha }_1$	${\alpha }_2$	${\alpha }_3$	AIC	BIC
	8.5670	0.1140	0.1379	0.1009	1524.638	1542.955

gives the estimates of SPMLE and the values of the Akaike information criterion (AIC) and Bayesian information criterion (BIC). According to Table 3, it is clear to see that the values of AIC and BIC of PMthINARCH$(3)$ and GMthINARCH$(3)$ are smaller than those of the PINAR$(3)$ and INARCH$(3)$ models, the values of AIC and BIC of INARCH$(3)$ are smaller than those of the PINAR$(3)$ model. Moreover, the values of AIC and BIC of PMthINARCH$(3)$ are smaller than those of GMthINARCH$(3)$. In summary, the INARCH model performed better than the PINAR model; meanwhile, the PMthINARCH model and GMthINARCH model performed better than the PINAR model and INARCH model.

According to Aknouche and Scotto (2022) [14], the two-stage weighted least squares estimation (2SWLSE) was used to estimate the parameters of the MthINGARCH model. Therefore, to compare the performance of 2SWLSE and SPMLE, and the performance of PMthINARCH, GMthINARCH, and PINAR models, to consider the in-sample and out-of-sample forecasts of these two estimation methods and the three models above, respectively. First, we considered the in-sample forecast. We used all of the observations to estimate the model, and then we could forecast the last 10 observations 711–720, the last 15 observations 706–720, and the last 20 observations 701–720; these three-time horizons of in-sample forecast are denoted by C1, C2, and C3, respectively. Similar to the in-sample forecast process, we also considered the out-of-sample forecast and divided all the observations into three-time horizons: the first one was 1–710 and 711–720, the second one was 1–705 and 706–720, and the third one was 1–700 and 701–720, which are denoted by D1, D2, and D3, respectively.

Here we illustrate the performance of the considered models by comparing the MADEs of each forecast. The MADEs of in-sample forecasts and out-of-sample forecasts for three models with SPMLE are shown in

Table 4. MADEs of in-sample forecasts and out-of-sample forecasts for PMthINARCH$(3)$, GMthINARCH$(3)$, and PINAR$(3)$ models with SPMLE.

Methods of Forecast		PMthINARCH	GMthINARCH	PINAR
In-sample	C1	15.30	16.80	17.40
	C2	15.87	17.67	18.40
	C3	16.65	20.70	21.90
Out-of-sample	D1	17.50	17.70	22.50
	D2	19.47	19.80	23.80
	D3	20.50	25.25	27.50

. The MADEs of the in-sample forecasts and out-of-sample forecasts for the PMthINARCH model with 2SWLSE and SPMLE are shown in

Table 5. MADEs of in-sample forecasts and out-of-sample forecasts for PMthINARCH$(3)$ model with SPMLE and 2SWLSE.

Methods of Forecast		SPMLE	2SWLSE
In-sample	C1	15.30	16.20
	C2	15.87	17.20
	C3	16.65	18.55
Out-of-sample	D1	17.50	18.60
	D2	19.47	21.67
	D3	20.50	22.70

, and the in-sample forecasts and out-of-sample forecasts for the GMthINARCH model with 2SWLSE and SPMLE are shown in

Table 6. MADEs of in-sample forecasts and out-of-sample forecasts for GMthINARCH$(3)$ model with SPMLE and 2SWLSE.

Methods of Forecast		SPMLE	2SWLSE
In-sample	C1	16.80	17.20
	C2	17.67	18.07
	C3	20.70	21.05
Out-of-sample	D1	17.70	19.90
	D2	19.80	22.87
	D3	25.25	26.50

. According to Table 4, the MADEs of PMthINARCH$(3)$ and GMthINARCH$(3)$ are smaller than those of PINAR$(3)$, Table 5 and Table 6 show that the MADEs of PMthINARCH$(3)$ and GMthINARCH$(3)$ of SPMLE are smaller than those of 2SWLSE; meanwhile, in these three Tables, the MADEs of in-sample forecasts were smaller than those of out-of-sample forecasts. In summary, the PMthINARCH model and GMthINARCH model were superior to the PINAR model in modeling this real data set, and the PMthINARCH model performed better than the GMthINARCH model. Meanwhile, the performance of SPMLE was better than 2SWLSE for MthINARCH models.

5. Conclusions

In this paper, we modified a multiplicative thinning-based INARCH model. The probability mass function of random variables is provided by saddlepoint approximation. We used the SPMLE to estimate the parameters and obtain the asymptotic distribution of the SPMLE. Moreover, to show the superiority of the MthINARCH models and the SPMLE, we used the PMthINARCH$(q)$ process and GMthINARCH$(q)$ process for discussion and comparison. The SPMLE performs well in the simulation studies. A real dataset indicates that the PMthINARCH model and the GMthINARCH model are able to describe the overdispersed integer-valued data, and the real data example leads to a superior performance of the MthINARCH models compared with the PINAR and INARCH models. In addition, the results also show a superior performance of SPMLE compared with 2SWLSE.

For further discussion, more research is needed for some aspects. Here we used the Poisson distribution and geometric distribution for ${\epsilon }_t$; however, we could use the negative binomial distribution or some zero-inflated distributions as well. Moreover, we just considered the INARCH model, so the corresponding INGARCH model should be considered as well.