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

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.


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 integervalued 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: where α i ∈ [0, 1) for i = 1, . . . , p, and {ε t } is a sequence of independent and identically distributed (i.i.d.) non-negative integer-valued random variables with E(ε t ) = µ and Var(ε t ) = σ 2 ε . The binomial thinning operator • is defined by Steutel and Van Harn (1979) [6] as: where Y i are i.i.d. Bernoulli random variables, independent of X, with a success probability are defined by α. 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: where α 0 > 0, α i ≥ 0, β j ≥ 0, i = 1, . . . , p, j = 1, . . . , q, p ≥ 1, q ≥ 0, and F t−1 is the σ-field generated by {X t−1 , X t−2 , . . . }. 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: where the symbol • stands for the binomial thinning operator, and 0 ≤ ω ≤ 1, 0 ≤ α i < 1 and 0 ≤ β j < 1 (i = 1, . . . , q, j = 1, . . . , 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 {X t }, 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.

A Multiplicative Thinning-Based INARCH Model
Note that N = {0, 1, 2, . . . } and Z = {. . . , −1, 0, 1, . . . } are the set of non-negative integers and integers, respectively. It can be supposed that {ε t , t ∈ Z} is a sequence of i.i.d. random variables with a mean of one and finite variance of σ 2 . The modified multiplicative thinning-based INARCH (denoted by the MthINARCH(q)) model, which we deal with in this paper, is defined by where 0 < ω ≤ 1, 0 ≤ α i < 1, i = 1, . . . , 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 ω • m and α i • X t−i are mutually independent and independent of the sequence {ε t , t ∈ Z}. The reason that we defined the new model in this way can be explained as follows. The additive term 1 in λ t and in (1) is unnatural, and is posed to ensure λ t > 0, but we can achieve this by adjusting the range of ω; therefore, we adopted a simple version of λ t in (2). Now that we discuss the conditional mean and conditional variance of X t . Note that Then we can obtain the conditional variance; first, let ν t := Var(λ t |F t−1 ) and σ 2 t := Var(X t |F t−1 ). For E(ε t ) = 1, Var(ε t ) = σ 2 , so E(ε 2 t ) = σ 2 + 1. Therefore, Proposition 1. The necessary and sufficient condition for the first-order stationarity of X t defined in (2) is that all roots of 1 − ∑ q i=1 α 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 ( 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 {ε t } in (2). First, we let ε t ∼ P(1), then E(ε t ) = Var(ε t ) = 1, and this model is denoted by PMthINARCH(q). It is easy to obtain Second, let ε t ∼ Ge(p * ). The mean of ε t is (1 − p * )/p * = 1, so we have p * = 0.5 and the variance is Var(ε t ) = 2. This model is denoted by GMthINARCH(q), then we have

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.

Saddlepoint Maximum Likelihood Estimation
Let θ = (ω, α 1 , . . . , α q ) T be the unknown parameter vector. Note that according to the condition on ε t , σ 2 is no longer an unknown parameter. The maximum likelihood estimator of θ was obtained by maximizing the conditional log-likelihood function givingθ = arg max θ l(θ). 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 . Notice that X t and λ t are correlated, it is difficult and complex to show the conditional cumulant-generating function.
Using the binomial theorem (a + b) n = ∑ n k=0 C k n a n−k b k , we have Similarly, we also have Therefore, for the PMthINARCH(q) model, we have while for the GMthINARCH(q) model, we have Thus the conditional cumulant-generating function of X t is: A highly accurate approximation to the conditional mass function of X t at x t is provided by the saddlepoint approximation: whereũ t is the unique value of u which satisfies the saddlepoint equation K t (u) = x t , with K t 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 (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: The value θ maximizing this expression is called the saddlepoint maximum likelihood estimator (SPMLE).

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 (u) at u = 0 yields, where µ t (θ) and σ 2 t (θ) are the conditional mean and conditional variance of X t . Notice that u t can be given by K t (ũ t ) = x t , so with the Taylor series expansion of K t (u) in (6), we have: Then, we can obtain the second-order Taylor expansion of K t (u) at u = 0, which is: Focusing on the exponent of the saddlepoint approximation (4), Equation (8) gives Then using Equation (7), we have Hence, we can derive from (8) and (9) that the first-order saddlepoint approximation to the conditional probability mass function is approximately: is the quasi-likelihood function for the estimation of θ. To establish the large-sample properties, we have which is the ergodic approximation ofL n (θ). The first and second derivatives of the quasilikelihood function are given in the Appendix A. The strong convergence and asymptotic normality for the SPMLEθ n are established in the following theorems.
First of all, the assumptions for Theorems 1 and 2 are listed as follows.

Theorem 2.
Under Assumptions 1 and 2, there exists a sequence of maximizersθ n ofL n (θ) such as that of n → ∞, ∂θ∂θ T , and Σ is positively definite.

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) 1 s s ∑ j=1 |θ j − θ 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 (ω, α 1 , α 2 ) T as the true values to generate the random sam-  Tables 1 and 2 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.

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) [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 gives the estimates of SPMLE and the values of the Akaike information criterion (AIC) and Bayesian information criterion (BIC). According to   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-ofsample 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. 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, and the in-sample forecasts and out-of-sample forecasts for the GMthI-NARCH model with 2SWLSE and SPMLE are shown in Table 6. According to Table 4, the MADEs of PMthINARCH(3) and GMthINARCH(3) are smaller than those of PINAR(3), Tables 5 and 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 insample 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 GMthI-NARCH model. Meanwhile, the performance of SPMLE was better than 2SWLSE for MthINARCH models.

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 ε 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.

Data Availability Statement:
The dataset is available at the website http://www.histdata.com/ (accessed on 17 January 2023).

Acknowledgments:
The authors are very grateful to three reviewers for their constructive suggestions and comments, leading to a substantial improvement in the presentation and contents.

Conflicts of Interest:
The authors declare no conflict of interest.

Appendix A Appendix A.1. Details of SPMLE
Here, we give the derivatives of K t (u) mentioned in Section 3.1 of PMthINARCH(q) and GMthINARCH(q). Now we give K t (u) and K t (u) of PMthINARCH(q). In Section 3.1, we have so the derivatives of K t (u) are given by Then we give K t (u) and K t (u) of GMthINARCH(q). In Section 3.1, we have so the derivatives of K t (u) are given by

Appendix A.2. Derivatives of the Quasi-Likelihood Function
The conditional log-quasi-likelihood function l t (θ) is continuous on Θ: for 1 ≤ t ≤ n, Then the first and second derivatives of µ t (θ) and σ 2 t (θ) can be easily expressed by Next, we check (iii). By Jensen's inequality, we have The equality holds if Then the proof of (iv) is similar to that in the Supplementary Material A.4 in Xu and Zhu (2022) [22]. Here we omit the details.
Using vector form of Suppose the left-hand side is 0 almost surely, then the right-hand side is also 0 almost surely, which can be written as So the coefficients of the above equation must satisfy ξ i (σ 2 + 1) = 0, 2σ 2 ξ i = 0, i = 0, . . . , q.
Appendix A.5. Lemmas for the Proof of Theorem 2 Similar to the proof of Theorem 1.2 in Hu (2016) [9], we give some related lemmas for the proof of Theorem 2. According to the derivatives of the quasi-likelihood function, we have Likewise for the other terms of parameters.
Similarly, we can prove other terms, thus, the first part of (A1) holds. The proof of the second part of (A1) is similar, here we omit the details.
holds by the central limit theorem for martingale difference sequence in Billingsley (1961).
Before showing Lemma A2, we have we use T n to derive the asymptotic distribution ofθ n . For any u ∈ R q+1 , the Taylor series expansion of T n (u) at θ 0 is where θ * = θ * n (u) is on the line segment connecting θ 0 and θ 0 + u √ n . For Euclidean distance · and any compact set K ⊂ R q+1 , sup u∈K θ * − θ 0 → 0, as n → ∞.
Lemma A3. For any compact set K ∈ R q+1 and any ε > 0, Because of Lemmas A1 and A2, we have which completes our proof.