# Tobit Models for Count Time Series

Christian H. Weiß  Helmut Schmidt University, Department of Mathematics and
Statistics, Hamburg, Germany. E-Mail<EMAIL_ADDRESS>ORCID:
0000-0001-8739-6631. Fukang Zhu  School of Mathematics, Jilin University, 2699
Qianjin, Changchun 130012, China. E-Mail<EMAIL_ADDRESS>ORCID:
0000-0002-8808-8179.Corresponding author.

###### Abstract

Several models for count time series have been developed during the last
decades, often inspired by traditional autoregressive moving average (ARMA)
models for real-valued time series, including integer-valued ARMA (INARMA) and
integer-valued generalized autoregressive conditional heteroscedasticity
(INGARCH) models. Both INARMA and INGARCH models exhibit an ARMA-like
autocorrelation function (ACF). To achieve negative ACF values within the
class of INGARCH models, log and softplus link functions are suggested in the
literature, where the softplus approach leads to conditional linearity in good
approximation. However, the softplus approach is limited to the INGARCH family
for unbounded counts, i. e., it can neither be used for bounded counts, nor
for count processes from the INARMA family. In this paper, we present an
alternative solution, named the Tobit approach, for achieving approximate
linearity together with negative ACF values, which is more generally
applicable than the softplus approach. A Skellam–Tobit INGARCH model for
unbounded counts is studied in detail, including stationarity, approximate
computation of moments, maximum likelihood and censored least absolute
deviations estimation for unknown parameters and corresponding simulations.
Extensions of the Tobit approach to other situations are also discussed,
including underlying discrete distributions, INAR models, and bounded counts.
Three real-data examples are considered to illustrate the usefulness of the
new approach.

Key words: count time series; INGARCH models; maximum likelihood estimation;
negative autocorrelation; Skellam distribution; Tobit model.

##  1 Introduction

Count time series $x_{1},\ldots,x_{n}$, $n\in\mathbb{N}=\\{1,2,\ldots\\}$, are
discrete-valued quantitative time series, which stem from an underlying count
process $\\{X_{t}\ |\ t\in\mathbb{Z}\\}$,
$\mathbb{Z}=\\{\ldots,-1,0,1,\ldots\\}$, with the range being contained in
$\mathbb{N}_{0}=\\{0,1,\ldots\\}$. More precisely, if the range of
$\\{X_{t}\\}$ is equal to $\mathbb{N}_{0}$, we are concerned with unbounded
counts, whereas bounded counts have the finite range $\\{0,\ldots,N\\}$ with
some upper bound $N\in\mathbb{N}$ (Weiß, 2018). Several models for count time
series have been developed during the last decades, often inspired by the
traditional autoregressive moving average (ARMA) models for real-valued time
series. Such proposals include the integer-valued ARMA (INARMA) models for
stationary count time series (among others by McKenzie, 1985; Alzaid & Al-Osh,
1990; Du & Li, 1991), where the ARMA recursion is adapted by using types of
“thinning operators”, and the integer-valued generalized autoregressive
conditional heteroscedasticity (INGARCH) models (among others by Ferland et
al., 2006; Fokianos et al., 2009; Xu et al., 2012; Zhu, 2012), which are
conditional regression models having a linear conditional mean. Both the
INARMA and the INGARCH models exhibit an ARMA-like autocorrelation structure,
in the sense that their autocorrelation function (ACF) satisfies a set of
Yule–Walker equations. However, by contrast to the ordinary ARMA models, the
attainable range of ACF values of these ARMA-like count models is often
restricted to only positive ACF values. This follows from the parameter
constraints being needed to ensure the non-negative outcomes for the process.

To achieve negative ACF values within the class of conditional regression
models, a traditional solution is to use a log link (Fokianos & Tjøstheim,
2011), but then the model behaves distinctly non-linear and the ARMA-like ACF
is lost. Linear models, by contrast, are popular among practitioners as they
are “simple, useful, and interpretable in a wide range of contexts” (Grunwald
et al., 2000, p. 479). They also allow to describe important stochastic
properties by closed formulae (e. g., the ACF by Yule–Walker equations) and to
derive simple moment estimators for initial model fit. Thus, Weiß et al.
(2022) recently proposed the family of softplus INGARCH models, which use the
softplus function as a smooth but nearly (piecewise) linear response function.
This model successfully combines an approximately ARMA-like ACF with the
ability to handle negative ACF values. However, the softplus approach is
limited to the INGARCH family for unbounded counts, i. e., it can neither be
used for bounded counts, nor for count processes from the INARMA family.

Therefore, in this paper, we present an alternative solution for achieving
approximate linearity together with negative ACF values, which is more
generally applicable than the softplus approach. For comparison purposes, our
novel approach is first exemplified and discussed within the INGARCH family
for unbounded counts, but later we show how to adapt it to other situations as
well. Thus, as the first step, we propose a novel family of conditional
regression models for unbounded counts, which is defined as a combination of
the INGARCH approach with the so-called Tobit model (_To_ bin’s pro _bit_
model, see Tobin (1958); Goldberger (1964) as well as Section 2). An
INGARCH$(p,q)$ model with $p\geq 1$ and $q\geq 0$ is characterized by
requiring the conditional mean $M_{t}=E(X_{t}|\mathcal{F}_{t-1})$, where
$\mathcal{F}_{t-1}$ is the $\sigma$-field generated by
$\\{X_{t-1},X_{t-2},\ldots\\}$, to be a linear expression in the last $p$
observations and the last $q$ conditional means (“feedback terms”), i. e.,

$M_{t}\ =\ a_{0}\ +\ \sum_{i=1}^{p}\,a_{i}\,X_{t-i}\ +\
\sum_{j=1}^{q}\,b_{j}\,M_{t-j}.$ (1.1)

In standard INGARCH models, see Appendix A.1 for details, $M_{t}$ is used as
the conditional mean of a count random variable and, thus, has to be a
positive real number. Therefore, recursion (1.1) is accompanied by the
constraints $a_{0}>0$ and $a_{1},\ldots,a_{p},b_{1},\ldots,b_{q}\geq 0$,
which, in turn, cause the ACF to be non-negative. So the first idea to allow
for negative ACF values might be to skip the aforementioned parameter
constraints and to use $\max\big{\\{}0,M_{t}\big{\\}}$ as the conditional mean
of the generated count. Since the response function $\max\\{0,x\\}$ is also
known as the rectified linear unit (ReLU) function, let us refer to this
approach as the _ReLU INGARCH model_. This model, however, has a serious
drawback: if $M_{t}\leq 0$ happens at time $t$, then the next count generated
by the ReLU INGARCH model has the mean zero, i. e., we have a degenerate count
distribution with all probability mass in zero. This, in turn, will cause
problems in likelihood computations, time series forecasting, and others.
Therefore, a different solution is required to allow $M_{t}$ according to
(1.1) to become negative.

The above drawbacks of the ReLU INGARCH model are avoided by our novel _Tobit
INGARCH model_. There, we first generate an integer-valued random variable
$X_{t}^{*}$ with range $\mathbb{Z}$ (so including negative integers!) and
having the conditional mean $M_{t}$ (which is possibly negative). $X_{t}^{*}$
might be understood as a kind of “auxiliary” or “intermediate” variable, which
is only required to catch the case where the negative dependence parameters in
(1.1) cause a negative outcome. Such negative outcomes are undesirable as we
aim at modeling a count process. Therefore, in a second step, $X_{t}^{*}$ is
left-censored at zero by computing $X_{t}=\max\big{\\{}0,X_{t}^{*}\big{\\}}$
(Tobit approach); details are given in Section 2. Note that our Tobit approach
differs from the aforementioned ReLU approach, because the counts generated by
$X_{t}=\max\big{\\{}0,X_{t}^{*}\big{\\}}$ always have a truly positive mean,
i. e., we never have a degenerate count distribution (provided that
$X_{t}^{*}$ takes positive values with non-zero probability). A detailed
discussion of this general model as well as of its stochastic properties is
provided by Section 2. Afterwards in Section 3, we consider the important
special case, where $X_{t}^{*}$ is generated from a conditional Skellam
distribution (Johnson et al., 2005, Section 4.12.3), leading to the
Skellam–Tobit INGARCH model. For this model, we also evaluate its abilities
for approximating linearity and for capturing negative ACF values. Section 4
then presents and investigates approaches for parameter estimation. Section 5
discusses a variety of possible extensions, i. e., how the Tobit approach
could be adapted to models outside the unbounded INGARCH class, and the
extensions are mainly about underlying discrete distributions, INAR models,
and bounded counts. Section 6 illustrates the application of the Tobit count
models to three real-world count time series: lottery winners counts, chemical
process yields, and daily air quality levels. Finally, Section 7 concludes and
outlines issues for future research.

##  2 The Class of Tobit INGARCH Models

The basic Tobit model, a censored probit regression model, was introduced by
Tobin (1958) and has been extended in several ways, see Chapter 6 in Maddala
(1983) for a comprehensive survey. Some of these extensions refer to time
series data. For example, a dynamic Tobit model for real-valued random
variables $Y_{t}$ (left-censored at zero) might be of the form

$Y_{t}\ =\ \max\big{\\{}0,\ Y_{t}^{*}\big{\\}}\quad\text{with}\quad Y_{t}^{*}\
=\ \sum_{i=1}^{p}\rho_{i}\,Y_{t-i}+\epsilon_{t}\
+\mbox{\boldmath$\gamma$}^{\top}\mbox{\boldmath$z$}_{t},$ (2.1)

see Maddala (1983); Honoré (1993); Honoré & Hu (2004); de Jong & Herrera
(2011). If $\epsilon_{t}$ is normally distributed, then the conditional
distribution of $Y_{t}^{*}$ given $Y_{t-1},\ldots$ and
$\mbox{\boldmath$z$}_{t}$ is also normal, i. e., $Y_{t}$ with range
$[0;\infty)$ is generated by censoring the $\mathbb{R}$-valued $Y_{t}^{*}$ at
zero. In what follows, this idea is transfered to the censoring of a
$\mathbb{Z}$-valued random variable $X_{t}^{*}$, leading to the count random
variable $X_{t}=\max\big{\\{}0,X_{t}^{*}\big{\\}}$ with range
$\mathbb{N}_{0}$.

To overcome the limitations of the ordinary INGARCH model while
(approximately) preserving its additive structure, we propose the novel class
of Tobit INGARCH models. For defining this model, let us introduce the
$\mathbb{Z}$-valued random operator
$\mathcal{S}_{\mbox{\scriptsize\boldmath$\nu$}}(\mu)$, where
$\mu\in\mathbb{R}$ denotes the mean of
$\mathcal{S}_{\mbox{\scriptsize\boldmath$\nu$}}(\mu)$, and $\nu$ comprises
some further distributional parameters (e. g., such as a dispersion
parameter). The letter “$\mathcal{S}$” might be understood as referring to _s_
igned integer outcomes from $\mathbb{Z}$. Given the values of $\mu$ and $\nu$,
$\mathcal{S}_{\mbox{\scriptsize\boldmath$\nu$}}(\mu)$ is a random variable on
$\mathbb{Z}$. The terminology “operator” is used in analogy to the thinning
operators in INARMA modeling, to express that
$\mathcal{S}_{\mbox{\scriptsize\boldmath$\nu$}}(\mu)$ acts on $\mu$ by
generating an integer value. It is assumed that the distribution of
$\mathcal{S}_{\mbox{\scriptsize\boldmath$\nu$}}(\mu)$ is such that
$P\big{(}\mathcal{S}_{\mbox{\scriptsize\boldmath$\nu$}}(\mu)>0\big{)}>0$ holds
for all $\mu$ and $\nu$. Later in Section 3, we consider the special case of
the operator following the Skellam distribution with mean $\mu$ and dispersion
parameter $\delta>0$, where this assumption is satisfied.

In analogy to (2.1), the _Tobit INGARCH model_ is defined by the equation

$X_{t}\ =\ \max\big{\\{}0,\ X_{t}^{*}\big{\\}}\quad\text{with}\quad X_{t}^{*}\
=\ \mathcal{S}_{\mbox{\scriptsize\boldmath$\nu$}}(M_{t}),$ (2.2)

where we assume the conditional mean $M_{t}=E(X_{t}^{*}|\mathcal{F}_{t-1})$ to
satisfy the linear INGARCH equation

$M_{t}\ =\
\alpha_{0}+\sum_{i=1}^{p}\alpha_{i}\,X_{t-i}+\sum_{j=1}^{q}\beta_{j}\,M_{t-j}$
(2.3)

with
$\alpha_{0},\alpha_{1},\ldots,\alpha_{p},\beta_{1},\ldots,\beta_{q}\in\mathbb{R}$,
recall (1.1). Note that (2.2) is equivalently stated as
$X_{t}=X_{t}^{*}\,\mathds{1}_{\mathbb{N}}(X_{t}^{*})$, because the integer
ReLU function $\max\\{0,x\\}$ can be rewritten as
$x\,\mathds{1}_{\mathbb{N}}(x)$, where the indicator function
$\mathds{1}_{A}(x)$ equals 1 (0) if $x\in A$ ($x\not\in A$).

According to Grunwald et al. (2000), a time series model for $\\{X_{t}\\}$ can
be classified as being linear (thus satisfying, e. g., a set of Yule–Walker
equations for the ACF) if its conditional mean $E(X_{t}|\mathcal{F}_{t-1})$ is
linear in past observations. From (2.2) and (2.3), it gets clear that the
$\mathbb{Z}$-valued $X_{t}^{*}$ has an exactly linear conditional mean and
would, thus, satisfy the linearity requirement of Grunwald et al. (2000).
However, it might happen that
$X_{t}^{*}=\mathcal{S}_{\mbox{\scriptsize\boldmath$\nu$}}(M_{t})$ takes a
negative value, which contradicts the count nature of the target process.
Therefore, (2.2) turns $X_{t}^{*}$ into the count $X_{t}$ that we are actually
interested in, at the price that $X_{t}$ is no longer linear in $X_{t-i}$ and
$M_{t-j}$ because of the censoring. Nevertheless, we shall see in Section 3.3
below that by an appropriate choice of
$\mathcal{S}_{\mbox{\scriptsize\boldmath$\nu$}}(\cdot)$, we can achieve that
$E(X_{t}|\mathcal{F}_{t-1})$ behaves at least _nearly_ linearly in past
observations, such that the corresponding Tobit INGARCH model can be applied
to linear count time series data in practice.

For being able to evaluate the extent of linearity, let us first derive
formulae for the conditional probabilities and moments of $\\{X_{t}\\}$. While
$X_{t}^{*}$ is a $\mathbb{Z}$-valued random variable, the $X_{t}$ generated by
(2.2) is a count random variable with conditional distribution

$P(X_{t}=x|\mathcal{F}_{t-1})\ =\
\left\\{\begin{array}[]{ll}P\big{(}\mathcal{S}_{\mbox{\scriptsize\boldmath$\nu$}}(M_{t})=x\big{)}&\text{if
}x>0,\\\\[2.15277pt]
P\big{(}\mathcal{S}_{\mbox{\scriptsize\boldmath$\nu$}}(M_{t})\leq
0\big{)}&\text{if }x=0.\end{array}\right.$

As a result, the conditional mean of $X_{t}$ given $\mathcal{F}_{t-1}$
computes as

$E(X_{t}|\mathcal{F}_{t-1})\ =\
E\big{(}X_{t}^{*}\,\mathds{1}_{\mathbb{N}}(X_{t}^{*})\big{|}\mathcal{F}_{t-1}\big{)}\
=\
\sum_{x=1}^{\infty}x\,P\big{(}\mathcal{S}_{\mbox{\scriptsize\boldmath$\nu$}}(M_{t})=x\big{)}.$
(2.4)

Higher-order conditional moments are computed analogously as

$E(X_{t}^{r}|\mathcal{F}_{t-1})=E\big{(}(X_{t}^{*})^{r}\,\mathds{1}_{\mathbb{N}}(X_{t}^{*})\big{|}\mathcal{F}_{t-1}\big{)}\quad\text{with}\quad
r\in\mathbb{N}.$ (2.5)

We note that these moments are equal to the corresponding partial moments of
$X_{t}^{*}$ (left-censored at zero). So having specified the
$\mathbb{Z}$-valued distribution of
$\mathcal{S}_{\mbox{\scriptsize\boldmath$\nu$}}(\mu)$, these partial moments
can be derived. This is illustrated in Section 3, where the special case of a
conditional Skellam distribution is discussed in detail.

##  3 Skellam–Tobit INGARCH Model

###  3.1 Skellam Distribution and Partial Moments

Probably the most well-known example of a parametric distribution for a
$\mathbb{Z}$-valued random variable is the _Skellam distribution_ ,
$\textup{Sk}(\lambda_{1},\lambda_{2})$ with $\lambda_{1},\lambda_{2}>0$, see
Section 4.12.3 in Johnson et al. (2005) and the references mentioned therein.
It is sometimes also referred to as the “Poisson-difference distribution”,
because $X^{*}\sim\textup{Sk}(\lambda_{1},\lambda_{2})$ can be generated as
the difference $X^{*}=X_{1}-X_{2}$ of the two independent Poisson variables
$X_{1}\sim\textup{Poi}(\lambda_{1})$ and $X_{2}\sim\textup{Poi}(\lambda_{2})$.
Its mean and variance are given by $\mu=\lambda_{1}-\lambda_{2}\in\mathbb{R}$
and $\sigma^{2}=\lambda_{1}+\lambda_{2}>|\mu|\geq 0$, respectively. Its
probability mass function (PMF) equals

$P(X^{*}=x)\ =\ {\rm
e}^{-\lambda_{1}-\lambda_{2}}\,\big{(}\lambda_{1}/\lambda_{2}\big{)}^{x/2}\,I_{x}\big{(}2\,\sqrt{\lambda_{1}\lambda_{2}}\big{)}\qquad\text{for
}x\in\mathbb{Z},$ (3.1)

where

$I_{n}(z)\>=\
\left(\frac{z}{2}\right)^{n}\sum_{k=0}^{\infty}\frac{(\frac{z}{2})^{2k}}{k!\cdot(k+n)!}\
\text{for }n\in\mathbb{N}_{0},\qquad I_{-n}(z)\>=\ I_{n}(z),$ (3.2)

denotes the modified Bessel function of the first kind, which is implemented
in R by besselI. This distribution has been used to analyze some independent
data (such as Karlis & Ntzoufras (2006, 2009)) and time series data (such as
Koopman et al. (2017)).

###### Remark 1

The modified Bessel function of the first kind satisfies

$I_{n+1}(z)\ =\ I_{n-1}(z)\ -\ \frac{2n}{z}\cdot I_{n}(z).$ (3.3)

In addition, we have

$\displaystyle
I_{n}^{\prime}(z)=\frac{1}{2}[I_{n-1}(z)+I_{n+1}(z)],~{}~{}I_{n}^{\prime}(z)=I_{n-1}(z)-\frac{n}{z}I_{n}(z),~{}~{}I_{n}^{\prime}(z)=\frac{n}{z}I_{n}(z)+I_{n+1}(z),$
(3.4)

see Jeffrey & Dai (2008, §§ 17) for these and further properties.

Equations (3.2) and (3.3) immediately imply important relations for the
Skellam’s PMF (3.1), namely

$P(X^{*}=-x)\ \overset{\eqref{bessel}}{=}\
e^{-\lambda_{1}-\lambda_{2}}\,\big{(}\lambda_{1}/\lambda_{2}\big{)}^{x/2-x}\,I_{x}\big{(}2\,\sqrt{\lambda_{1}\lambda_{2}}\big{)}\
=\ \big{(}\lambda_{2}/\lambda_{1}\big{)}^{x}\,P(X^{*}=x),$ (3.5)

and

$\begin{array}[]{@{}rl}P(X^{*}=x+1)\
\overset{\eqref{bessel2}}{=}&e^{-\lambda_{1}-\lambda_{2}}\,\big{(}\lambda_{1}/\lambda_{2}\big{)}^{(x+1)/2}\,\Big{(}I_{x-1}\big{(}2\,\sqrt{\lambda_{1}\lambda_{2}}\big{)}\
-\
\frac{x}{\sqrt{\lambda_{1}\lambda_{2}}}\,I_{x}\big{(}2\,\sqrt{\lambda_{1}\lambda_{2}}\big{)}\Big{)}\\\\[4.30554pt]
\ =&\big{(}\lambda_{1}/\lambda_{2}\big{)}\,P(X^{*}=x-1)\ -\
\big{(}x/\lambda_{2}\big{)}\,P(X^{*}=x).\end{array}$ (3.6)

These results shall be useful for computing partial moments, see the details
below.

The cumulative distribution function (CDF) of
$\textup{Sk}(\lambda_{1},\lambda_{2})$ is easily computed in practice, by
utilizing Equation (4) in Johnson (1959) on a connection to the non-central
$\chi^{2}$-distribution. Let $\chi^{2}(\nu,\tau)$ denote the non-central
$\chi^{2}$-distribution with $\nu$ degrees of freedom and non-centrality
parameter $\tau$ (whose CDF is provided by pchisq in R). Then, Equation (4) in
Johnson (1959) implies that

$P(X^{*}\leq x)\ =\ \left\\{\begin{array}[]{ll}P(Q\leq 2\lambda_{2})\text{
with }Q\sim\chi^{2}(-2x,2\lambda_{1})&\text{if }x\leq 0,\\\\[4.30554pt]
1-P(Q\leq 2\lambda_{1})\text{ with
}Q\sim\chi^{2}\big{(}2(x+1),2\lambda_{2}\big{)}&\text{if
}x>0.\end{array}\right.$ (3.7)

Finally, let us turn to the computation of partial moments. In view of (2.4),
for example, the partial moment
$E\big{(}X^{*}\,\mathds{1}_{\mathbb{N}}(X^{*})\big{)}$ is relevant to us,
while the second-order moment
$E\big{(}(X^{*})^{2}\,\mathds{1}_{\mathbb{N}}(X^{*})\big{)}$, also see (2.5),
is used for computing the partial variance. Note that the Skellam distribution
is characterized by the Stein equation

$E\big{(}X^{*}\,f(X^{*})\big{)}\ =\ \lambda_{1}\,E\big{(}f(X^{*}+1)\big{)}\ -\
\lambda_{2}\,E\big{(}f(X^{*}-1)\big{)}$ (3.8)

for any bounded function $f:\mathbb{Z}\to\mathbb{R}$, see Theorem 1 in
Balachandran et al. (2017). Such a Stein equation can be utilized to computing
sophisticated moment expressions, see Weiß & Aleksandrov (2022), and it is now
applied for deriving the relevant partial moments.

###### Proposition 1

Let $X^{*}\sim\textup{Sk}(\lambda_{1},\lambda_{2})$, then it has the partial
mean

$E\big{(}X^{*}\,\mathds{1}_{\mathbb{N}}(X^{*})\big{)}\ =\ \mu\,P(X^{*}\geq 0)\
+\ \lambda_{2}\,\big{(}P(X^{*}=0)+P(X^{*}=1)\big{)}$ (3.9)

and the second-order partial moment

$E\big{(}(X^{*})^{2}\,\mathds{1}_{\mathbb{N}}(X^{*})\big{)}\ =\
(\sigma^{2}+\mu^{2})\,P(X^{*}\geq 1)\ +\ \lambda_{2}\,\mu\,P(X^{*}=1)\ +\
\lambda_{1}\,(1+\mu)\,P(X^{*}=0).$ (3.10)

The proof of Proposition 1 is provided by Appendix A.2. The partial-moment
expressions in Proposition 1 are evaluated by using Equations (3.1) and (3.7)
on the Skellam’s PMF and CDF, respectively. Note that (3.9) and (3.10) imply
the variance of $X^{*}\,\mathds{1}_{\mathbb{N}}(X^{*})$.

In practice, it may be advantageous to reparametrize
$\textup{Sk}(\lambda_{1},\lambda_{2})$ by the mean $\mu$ and an additional
dispersion parameter $\delta>0$, defined by the additive decomposition
$\sigma^{2}=|\mu|+\delta$. Thus, the reparametrized version
$\textup{Sk}^{*}(\mu,\delta)$ is related to the original version
$\textup{Sk}(\lambda_{1},\lambda_{2})$ by

$\lambda_{1}\ =\
\tfrac{1}{2}\,\big{(}|\mu|+\mu+\delta\big{)},\qquad\lambda_{2}\ =\
\tfrac{1}{2}\,\big{(}|\mu|-\mu+\delta\big{)}.$ (3.11)

Note that the limit $\delta\to 0$ leads to (the negative of) a
$\textup{Poi}\big{(}|\mu|\big{)}$-variate if $\mu>0$ ($\mu<0$), i. e., the
Poisson distribution is a boundary case of the Skellam distribution.

###### Example 1

In the case $\mu=0$ (i. e., $\lambda_{1}=\lambda_{2}=\delta/2$), the Skellam
distribution is symmetric around 0 with
$P(X^{*}=x)=e^{-\delta}\,I_{x}(\delta)$ and $P(X^{*}\geq
1)=\big{(}1-P(X^{*}=0)\big{)}/2$. Equations (3.9) and (3.10) simplify to

$\begin{array}[]{rl}E\big{(}X^{*}\,\mathds{1}_{\mathbb{N}}(X^{*})\big{)}\
=&\tfrac{1}{2}\,\delta\,\big{(}P(X^{*}=0)+P(X^{*}=1)\big{)}\ =\
\tfrac{1}{2}\,\delta\,e^{-\delta}\,\big{(}I_{0}(\delta)+I_{1}(\delta)\big{)}.\\\\[4.30554pt]
E\big{(}(X^{*})^{2}\,\mathds{1}_{\mathbb{N}}(X^{*})\big{)}\
=&\delta\,P(X^{*}\geq 1)\ +\ \tfrac{1}{2}\,\delta\,P(X^{*}=0)\ =\
\tfrac{1}{2}\,\delta.\end{array}$

Thus, we have

$\begin{array}[]{rl}V\big{(}X^{*}\,\mathds{1}_{\mathbb{N}}(X^{*})\big{)}\
=&\tfrac{1}{2}\,\delta-\tfrac{1}{4}\,\delta^{2}\,e^{-2\,\delta}\,\big{(}I_{0}(\delta)+I_{1}(\delta)\big{)}^{2}.\end{array}$

###  3.2 Skellam–Tobit INGARCH Model: Definition and Properties

Let us define the Skellam operator
$\mathcal{S}_{\mbox{\scriptsize\boldmath$\nu$}}(\mu)$ to generate a random
integer number according to $\textup{Sk}(\lambda_{1},\lambda_{2})$, where
$\mbox{\boldmath$\nu$}=(\lambda_{1},\lambda_{2}),\mu=\lambda_{1}-\lambda_{2}$.
Then, adapting the general Tobit INGARCH approach (2.2), we define the
_Skellam–Tobit INGARCH (STINGARCH) model_ by

$X_{t}\ =\ \max\big{\\{}0,\ X_{t}^{*}\big{\\}}\quad\text{with}\quad X_{t}^{*}\
=\ \mathcal{S}_{\mbox{\scriptsize\boldmath$\nu$}}(M_{t}),$ (3.12)

where $M_{t}=E(X_{t}^{*}|\mathcal{F}_{t-1})$ satisfies the linear INGARCH
equation (2.3). Note that Doukhan et al. (2021) considered a Skellam GARCH
model for $\mathbb{Z}$-valued time series based on modelling the conditional
variance, but here we will provide a method for analyzing non-negative count
time series based on modeling the conditional mean.

###### Theorem 1

Consider the STINGARCH process $\\{X_{t}\\}$ defined in (3.12). If

$\sum_{i=1}^{p}\max\\{0,\alpha_{i}\\}+\sum_{j=1}^{q}|\beta_{j}|<1,$ (3.13)

then $\\{X_{t}\\}$ exists and is stationary.

The proof of Theorem 1 is provided by Appendix A.3.

###### Remark 2

The above condition is comparable to the ones in other real-valued dynamic
Tobit models, such as de Jong & Herrera (2011), Michel & de Jong (2018), and
Bykhovskaya (2023). This condition is also comparable to the one for the
INGARCH model in Weiß et al. (2022). Due to the censoring mechanism, we do not
need to impose lower bounds on the parameters $\alpha_{i}$, $i=1,\ldots,p$.

In what follows, we consider the STINGARCH model based on the additive
reparametrization (3.11), i. e.,

$X_{t}\ =\ \max\big{\\{}0,\ X_{t}^{*}\big{\\}}\quad\text{with}\quad X_{t}^{*}\
=\ \mathcal{S}_{\delta}(M_{t}).$ (3.14)

Then, the conditional probability of $X_{t}^{\ast}$ given $M_{t}$ is

$\displaystyle P(X_{t}^{\ast}=x|M_{t})$ $\displaystyle={\rm
e}^{-|M_{t}|-\delta}\left(\frac{|M_{t}|+M_{t}+\delta}{|M_{t}|-M_{t}+\delta}\right)^{x/2}I_{x}\Big{(}\sqrt{\big{(}|M_{t}|+M_{t}+\delta\big{)}\big{(}|M_{t}|-M_{t}+\delta\big{)}}\Big{)}$
$\displaystyle=\left\\{\begin{array}[]{ll}\displaystyle{\rm
e}^{-M_{t}-\delta}\left(\frac{2M_{t}+\delta}{\delta}\right)^{x/2}\,I_{x}\Big{(}\sqrt{\delta\,(2M_{t}+\delta)}\Big{)}&\text{if
}M_{t}\geq 0,\\\ \displaystyle{\rm
e}^{M_{t}-\delta}\left(\frac{\delta}{-2M_{t}+\delta}\right)^{x/2}\,I_{x}\Big{(}\sqrt{\delta\,(-2M_{t}+\delta)}\Big{)}&\text{if
}M_{t}<0.\end{array}\right.$ (3.17)

For the reparametrized STINGARCH model (3.14), the limit $\delta\to 0$ leads
to the ReLU INGARCH model with conditional Poisson distribution
$\textup{Poi}\big{(}\max\\{0,M_{t}\\}\big{)}$.

###### Remark 3

As explained before, we use an additive reparametrization for the Skellam’s
relation between variance and mean. This differs from Carallo et al. (2024),
who use the multiplicative decomposition $\sigma^{2}=\tilde{\delta}\,|\mu|$
with $\tilde{\delta}>1$ (as $|\mu|<\sigma^{2}$ always holds). However, to
guarantee $\sigma^{2}>0$, such a multiplicative parametrization requires to
assume that $\lambda_{1}\neq\lambda_{2}$, which is unnecessarily restrictive.
It would also lead to rather different moment plots than those in Figure 1 in
Section 3.3 below.

It should be noted that the considered reparametrization does not affect
Theorem 1, because we established the existence and stationarity for the
original version of Skellam.

Figure 1: Plots (against $\mu$) of mean, variance, and dispersion ratio (from
left to right) of left-censored $\textup{Sk}^{*}(\mu,\delta)$-variate.

###  3.3 On Moments and Approximate Linearity

The Tobit INGARCH approach was proposed as an approximate solution for
achieving a conditionally linear count model allowing for negative ACF values.
This conjecture shall now to be verified. To compute the conditional moments
(2.4) and (2.5) for the STINGARCH model (3.14), the Skellam’s partial moments
in Proposition 1 are required. The resulting mean, variance, and dispersion
ratio of such left-censored $\textup{Sk}^{*}(\mu,\delta)$-variates are
illustrated by Figure 1. It can be seen that both the mean and the variance
behave roughly linear outside $\mu=0$, and the extent of linearity increases
with decreasing $\delta$. While the mean is rather stable for varying values
of $\delta$, the dispersion ratio increases with increasing $\delta$ (but is
roughly constant for varying $\mu$). For $\delta\to 0$, the dispersion ratio
converges to the constant 1 (for all $\mu$), whereas the limiting mean
(=limiting variance) equals $\max\\{0,\mu\\}$. This goes along with our above
statement that the limit $\delta\to 0$ leads to the Poisson ReLU INGARCH
model.

The approximate (piecewise) linearity of the left-censored Skellam’s mean in
Figure 1 (for small $\delta$ such as 0.25) leads us to conjecture that also
the STINGARCH model behaves similar to a truly linear model, because it
implies that the conditional mean $E(X_{t}|\mathcal{F}_{t-1})$ is
approximately linear in past observations, also recall Grunwald et al. (2000).
More precisely, we conjecture that the marginal mean $\mu$, variance
$\sigma^{2}$, and ACF $\rho(h)$ with time lag $h\in\mathbb{N}$ of a STINGARCH
process $\\{X_{t}\\}$ can be computed in good approximation according to the
same formula as for the truly linear, so-called “dispersed INGARCH” (DINGARCH)
model of Xu et al. (2012), see Appendix A.1 for details, but by also allowing
for negative parameter values (which are impossible for the DINGARCH model).
Here, the DINGARCH’s dispersion parameter $\eta$ is chosen to equal the
dispersion ratio
$I_{\mu,\delta}^{*}:=V\big{(}X^{*}\,\mathds{1}_{\mathbb{N}}(X^{*})\big{)}\,\big{/}\,E\big{(}X^{*}\,\mathds{1}_{\mathbb{N}}(X^{*})\big{)}$
for $X^{*}\sim\textup{Sk}^{*}(\mu,\delta)$, which is computed according to
Proposition 1. Let us illustrate this approximate approach for the
STINARCH$(1)$ and STINGARCH$(1,1)$ model. We propose the approximate scheme

$\begin{array}[]{l@{\qquad}l}(p,q)=(1,0):&\displaystyle\mu\ \approx\
\frac{\alpha_{0}}{1-\alpha_{1}},\quad\frac{\sigma^{2}}{\mu}\ \approx\
I_{\mu,\delta}^{*}\cdot\frac{1}{1-\alpha_{1}^{2}},\quad\rho(h)\ \approx\
\alpha_{1}^{h}\text{ for }h\geq 1,\\\\[12.91663pt]
(p,q)=(1,1):&\displaystyle\mu\ \approx\
\frac{\alpha_{0}}{1-\alpha_{1}-\beta_{1}},\quad\frac{\sigma^{2}}{\mu}\
\approx\
I_{\mu,\delta}^{*}\cdot\frac{1-(\alpha_{1}+\beta_{1})^{2}+\alpha_{1}^{2}}{1-(\alpha_{1}+\beta_{1})^{2}},\\\\[8.61108pt]
&\displaystyle\rho(h)\ \approx\
(\alpha_{1}+\beta_{1})^{h-1}\,\frac{\alpha_{1}\,\big{(}1-\beta_{1}(\alpha_{1}+\beta_{1})\big{)}}{1-(\alpha_{1}+\beta_{1})^{2}+\alpha_{1}^{2}}\text{
for }h\geq 1,\end{array}$ (3.18)

for moment computation in these cases. To check the quality of this
approximation and, thus, the extent of linearity of the STINGARCH models, we
compare the linear moment values being computed according to (3.18) with their
true values, which are computed numerically exactly for the STINARCH$(1)$
model (by utilizing its Markov property), and approximately for
STINGARCH$(1,1)$ (based on a simulated sample path of length $10^{6}$).
Results are tabulated in Appendix B; let us start with the STINARCH$(1)$ case.
For different choices of $\alpha_{0},\alpha_{1},\delta$, the true values
(columns “tob”) and anticipated linear values (columns “lin”, computed via
(3.18)) for marginal mean, dispersion ratio, and partial ACF (PACF) at lags
$h=1,2,3$ have been computed, where $\delta=0$ corresponds to the limiting
case of the Poi-ReLU INARCH$(1)$ model. The latter just serves as the
benchmark for the maximal extent of linearity within the INARCH$(1)$ approach,
but it is not relevant for practice because of the drawbacks outlined in
Section 1.

Table B.1 provides results for the case, where the marginal mean is close to 5
(“low counts”). Comparing the “tob-columns” to the “lin-columns”, we generally
see a good agreement, which further improves for decreasing $\delta$ (and
which is even better for larger mean, see Table B.2 for comparison). Except
for strong negative autocorrelation (as implied by $\alpha_{1}=-0.75$), the
choice $\delta=0.25$, for example, leads to linearity in close approximation.
For $\alpha_{1}=-0.75$, it should be noted that even the exactly piecewise-
linear Poi-ReLU INARCH$(1)$ model causes departures from linearity, so some
extent of non-linearity cannot be avoided in this extreme situation. It is
also worth pointing out that the tob- and lin-values of the dispersion ratio
$\sigma^{2}/\mu$ agree quite well even for larger $\delta$ (except the extreme
scenario $\alpha_{1}=-0.75$), i. e., the approximation
$I_{\mu,\delta}^{*}/(1-\alpha_{1}^{2})$ generally shows a good performance.

If we now turn to the corresponding STINGARCH$(1,1)$ results in Tables
B.3–B.4, where the ACF instead of the PACF is considered, we again see a very
close agreement throughout. The additional feedback term did not affect the
closeness to linearity compared to our previous STINARCH$(1)$ conclusions.
Thus, to sum up, if $\delta$ is chosen sufficiently small such as
$\delta=0.25$, the novel Tobit approach leads to approximately linear count
time series models that capture a wide range of negative ACF values.
Therefore, STINGARCH models are recommended as a kind of “workhorse” for
linear count processes with signed ACF values.

##  4 Estimation

Let $X_{1},\ldots,X_{n}$ be observations from model (3.12). Define the
parameter vector $\mbox{\boldmath$\theta$}=(\alpha_{0},\alpha_{1},\ldots$,
$\alpha_{p},\beta_{1},\ldots,\beta_{q},\delta)^{\top}$, and let
$\mbox{\boldmath$\vartheta$}=(\alpha_{0},\alpha_{1},\ldots,\alpha_{p},\beta_{1},\ldots,\beta_{q})^{\top}$,
i. e.,
$\mbox{\boldmath$\theta$}=(\mbox{\boldmath$\vartheta$}^{\top},\delta)^{\top}$.
Let
$\mbox{\boldmath$\vartheta$}_{0}=(\alpha_{00},\alpha_{10},\ldots,\alpha_{p0},\beta_{10},\ldots,\beta_{q0})^{\top}$
be the true value of $\vartheta$, and
$\mbox{\boldmath$\theta$}_{0}=(\mbox{\boldmath$\vartheta$}_{0}^{\top},\delta_{0})^{\top}$
be the one of $\theta$. The parameter space of $\theta$ is $\Theta$ and the
parameter space of $\vartheta$ is $\mbox{\boldmath$\Theta$}^{\ast}$, which are
compact subsets of $\mathbb{R}^{p+q+2}$ and $\mathbb{R}^{p+q+1}$,
respectively, and can be further refined when deriving large-sample properties
of the following considered estimators. We have to distinguish two general
estimation scenarios:

Scenario 1:

$\delta$ is considered as a tuning parameter that is specified in advance to
achieve approximate linearity. Recalling the discussion in Section 3.3,
setting $\delta:=0.25$ or smaller might be such a choice. Then, parameter
estimation refers to solely $\vartheta$.

Scenario 2:

$\delta$ is considered as a parameter of the data-generating process (DGP).
Then, parameter estimation refers to whole $\theta$.

In our opinion, scenario 1 is more relevant for practice. In view of George
Box’ famous words “all models are wrong, but some are useful”, we consider the
STINGARCH model as tool (workhorse) for describing count time series
exhibiting linear sample properties and negative dependencies. Here, the
parameter $\delta$ is used as a tuning parameter set by the user to control
the extent of linearity, similar to the parameter $c$ of the softplus INGARCH
model of Weiß et al. (2022) mentioned in Section 1. The user may also define
multiple candidate models with different extents of linearity. But for
completeness, we also investigate scenario 2 in some detail.

We consider two estimation methods: the maximum likelihood estimator (MLE) in
Section 4.1, and the censored least absolute deviations estimator (CLADE) in
Section 4.2. The second method is expected to be robust if the underlying
distribution is mis-specified, and it can be applied to estimation scenario 1.
The first method can be used for both estimation scenarios.

###  4.1 Maximum Likelihood Estimation

First, we consider the MLE of $\mbox{\boldmath$\theta$}_{0}$ (i. e.,
estimation scenario 2). Based on (3.2), define
$p_{t}(x)=P(X_{t}^{\ast}=x|M_{t})$ and $f_{t}(x)=P(X_{t}^{\ast}\leq x|M_{t})$.
Then, the conditional log-likelihood function is

$L(\mbox{\boldmath$\theta$})=\sum_{t=1}^{n}l_{t}(\mbox{\boldmath$\theta$})\quad\text{with}\quad
l_{t}(\mbox{\boldmath$\theta$})\ =\
\mathds{1}_{\mathbb{N}}(X_{t})\,\log{p_{t}(X_{t})}\ +\
\mathds{1}_{\\{0\\}}(X_{t})\,\log{f_{t}(0)},$ (4.1)

where Equations (3.1) and (3.7) are used for the PMF and CDF computations,
respectively. If $X_{t}>0$, then

$\displaystyle l_{t}(\mbox{\boldmath$\theta$})\ =\
\left\\{\begin{array}[]{ll}\displaystyle-
M_{t}-\delta+\tfrac{1}{2}\,X_{t}\,\Big{[}\log(2M_{t}+\delta)-\log\delta\Big{]}+\log
I_{X_{t}}\Big{(}\sqrt{\delta(2M_{t}+\delta)}\Big{)}&\text{if }M_{t}\geq
0,\\\\[4.30554pt] \displaystyle
M_{t}-\delta-\tfrac{1}{2}\,X_{t}\,\Big{[}\log(-2M_{t}+\delta)-\log\delta\Big{]}+\log
I_{X_{t}}\Big{(}\sqrt{\delta(-2M_{t}+\delta)}\Big{)}&\text{if
}M_{t}<0.\end{array}\right.$

If $X_{t}=0$, then $l_{t}(\mbox{\boldmath$\theta$})=\log{f_{t}(0)}$. The MLE
$\hat{}\mbox{\boldmath$\theta$}$ is defined as
$\hat{}\mbox{\boldmath$\theta$}={\arg\max}_{\mbox{\scriptsize\boldmath$\theta$}}L(\mbox{\boldmath$\theta$})$.
If MLE is applied to estimation scenario 1, then $\delta$ is specified in
advance, and the maximization of (4.1) is done with respect to $\vartheta$.

The above likelihood function is cumbersome (for both scenarios), so we need
some numerical algorithms to find the MLE $\hat{}\mbox{\boldmath$\theta$}$.
The score function and the Hessian matrix are defined as

$\mbox{\boldmath$S$}_{n}(\mbox{\boldmath$\theta$})\ =\
\sum_{t=1}^{n}\frac{\partial
l_{t}(\mbox{\boldmath$\theta$})}{\partial\mbox{\boldmath$\theta$}}\qquad\text{and}\qquad\mbox{{{H}}}_{n}(\mbox{\boldmath$\theta$})\
=\
-\sum_{t=1}^{n}\frac{\partial^{2}l_{t}(\mbox{\boldmath$\theta$})}{\partial\mbox{\boldmath$\theta$}\partial\mbox{\boldmath$\theta$}^{\top}},$

respectively. The MLE $\hat{\mbox{\boldmath$\theta$}}$ also is a solution to
the equation
$\mbox{\boldmath$S$}_{n}(\mbox{\boldmath$\theta$})=\mbox{\boldmath$0$}$.
Detailed expressions for the partial derivatives are discussed in Appendix
A.4.

For establishing consistency and asymptotic normality of MLE, one could verify
the regularity conditions in Theorems 2.7 and 3.1 in Newey & McFadden (1994)
or Theorems 4.1.1 and 4.1.3 in Amemiya (1985). But because of the complex
expressions for the modified Bessel function of the first kind in the
likelihood function, it is hardly possible to do so. Note that the
aforementioned theorems would imply that the asymptotic variance of the MLE is
$\mbox{{{U}}}^{-1}\mbox{{{V}}}\mbox{{{U}}}^{-1}/n$, where

$\mbox{{{V}}}\ =\ E\left(\frac{\partial
l_{t}(\mbox{\boldmath$\theta$})}{\partial\mbox{\boldmath$\theta$}}\,\frac{\partial
l_{t}(\mbox{\boldmath$\theta$})}{\partial\mbox{\boldmath$\theta$}^{\top}}\right),\qquad\mbox{{{U}}}\
=\
-E\left(\frac{\partial^{2}l_{t}(\mbox{\boldmath$\theta$})}{\partial\mbox{\boldmath$\theta$}\partial\mbox{\boldmath$\theta$}^{\top}}\right),$

and where these two matrices are evaluated at the true value
$\mbox{\boldmath$\theta$}=\mbox{\boldmath$\theta$}_{0}$. This implies to
approximate the asymptotic variance by
$(\hat{\mbox{{{U}}}}_{n}\hat{\mbox{{{V}}}}^{-1}_{n}\hat{\mbox{{{U}}}}_{n})^{-1}/n$,
where

$\hat{\mbox{{{V}}}}_{n}\ =\ \left.\frac{1}{n}\sum_{t=1}^{n}\frac{\partial
l_{t}(\mbox{\boldmath$\theta$})}{\partial\mbox{\boldmath$\theta$}}\,\frac{\partial
l_{t}(\mbox{\boldmath$\theta$})}{\partial\mbox{\boldmath$\theta$}^{\top}}\right|_{\mbox{\scriptsize\boldmath$\theta$}=\hat{}\mbox{\scriptsize\boldmath$\theta$}},\qquad\hat{\mbox{{{U}}}}_{n}\
=\
-\left.\frac{1}{n}\sum_{t=1}^{n}\frac{\partial^{2}l_{t}(\mbox{\boldmath$\theta$})}{\partial\mbox{\boldmath$\theta$}\partial\mbox{\boldmath$\theta$}^{\top}}\right|_{\mbox{\scriptsize\boldmath$\theta$}=\hat{}\mbox{\scriptsize\boldmath$\theta$}}.$

As a practically feasible solution, we recommend to approximate the MLE’s
standard errors by using the numerical Hessian of the maximized log-likelihood
function. More precisely, we use the square-roots of the diagonal of the
inverse Hessian, as it is common practice in applications of the ML approach.
The performance of these approximate standard errors is later investigated by
simulations in Section 4.3.

###  4.2 Censored Least Absolute Deviations Estimation

CLADE is a popular estimation method for censored models since Powell (1984),
also see de Jong & Herrera (2011), Bilias et al. (2019) and Bykhovskaya
(2023). As before, let
$\mbox{\boldmath$\vartheta$}=(\alpha_{0},\alpha_{1},\ldots$,
$\alpha_{p},\beta_{1},\ldots,\beta_{q})^{\top}$, and let
$\mbox{\boldmath$\vartheta$}_{0}=(\alpha_{00},\alpha_{10},\ldots,\alpha_{p0},\beta_{10},\ldots,\beta_{q0})^{\top}$
be its true value. The CLADE $\check{}\mbox{\boldmath$\vartheta$}$ is
formalized as

$\check{}\mbox{\boldmath$\vartheta$}=\arg\min_{\mbox{\scriptsize\boldmath$\vartheta$}\in\mbox{\scriptsize\boldmath$\Theta$}^{\ast}}\sum_{t=1}^{n}\big{|}X_{t}-\max\\{0,M_{t}\\}\big{|},$
(4.2)

where $M_{t}$ is defined in (2.3). Here, we consider the CLADE for
$\mbox{\boldmath$\vartheta$}_{0}$, not for $\mbox{\boldmath$\theta$}_{0}$,
because this estimator does not provide a first-round estimate for $\delta$.

For real-valued dynamic Tobit models, de Jong & Herrera (2011) established
consistency and asymptotic normality of CLADE based on Powell (1984), and
Bykhovskaya (2023) considered CLADE in a more complex case. But their analyses
cannot be generalized here, because their techniques use a error density,
which is absent for our model. This problem deserves further study in future
research.

###### Remark 4

In general, the minimization of absolute deviations in models with a positive
part is a non-convex problem, and designing numerical algorithms requires
special care. The method in Bilias et al. (2019) allows to find an exact
computation of $\check{}\mbox{\boldmath$\vartheta$}$, i. e., the CLADE
minimization problem (4.2) can be reformulated as the following mixed integer
programming problem:

$\displaystyle\min_{\\{\mbox{\scriptsize\boldmath$\vartheta$},(\gamma_{t},\phi_{t},sm_{t},sp_{t},t=1,\ldots,n)\\}}\sum_{t=1}^{n}(sm_{t}+sp_{t})\quad\text{subject
to}$ $\displaystyle X_{t}-\phi_{t}+sm_{t}-sp_{t}=0$
$\displaystyle\phi_{t}\geq\max\\{0,\ M_{t}\\},$
$\displaystyle\phi_{t}\leq\min\big{\\{}M_{t}+\omega_{t}(1-\gamma_{t}),\
\omega_{t}\gamma_{t}\big{\\}},$ $\displaystyle sm_{t}\geq 0,\quad sp_{t}\geq
0,\quad\gamma_{t}\in\\{0,1\\},$

where $t=1,\ldots,n$. Here, $\omega_{t}$ is an arbitrary parameter that must
be sufficiently large such that a solution to the system of inequalities
exists. In other words, the inequalities only represent the max-operator
exactly if $\omega_{t}$ is large enough, specifically if $\omega_{t}>|M_{t}|$.
The auxiliary choice variables $sm_{t}$ and $sp_{t}$ are the under-achievement
and over-achievement positive slack variables for the CLAD model, and the
positive variable $\phi_{t}$ is equal to $\max\\{0,M_{t}\\}$. The optimization
is done with respect to the binary choice variables $\gamma_{t}$ and the
continuous choice variables $sm_{t}$, $sp_{t}$, $\phi_{t}$, and $\vartheta$.
See Bilias et al. (2019) for more details about the method.

For the STINARCH case (i. e., all $\beta_{j}$s in (2.3) are zero, so $M_{t}\
=\ \alpha_{0}+\sum_{i=1}^{p}\alpha_{i}\,X_{t-i}$ and
$\mbox{\boldmath$\vartheta$}=(\alpha_{0},\alpha_{1},\ldots$,
$\alpha_{p})^{\top}$), a possible alternative to (4.2) could be to use squared
deviations instead of absolute deviations, thus leading to a conditional least
squares (CLS) estimator. The ordinary CLS estimator, however, obtained by
minimizing $\sum_{t=1}^{n}\big{(}X_{t}-M_{t}\big{)}^{2}$ over
$\mbox{\boldmath$\vartheta$}\in\mbox{\boldmath$\Theta$}^{\ast}$, will not lead
a consistent estimator, which can be seen from a simple example: If $p=q=0$
with $M_{t}=\alpha_{0}$, then we have the estimator
$\tilde{\alpha}_{0}=\sum_{t=1}^{n}X_{t}/n$. Hence, $\tilde{\alpha}_{0}$
converges to $E(X_{t})\geq 0$ using the law of large numbers, but
$\alpha_{0}\in\mathbb{R}$, so $\tilde{\alpha}_{0}$ is inconsistent. A similar
discussion for other models can also be found in Bykhovskaya (2023).
Therefore, we consider the following CLS estimator
$\breve{\mbox{\boldmath$\vartheta$}}$:

$\displaystyle\breve{\mbox{\boldmath$\vartheta$}}=\arg\min_{\mbox{\scriptsize\boldmath$\vartheta$}\in\mbox{\scriptsize\boldmath$\Theta$}^{\ast}}\sum_{t=1}^{n}\big{(}X_{t}-\max\\{0,M_{t}\\}\big{)}^{2}.$

According to Proposition 1, we know that $E(X_{t}|\mathcal{F}_{t-1})\neq
M_{t}$ or $\max\\{0,M_{t}\\}$, so the above CLS estimator is not the common
one (where the conditional expectation is used to construct the objective
function). Note that $\max\\{0,M_{t}\\}\neq M_{t}\mathds{1}_{\\{X_{t}>0\\}}$,
and the objective function is discontinuous with respect to the parameters. So
we cannot obtain an explicit expression for the CLS estimator, the theoretic
properties of which need to be further studied. From the simulation experiment
in Table B.5 below, we find that the CLS estimator
$\breve{\mbox{\boldmath$\vartheta$}}$ shows a satisfactory finite-sample
performance in practice.

###  4.3 Results from Simulation Study

To analyze the performance of the different estimation approaches, we did a
comprehensive simulation study. For various parametrizations of STINARCH$(1)$
and STINGARCH$(1,1)$ processes, for different sample sizes $n$, and for both
estimation scenarios, 1,000 stationary count time series were generated. The
estimates were determined via numerical optimization (by using R’s optim),
where for the STINGARCH$(1,1)$ model, we used the initialization
$M_{0}:=\alpha_{0}$. Relevant results are tabulated in Appendix B. Table B.5
for the STINARCH$(1)$ process refers to estimation scenario 1, i. e., $\delta$
is considered as a tuning parameter and specified in advance as $0.25$ such
that the model behaves nearly linear. In estimation scenario 1, the model is
taken as a working tool to handle linear counts, i. e., one does not
necessarily expect that the DGP is also STINARCH$(1)$ with exactly
$\delta=0.25$. Therefore, we considered STINARCH$(1)$ DGPs with true
$\delta\in\\{0,0.25,1\\}$, to check the robustness of estimates with respect
to a possible mis-specification of $\delta$. Table B.5 shows that in any case,
bias and standard errors (SEs) decrease with increasing $n$, confirming the
consistency of the estimates in all cases. But there are also notable
differences, especially for low $n$. For negative $\alpha_{1}$, the CLADE
$\check{\alpha}_{0}$ shows a clear negative bias, whereas the MLE
$\hat{\alpha}_{0}$ and CLS $\breve{\alpha}_{0}$ exhibit an only mild bias. For
positive $\alpha_{1}$, the situation is reversed, the MLE $\hat{\alpha}_{0}$
and CLS $\breve{\alpha}_{0}$ have a notable positive bias. The
$\alpha_{1}$-estimates, by contrast, have quite similar properties for all
estimation approaches and parametrizations. Table B.5 is complemented by Table
B.6, where the performance of the MLE’s approximative SE (derived from the
log-likelihood’s Hessian) is investigated. Especially for positive
$\alpha_{1}$, we see a rather good agreement between the simulated SEs and the
mean of the approximated SEs. For negative $\alpha_{1}$ and small $n$, we
observe a slight exceedance of the true SE, but the discrepancy vanishes with
increasing $n$. It is also worth noting that the computation of the
approximate SEs was always possible in the case of positive $\alpha_{1}$,
while the Hessian was not invertible for negative $\alpha_{1}$ in 1–2 % of all
simulation runs.

Table B.7 again considers the STINARCH$(1)$ model, but now regarding
estimation scenario 2, i. e., $\delta$ is estimated together with
$(\alpha_{0},\alpha_{1})$. In this case, only MLE is applicable. It can be
seen that the bias and SE properties of $(\hat{\alpha}_{0},\hat{\alpha}_{1})$
are essentially the same as if $\delta$ is specified (recall Table B.5), but
the performance of $\hat{\delta}$ is interesting. Although $\hat{\delta}$ is
visibly biased for, say, $n\leq 250$, the bias decreases and becomes
negligible for $n\geq 500$. Also the standard errors of $\hat{\delta}$
decrease with increasing $n$, but they are still rather larger even for
$n=1000$ (compared to the true value of $\delta$). These difficulties in
estimating $\delta$ are plausible in view of Section 3.3, where we recognized
the nearly linear properties of STINGARCH models for a range of
$\delta$-values (impeding the distinguishability of different $\delta$). As
before, we also analyzed the performance of the approximate SEs. While the
approximate SEs still agree reasonably well with the simulated SEs in the mean
(also the ones for $\delta$), the inclusion of $\delta$ into parameter
estimation causes another problem: the percentage of simulation runs where the
Hessian cannot be inverted increases considerably in some cases, namely with
decreasing $n$ and $\delta$. For $\delta=0.25$, it takes values about 40 % if
$n=100$, but decreases below 10 % for $n=1000$. For $\delta=1$, the
percentages of non-invertibility are much smaller, namely around 12 % for
$n=100$ and between 0–2 % for $n=1000$.

Let us return to estimation scenario 1 (using $\delta$ as a tuning parameter),
which is more relevant for practice in our opinion. In Tables B.8–Table B.9,
estimation results for STINGARCH$(1,1)$ models are shown. Since MLE and CLS
again showed analogous bias and SE properties, we omit the CLS results this
time to save some space. Furthermore, we now focus on sample sizes $n\geq
250$, because it is well known for INGARCH-type models that the inclusion of
the feedback term $M_{t-1}$ increases the sample-size requirement for valid
estimation. For $(\alpha_{1},\beta_{1})=(-0.25,0.45)$, the CLADE of
$\alpha_{0}$ succeeds in terms of bias and SE, but otherwise, there is hardly
any clear advantage but often a deterioration compared to MLE. Therefore,
taking all results together, the MLE appears to be the preferable estimator
for STINGARCH processes. In addition, it allows to approximate the SEs from
the log-likelihood’s Hessian. As can be seen from Table B.10, the simulated SE
and the mean of the approximate ones agree quite well in nearly any case
(sometimes slight deviations if $n=250$). It is also worth noting that we
hardly experienced any problems in inverting the Hessian; if there was a non-
zero percentage at all, then it was usually much below 1 %.

##  5 Possible Extensions

First, we state some possible extensions based on the unbounded INGARCH
modeling framework developed in Section 2. In Section 3, the Skellam
distribution is used as the underlying conditional distribution for
illustration, but other choices are possible as well. For examples, if
stronger dispersion is required, one my use the generalized Poisson difference
distribution discussed in Carallo et al. (2024) instead, which includes the
Skellam distribution as a special case. In addition, unlike Agosto et al.
(2016) where a non-negative link function is used for possible covariates,
such covariates could be easily and naturally included in (2.3).

Second, the Tobit approach proposed in Section 2 is not only relevant for
unbounded INGARCH models, it might also be used for different integer-valued
time series models such as those from the INAR (integer-valued autoregressive)
family. While the ordinary INAR models for count time series (see Weiß, 2018,
Chapters 2–3) are well-defined only for positive parameter values and can,
thus, only describe positive ACF values, also extensions to
$\mathbb{Z}$-valued INAR-type models have been proposed in the literature. The
first contribution appears to be the INARS$(p)$ model (S=signed) by Kim & Park
(2008), which relies on a model recursion of the form

$X_{t}\ =\ \alpha_{1}\odot X_{t-1}\ +\ldots+\ \alpha_{p}\odot X_{t-{p}}\ +\
\epsilon_{t},$ (5.1)

where “$\odot$” denotes the signed binomial thinning operator, and the
$(\epsilon_{t})$ are i. i. d. $\mathbb{Z}$-valued innovations. “$\odot$” is
defined by requiring the conditional distribution $\alpha\odot X\big{|}X$ to
be
$\textup{sgn}(\alpha)\cdot\textup{sgn}(X)\cdot\textup{Bin}\big{(}|X|,|\alpha|\big{)}$,
where $-1<\alpha<1$, where $X$ is a $\mathbb{Z}$-valued random variable, and
where $\textup{sgn}(z)=1$ for $z\geq 0$ and $-1$ otherwise. The special case
of the INARS$(1)$ model according to (5.1) has been investigated by Andersson
& Karlis (2014) together with Skellam-distributed innovations, whereas Alzaid
& Omair (2014) discuss a modification allowing for a Skellam marginal
distribution.

Applying the Tobit approach of Section 2 to the INARS$(p)$ model (5.1), we
propose the $p$th-order _Tobit INAR model_ defined by

$X_{t}\ =\ \max\big{\\{}0,\ X_{t}^{*}\big{\\}},\qquad X_{t}^{*}\ =\
\alpha_{1}\odot X_{t-1}\ +\ldots+\ \alpha_{p}\odot X_{t-{p}}\ +\
\epsilon_{t}.$ (5.2)

Unlike the ordinary INAR$(p)$ model for count time series (see Alzaid & Al-
Osh, 1990; Du & Li, 1991), the AR-parameters
$\alpha_{1},\ldots,\alpha_{p}\in(-1,1)$ are allowed to take negative values
(obviously, if all $\alpha_{1},\ldots,\alpha_{p}>0$, then (5.2) reduces to the
ordinary INAR$(p)$ recursion). The observations $\\{X_{t}\\}$ in (5.2) are
counts, and $X_{t}^{*}$ is allowed to take negative values only to capture the
effect of negative AR-parameters. But as we are actually interested in count
time series modeling, we assume the i. i. d. innovations $(\epsilon_{t})$ to
be counts as well. A more detailed investigation of such Tobit INAR$(p)$
models and their abilities to handle negative autocorrelations appears to be
an interesting direction for future research.

Third, another possible extension refers to the case of bounded counts, i. e.,
where the range of $\\{X_{t}\\}$ is bounded from above by some
$N\in\mathbb{N}$. Both INAR- and INGARCH-type models have been proposed for
such a case, see McKenzie (1985); Ristić et al. (2016) as well as the survey
in Weiß (2018, Chapters 3–4). Here, the additional difficulty arises from the
fact that the range $\\{0,\ldots,N\\}$ is bounded from two sides. Thus,
instead of using the ReLU function for constructing a Tobit model, it suggests
itself to use a type of clipped ReLU (cReLU) function (see Cai et al., 2017)
instead. Let us focus on the BINGARCH case (B=bounded) for illustration. Using
again the recursive scheme (2.3) for the conditional mean
$M_{t}=E(X_{t}^{*}|\mathcal{F}_{t-1})$, the _Tobit BINGARCH model_ is defined
by

$X_{t}\ =\ \max\big{\\{}0,\ \min\\{N,\
X_{t}^{*}\\}\big{\\}}\quad\text{with}\quad X_{t}^{*}\ =\
\mathcal{S}_{\mbox{\scriptsize\boldmath$\nu$}}(M_{t}),$ (5.3)

where $\mathcal{S}_{\mbox{\scriptsize\boldmath$\nu$}}(\mu)$,
$\mu\in\mathbb{R}$, again denotes a $\mathbb{Z}$-valued random operator. An
analogous construction might be used for adapting the Tobit INAR model (5.2)
to the case of bounded counts. Again, future research should study the
properties and abilities of these models.

##  6 Real Applications

In what follows, we illustrate the application of the novel STINGARCH model by
a couple of real-world data examples. In addition, these data examples serve
for discussing some of the potential extensions sketched in Section 5.

Figure 2: Plot of lottery winners $x_{1},\ldots,x_{666}$ and their sample ACF
$\hat{\rho}(k)$, see Section 6.1.

###  6.1 Lottery Winners Data

The main lottery game in Germany is “Lotto 6 aus 49”, where 6 out of 49
numbers (plus one bonus number) are randomly drawn twice a week (on Wednesdays
and Saturdays). Information on past lottery draws is available on
https://www.westlotto.de/lotto-6aus49/gewinnzahlen/gewinnzahlen.html, from
which we obtained the number of lottery wins in the second prize category (“6
numbers”) for each lottery draw between Saturday, May 4, 2013, and Wednesday,
September 18, 2019. The resulting count time series of length $n=666$ is
plotted in Figure 2. The corresponding sample ACF shows an alternating pattern
of positive and negative ACF values, which decrease only slowly with
increasing time lag $k$. Because of this kind of a “long memory”, the
STINGARCH$(1,1)$ model appears to be a plausible candidate for these data.
Since it is not immediately clear that we are dealing with linear time series
data here, we consider estimation scenario 2 (recall Section 4), i. e., the
parameter $\delta$ is included into ML estimation. The resulting ML estimates
are shown in the first block of Table 2. $\hat{\delta}$ takes a rather large
value (similar for the remaining candidate models), which indicates a non-
linear dependence structure (recall Figure 2). The strong value for
$\hat{\beta}_{1}$, in turn, appears plausible in view of the apparent “long
memory”. Since $|\hat{\beta}_{1}|$ is rather close to 1, we also fitted the
STINGARCH$(1,2)$ model to the data. Its estimates are more balanced, but
Akaike’s and Bayesian information criterion (AIC and BIC, respectively) are
worse. We also considered the Poisson softplus INGARCH$(1,1)$ (Poi-
spINGARCH$(1,1)$) model of Weiß et al. (2022) with $c=1$ as a further
competitor (which is discussed in some more detail later in Section 6.2), but
it clearly performs worse. Computing the standardized Pearson residuals (Weiß,
2018, Sect. 2.4) for checking the adequacy of the ML-fitted models, their
sample ACF values in Table 2 are not significant, indicating that both fitted
models capture the ACF structure very well. Also their means are very close to
zero, i. e., the fitted models well describe the mean of the data (sample mean
2.443). The variances of the residuals in Table 2 are only a little larger
than one, and this slight exceedance is caused by a few extreme counts in
Figure 2: excluding, for example, the residuals corresponding to $x_{405}=21$
and $x_{609}=20$, the variances reduce to 1.180 and 1.215. Thus, also the
dispersion of the data appears to be covered quite well by the fitted models.

Table 1: ML estimation for lottery winners data: estimates and approximate
standard errors (SEs), information criteria.

Model | $\hat{\delta}$ | SE | $\hat{\alpha}_{0}$ | SE | $\hat{\alpha}_{1}$ | SE | $\hat{\beta}_{1}$ | SE | $\hat{\beta}_{2}$ | SE | AIC | BIC
---|---|---|---|---|---|---|---|---|---|---|---|---
STINGARCH$(1,1)$ | 3.619 | (0.471) | 4.442 | (0.197) | -0.008 | (0.008) | -0.990 | (0.007) |  |  | 2782.7 | 2800.7
STINGARCH$(1,2)$ | 3.595 | (0.483) | 1.630 | (0.617) | -0.044 | (0.017) | -0.316 | (0.137) | 0.642 | (0.136) | 2800.1 | 2822.6
Poi-spINGARCH$(1,1)$ |  |  | 5.023 | (0.118) | -0.016 | (0.011) | -1.099 | (0.012) |  |  | 2923.3 | 2936.8
Model | $\hat{\delta}$ | SE | $\hat{\alpha}_{0}$ | SE | $\hat{\gamma}_{1}$ | SE |  |  |  |  | AIC | BIC
ST-regression | 3.458 | (0.124) | 1.239 | (0.185) | 2.031 | (0.471) |  |  |  |  | 2783.3 | 2796.8

Table 2: Model diagnostics by Pearson residuals for lottery winners data.

Model | Mean | Variance | ACF(1) | ACF(2) | ACF(3) | ACF(4) | ACF(5)
---|---|---|---|---|---|---|---
STINGARCH$(1,1)$ | -0.007 | 1.351 | 0.012 | 0.008 | -0.043 | 0.042 | -0.009
STINGARCH$(1,2)$ | -0.021 | 1.368 | 0.025 | 0.006 | -0.036 | 0.050 | -0.005
ST-regression | -0.019 | 1.358 | -0.009 | 0.010 | -0.060 | 0.057 | -0.026

Revisiting the ACF in Figure 2, another explanation for the alternating
pattern could be periodicity with period 2, which indeed appears plausible in
view of the draws on Wednesdays and Saturdays. Hence, as a further candidate
for modelling the lottery winners data, let us use the weekday information as
a covariate, say, $z_{t}=1$ ($=0$) if the $t$th lottery draw is on Saturday
(Wednesday). As already mentioned in Section 5, covariates
$z_{t,1},\ldots,z_{t,r}$ are easily included into the model by just adding the
sum $\sum_{k=1}^{r}\gamma_{k}\,z_{t,k}$ on the right-hand side of Equation
(2.3). For the ML-fit including the covariate term $+\gamma_{1}\,z_{t}$
referring to the weekday, the autoregressive parameters are not significant
anymore, and the model can be reduced to the equation
$M_{t}=\alpha_{0}+\gamma_{1}\,z_{t}$. Its estimation results (including
$\hat{\delta}$) are shown in the second block (“ST-regression”) of Table 2.
The positive value for $\hat{\gamma}_{1}$ implies that the event “6 numbers”
is more likely on Saturdays than on Wednesdays, which can be explained by the
fact that more people participate in Saturday’s lottery. While the AIC is
slightly worse than for the STINGARCH$(1,1)$ model, the BIC is better. The
checks for model adequacy based on the Pearson residuals, see Table 2,
indicate a similarly well performance: their mean is close to zero, their
variance (without $x_{405},x_{609}$: 1.209) is close to one, and their ACF is
not significant. Hence, the user may decide between the candidate models, e.
g., based on their interpretability.

Figure 3: Plot of process yields $x_{1},\ldots,x_{70}$ and their sample PACF
$\hat{\rho}_{\textup{p}}(k)$, see Section 6.2.

###  6.2 Chemical Process Data

As our second example, let us revisit the chemical process data in Appendix
A.3 of O’Donovan (1983, “Time series 4.1”), which consist of $n=70$
consecutive yields from a batch chemical process. The time series plot in
Figure 3 shows a “spiky behaviour”. Because of the apparent AR$(1)$-like PACF
with strongly negative lag-1 autocorrelation, $\hat{\rho}(1)=-0.588$, these
data were considered by Weiß et al. (2022, Section 5.2) as an illustrative
example for their Poi-spINARCH$(1)$ model. As explained in Section 1, our
STINARCH$(1)$ constitutes a natural competitor to that model, so we shall now
compare their model fits. In addition, recall that our Tobit approach is not
limited to unbounded INGARCH-like counts, as it is the case for the softplus
INGARCH model, but can also be adapted to different model families, see
Section 5. Also this topic shall be made a subject of discussion in this
section.

Table 3: ML estimation for process yields data: estimates and approximate
standard errors (SEs), information criteria.

Model |  $c$ or $\delta$ | $\hat{\alpha}_{0}$ | SE | $\hat{\alpha}_{1}$ | SE | AIC | BIC
---|---|---|---|---|---|---|---
Poi-spINARCH$(1)$ | 1 | 79.783 | (4.820) | -0.603 | (0.094) | 486.94 | 491.44
STINARCH$(1)$ | 0.25 | 79.767 | (4.833) | -0.602 | (0.094) | 486.81 | 491.31
Poi-TINARS$(1)$ | — | 73.792 | (6.066) | -0.482 | (0.121) | 489.18 | 493.68

Table 4: Model diagnostics by Pearson residuals for process yields data.

Model | Mean | Variance | ACF(1) | ACF(2) | ACF(3) | ACF(4) | ACF(5)
---|---|---|---|---|---|---|---
Poi-spINARCH$(1)$ | 0.000 | 1.142 | -0.064 | -0.054 | 0.119 | 0.044 | 0.082
STINARCH$(1)$ | 0.000 | 1.136 | -0.065 | -0.054 | 0.118 | 0.044 | 0.082

Comparing the ML-estimation results for the Poi-spINARCH$(1)$ and
STINARCH$(1)$ models, see rows 1–2 in Table 4, we recognize a nearly identical
performance (with slightly better AIC and BIC for STINARCH$(1)$). This also
applies to the stochastic properties of both model fits. Regarding the Pearson
residuals in Table 4, both have a mean about zero, i. e., they nearly
perfectly capture the data’s mean (sample mean 49.686). The Pearson residuals’
variances only slightly exceed one, see Table 4, so both models only slightly
fall below the sample dispersion 1.706. Also the data’s autocorrelation
structure is described quite well by both model fits, because none of the
Pearson residuals’ ACF value is significantly different from zero. Note that
the STINARCH$(1)$’s CLADE are given by $\check{\alpha}_{0}=72.892$ and
$\check{\alpha}_{1}=-0.461$, whereas CLS leads to $\breve{\alpha}_{0}=79.071$
and $\breve{\alpha}_{1}=-0.588$. Hence, like in the simulation study of
Section 4.3, the CLS estimates are more close to the MLE. Finally, like the
Poi-spINARCH$(1)$ model’s tuning parameter $c$ was specified as 1, to
guarantee conditional linearity in close approximation, we specified the
STINARCH$(1)$’s tuning parameter as $\delta=0.25$ (as recommended in Section
3.3; estimation scenario 1). Actually, because of the small sample size
$n=70$, the estimation of $\delta$ does not work reliably here because of huge
SEs. But we noted that a variation of $\delta$ had only negligible effects on
the estimation results anyway.

In the previous paragraph, we recognized that Poi-spINARCH$(1)$ and
STINARCH$(1)$ lead to nearly equally well model fits. But as shown in Section
5, the Tobit approach is more flexible and might be applied to different model
families as well. As we are concerned with AR$(1)$-like unbounded counts, also
the Tobit INARS$(1)$ model (5.2) would be a plausible further candidate for
the chemical process data. Thus, let us also fit the Tobit INARS$(1)$ model
$X_{t}\,=\,\max\big{\\{}0,\,\alpha_{1}\odot X_{t-1}+\epsilon_{t}\big{\\}}$
with i. i. d. $\textup{Poi}(\alpha_{0})$-innovations (Poi-TINARS$(1)$) to the
data. Its transition probabilities, as required for likelihood computation,
are

$\begin{array}[]{rl}P(X_{t}=x>0\ |\ X_{t-1})\
=&\sum_{j=0}^{X_{t-1}}P\Big{(}\textup{Bin}\big{(}X_{t-1},|\alpha_{1}|\big{)}=j\Big{)}\cdot
P\Big{(}\epsilon_{t}=x-\textup{sgn}(\alpha_{1})\,j\Big{)},\\\\[4.30554pt]
P(X_{t}=0\ |\ X_{t-1})\
=&\sum_{j=0}^{X_{t-1}}P\Big{(}\textup{Bin}\big{(}X_{t-1},|\alpha_{1}|\big{)}=j\big{)}\cdot
P\Big{(}\epsilon_{t}\leq-\textup{sgn}(\alpha_{1})\,j\Big{)},\end{array}$

where “$\textup{Bin}(X_{t-1},|\alpha_{1}|)$” is used as an abbreviation for a
binomial random variable. But as can be seen from the third row in Table 4,
the fitted Poi-TINARS$(1)$ model performs worst. It does not capture the ACF
to full extent ($\hat{\alpha}_{1}=-0.482$ vs. $\hat{\rho}(1)=-0.588$), and the
fitted models dispersion ratio 2.257 is too large. So taking all results
together, the chemical process data are best explained by the STINARCH$(1)$
model.

Figure 4: Plot of first 365 air quality counts, and sample ACF $\hat{\rho}(k)$
of full data, see Section 6.3.

###  6.3 Air Quality Data

Also the third example serves to demonstrate the great flexibility of the
Tobit approach. For this purpose, let us revisit one of the time series
discussed by Liu et al. (2022), namely the daily air quality level (Dec.
2013–July 2019, so length $n=2\,068$) in Beijing. While air quality is
actually measured on an ordinal scale (from $s_{0}=\text{``excellent''}$ to
$s_{5}=\text{``severely polluted''}$), Liu et al. (2022) used a “rank-count
approach” for modelling these data, i. e., the ordinal $Y_{t}$ at time $t$ is
substituted by $s_{X_{t}}$ with the count $X_{t}$ having the upper bound
$N=5$. For the sake of readability, the plot in Figure 4 shows only the first
365 observations, as well as the sample ACF $\hat{\rho}(k)$ of the full data.
Liu et al. (2022) fitted an (exactly linear) INGARCH$(1,1)$ model to the data,
referred to as “ZOBPAR$(1,1)$”, where the conditional distribution is a
truncated Poisson distribution with additional one-inflation. Because of being
exactly linear, the ZOBPAR’s model parameters are not allowed to become
negative. However, looking at the sample ACF in Figure 4, we observe a rather
strong decline after lag 1, and this might be better explained if negative
parameter values are permitted. Thus, we use the Skellam Tobit BINGARCH$(1,1)$
(STBINGARCH$(1,1)$) model (5.3) with $N=5$ from Section 5 as a competitor,
which we also equip with a further parameter $\kappa\in[0;1)$ for capturing
one-inflation:

$P(X_{t}=x\ |\ X_{t-1},\ldots)\ =\ \kappa\cdot\mathds{1}_{\\{1\\}}(x)\ +\
(1-\kappa)\cdot\left\\{\begin{array}[]{ll}P\big{(}\mathcal{S}_{\delta}(M_{t})\leq
0\big{)}&\text{if }x=0,\\\\[2.15277pt]
P\big{(}\mathcal{S}_{\delta}(M_{t})=x\big{)}&\text{if }1\leq x\leq
N-1,\\\\[2.15277pt] P\big{(}\mathcal{S}_{\delta}(M_{t})\geq N\big{)}&\text{if
}x=N.\end{array}\right.$

When specifying $\delta$ in our data applications, we routinely tried out
candidate models with different levels $\delta$ of linearity. In this case, we
recognized that we can further improve upon our default choice $\delta=0.25$
by choosing the tuning parameter very close to 0; we used $\delta=0.01$. The
results of ML estimation are summarized in Table 6.

Table 5: ML estimation for air quality data: estimates and approximate
standard errors (SEs), information criteria.

Model | $\delta$ | $\hat{\alpha}_{0}$ | SE | $\hat{\alpha}_{1}$ | SE | $\hat{\beta}_{1}$ | SE | $\kappa$ | SE | AIC | BIC
---|---|---|---|---|---|---|---|---|---|---|---
ZOBPAR$(1,1)$ | — | 0.592 | (0.059) | 0.756 | (0.038) | 0.000 | (0.040) | 0.128 | (0.015) | 5721.6 | 5744.1
STBINGARCH$(1,1)$ | 0.01 | 0.787 | (0.067) | 0.699 | (0.030) | -0.127 | (0.037) | 0.118 | (0.015) | 5745.5 | 5768.0

Table 6: Model diagnostics by Pearson residuals for air quality data.

Model | Mean | Variance | ACF(1) | ACF(2) | ACF(3) | ACF(4) | ACF(5)
---|---|---|---|---|---|---|---
ZOBPAR$(1,1)$ | 0.010 | 0.801 | 0.060 | -0.059 | -0.028 | -0.023 | 0.014
STBINGARCH$(1,1)$ | -0.034 | 0.752 | 0.018 | -0.087 | -0.043 | -0.031 | 0.011

If solely looking at the information criteria, then Table 6 clearly votes for
the ZOBPAR$(1,1)$ model of Liu et al. (2022). But in fact, a more in-depth
analysis is advisable. It can be seen that the $\beta_{1}$-estimate of the
ZOBPAR model falls on the lower bound 0 of its admissible range, whereas
STBINGARCH leads to a significantly negative $\hat{\beta}_{1}$. If comparing
the ACFs of the fitted models to the sample ACF, namely
$0.521,0.272,0.144,\ldots$ (ZOBPAR) vs. $0.527,0.228,0.100,\ldots$
(STBINGARCH) vs. $0.539,0.214,0.073,\ldots$ (data), we recognize that
STBINGARCH better explains the strong decline after lag 1. This is also
confirmed by an analysis of the Pearson residuals’ ACF in Table 6, where the
ZOBPAR model has significant values at lags 1 and 2, and STBINGARCH only at
lag 2. Furthermore, both models lead to a Pearson residuals’ variance below
one this time, i. e., both models include too much dispersion, and now
STBINGARCH does worse. The reason becomes clear by comparing the full marginal
distributions $p(0),\ldots,p(5)$:

$\begin{array}[]{rlllllll}&0.220,&0.375,&0.183,&0.117,&0.068,&0.037&\text{(ZOBPAR)}\\\
\text{vs.}&0.229,&0.374,&0.182,&0.107,&0.057,&0.051&\text{(STBINGARCH)}\\\
\text{vs.}&0.174,&0.383,&0.222,&0.131,&0.065,&0.025&\text{(data).}\end{array}$

Both models over-estimate the zero probability at the cost of under-estimating
some central probabilities, but STBINGARCH also over-estimates the probability
for reaching the upper bound $N=5$. Thus, the STBINGARCH’s disadvantage
regarding the information criteria seems to be caused by this deviation of the
marginal distribution, whereas it performs better regarding the
autocorrelation structure.

##  7 Conclusions

To give an alternative solution to the log and softplus approaches for
achieving negative ACF values within the class of INGARCH models, the Tobit
approach is proposed, and the Skellam–Tobit INGARCH model for unbounded counts
is studied in detail. Several extensions of this approach to other situations
are discussed as well. Three real-data examples are analyzed to illustrate the
usefulness of the new approach.

There are some suggestions for possible future research. First, one could work
out the TINGARCH model based on the generalized Poisson difference
distribution (recall Sections 5) or other distributions. In particular,
closed-form expressions for distribution and moments of the corresponding
censored random variable $X^{*}\,\mathds{1}_{\mathbb{N}}(X^{*})$ should be
derived. Second, more studies regarding the Tobit INAR models for unbounded or
bounded counts mentioned in Sections 5 are needed. Third, similar to the
arguments in Weiß et al. (2022), we can apply the new approach to many real
datasets and find some causal relationship like Chen et al. (2018), who
investigated the relation between human influenza cases and air pollution.

#### Acknowledgements

The authors thank Dr. Tobias Möller (Helmut Schmidt University, Hamburg) for
the collection and preparation of the lottery winners data used in Section
6.1. Zhu’s work is supported by National Natural Science Foundation of China
(No. 12271206), and Science and Technology Research Planning Project of Jilin
Provincial Department of Education (No. JJKH20231122KJ).

## References

  * Agosto et al. (2016) Agosto, A., Cavaliere, G., Kristensen, D., Rahbek, A. (2016) Modeling corporate defaults: Poisson autoregressions with exogenous covariates (PARX). Journal of Empirical Finance 38, 640–663.
  * Alzaid & Al-Osh (1990) Alzaid, A.A., Al-Osh, M.A. (1990) An integer-valued pth-order autoregressive structure (INAR($p$)) process. Journal of Applied Probability 27, 314–324.
  * Alzaid & Omair (2014) Alzaid, A.A., Omair, M.A. (2014) Poisson difference integer valued autoregressive model of order one. Bulletin of the Malaysian Mathematical Science Society 37, 465–485.
  * Amemiya (1985) Amemiya, T. (1985) Advanced Econometrics. Harvard University Press, Cambridge.
  * Andersson & Karlis (2014) Andersson, J., Karlis, D. (2014) A parametric time series model with covariates for integers in $\mathbb{Z}$. Statistical Modelling 14, 135–156.
  * Balachandran et al. (2017) Balachandran, P., Kolaczyk, E.D., Viles, W.D. (2017) On the propagation of low-rate measurement error to subgraph counts in large networks. Journal of Machine Learning Research 18, 1–33.
  * Bilias et al. (2019) Bilias, Y., Florios, K., Skouras, S. (2019) Exact computation of Censored Least Absolute Deviations estimator. Journal of Econometrics 212, 584–606.
  * Bykhovskaya (2023) Bykhovskaya, A. (2023) Time series approach to the evolution of networks: prediction and estimation. Journal of Business & Economic Statistics 41, 170–183.
  * Cai et al. (2017) Cai, Z., He, X., Sun, J., Vasconcelos, N. (2017) Deep learning with low precision by half-wave Gaussian quantization. 2017 IEEE Conference on Computer Vision and Pattern Recognition (CVPR), Honolulu, HI, 5406–5414.
  * Carallo et al. (2024) Carallo, G., Casarin, R., Robert, C.P. (2024) Generalized Poisson difference autoregressive processes. International Journal of Forecasting, forthcoming, https://doi.org/10.1016/j.ijforecast.2023.11.009
  * Chen et al. (2018) Chen, C.W.S., Hsieh, Y., Su, H., Wu, J.J. (2018) Causality test of ambient fine particles and human influenza in Taiwan: Age group-specific disparity and geographic heterogeneity. Environment International 111, 354–361.
  * de Jong & Herrera (2011) de Jong, R., Herrera, A.M. (2011) Dynamic censored regression and the open market desk reaction function. Journal of Business & Economic Statistics 29, 228–237.
  * Doukhan et al. (2021) Doukhan, P., Mamode Khan, N., Neumann, M.H. (2021) Mixing properties of integer-valued GARCH processes. ALEA 18, 401–420.
  * Du & Li (1991) Du, J.-G., Li, Y. (1991) The integer-valued autoregressive (INAR($p$)) model. Journal of Time Series Analysis 12, 129–142.
  * Ferland et al. (2006) Ferland, R., Latour, A., Oraichi, D. (2006) Integer-valued GARCH processes. Journal of Time Series Analysis 27, 923–942.
  * Fokianos et al. (2009) Fokianos, K., Rahbek, A., Tjøstheim, D. (2009) Poisson autoregression. Journal of the American Statistical Association 104, 1430–1439.
  * Fokianos & Tjøstheim (2011) Fokianos, K., Tjøstheim, D. (2011) Log-linear Poisson autoregression. Journal of Multivariate Analysis 102, 563–578.
  * Goldberger (1964) Goldberger, A.S. (1964) Econometric Theory. John Wiley & Sons, New York.
  * Grunwald et al. (2000) Grunwald, G., Hyndman, R.J., Tedesco, L., Tweedie, R.L. (2000) Non-Gaussian conditional linear AR(1) models. Australian and New Zealand Journal of Statistics 42, 479–495.
  * Honoré (1993) Honoré, B.E. (1993) Orthogonality conditions for Tobit models with fixed effects and lagged dependent variables. Journal of Econometrics 59, 35–61.
  * Honoré & Hu (2004) Honoré, B.E., Hu, L. (2004) Estimation of cross sectional and panel data censored regression models with endogeneity. Journal of Econometrics 122, 293–316.
  * Jeffrey & Dai (2008) Jeffrey, A., Dai, H. (2008) Handbook of Mathematical Formulas and Integrals. 4th edition, Academic Press, San Diego.
  * Johnson (1959) Johnson, N.L. (1959) On an extension of the connexion between Poisson and $\chi^{2}$ distributions. Biometrika 46, 352–363.
  * Johnson et al. (2005) Johnson, N.L., Kemp, A.W., Kotz, S. (2005) Univariate Discrete Distributions. 3rd edition, John Wiley & Sons, New Jersey.
  * Karlis & Ntzoufras (2006) Karlis, D., Ntzoufras, I. (2006) Bayesian analysis of the differences of count data. Statistics in Medicine 25, 1885–1905.
  * Karlis & Ntzoufras (2009) Karlis, D., Ntzoufras, I. (2009) Bayesian modelling of football outcomes: using the Skellam’s distribution for the goal difference. IMA Journal of Management Mathematics 20, 133–145.
  * Kim & Park (2008) Kim, H.-Y., Park, Y. (2008) A non-stationary integer-valued autoregressive model. Statistical Papers 49, 485–502.
  * Koopman et al. (2017) Koopman, S.J., Lit, R., Lucas, A. (2017) Intraday stochastic volatility in discrete price changes: The dynamic skellam model. Journal of the American Statistical Association 112, 1490–1503.
  * Liu et al. (2022) Liu, M., Zhu, F., Zhu, K. (2022) Modeling normalcy-dominant ordinal time series: An application to air quality level. Journal of Time Series Analysis 43, 460–478.
  * Maddala (1983) Maddala, G.S. (1983) Limited-Dependent and Qualitative Variables in Econometrics. Cambridge University Press, Cambridge.
  * McKenzie (1985) McKenzie, E. (1985) Some simple models for discrete variate time series. Water Resources Bulletin 21, 645–650.
  * Michel & de Jong (2018) Michel, J., de Jong, R.M. (2018) Mixing properties of the dynamic Tobit model with mixing errors. Economics Letters 162, 112–115.
  * Newey & McFadden (1994) Newey, W.K., McFadden, D. (1994) Large sample estimation and hypothesis testing. In R.F. Engle and D. MacFadden (eds): Handbook of Econometrics, Vol. 4, North-Holland, Amsterdam, 2111–2245.
  * O’Donovan (1983) O’Donovan, T.M. (1983) Short Term Forecasting: An Introduction to the Box–Jenkins Approach. John Wiley & Sons, Chichester.
  * Powell (1984) Powell, J.L. (1984) Least absolute deviations estimation for the censored regression model. Journal of Econometrics 25, 303–325.
  * Ristić et al. (2016) Ristić, M.M., Weiß, C.H., Janjić, A.D. (2016) A binomial integer-valued ARCH model. International Journal of Biostatistics 12, 20150051.
  * Tobin (1958) Tobin, J. (1958) Estimation of relationships for limited dependent variables. Econometrica 26, 24–36.
  * Weiß (2009) Weiß, C.H. (2009) Modelling time series of counts with overdispersion. Statistical Methods and Applications 18, 507–519.
  * Weiß (2018) Weiß, C.H. (2018) An Introduction to Discrete-valued Time Series. John Wiley & Sons, Chichester.
  * Weiß & Aleksandrov (2022) Weiß, C.H., Aleksandrov, B. (2022) Computing (bivariate) Poisson moments using Stein–Chen identities. The American Statistician 76, 10–15.
  * Weiß et al. (2022) Weiß, C.H., Zhu, F., Hoshiyar, A. (2022) Softplus INGARCH models. Statistica Sinica 32, 1099–1120.
  * Xu et al. (2012) Xu, H.-Y., Xie, M., Goh, T.N., Fu, X. (2012) A model for integer-valued time series with conditional overdispersion. Computational Statistics and Data Analysis 56, 4229–4242.
  * Zhu (2012) Zhu, F. (2012) Modeling overdispersed or underdispersed count data with generalized Poisson integer-valued GARCH models. Journal of Mathematical Analysis and Applications 389, 58–71.

##  Appendix A Appendices

###  A.1 About INGARCH models

For INGARCH$(p,q)$ models, the count $X_{t}$ being generated at time $t$ has
the conditional mean $M_{t}$ satisfying the linear recursive scheme (1.1) with
$a_{0}>0$ and $a_{1},\ldots,a_{p},b_{1},\ldots,b_{q}\geq 0$. The default
choice is a conditional Poisson distribution, i. e., $X_{t}$, conditioned on
$X_{t-1},\ldots$, is $\textup{Poi}(M_{t})$-distributed. Using the linearity of
the conditional mean, the unconditional mean equals (Ferland et al., 2006)

$\mu\ =\ \frac{a_{0}}{1-\sum_{i=1}^{p}a_{i}-\sum_{j=1}^{q}b_{j}},$ (A.1)

and variance and autocovariances satisfy a set of Yule–Walker equations (Weiß,
2009):

$\begin{array}[]{@{}rl}\gamma(0)\ =&\mu\ +\ \gamma_{M}(0),\qquad\gamma_{M}(0)\
=\ \sum\limits_{i=1}^{p}a_{i}\,\gamma(i)\ +\
\sum\limits_{j=1}^{q}b_{j}\,\gamma_{M}(j),\\\\[11.62495pt] \gamma(k)\
=&\sum\limits_{i=1}^{p}a_{i}\,\gamma(|k-i|)\ +\
\sum\limits_{j=1}^{\min{\\{k-1,q\\}}}b_{j}\,\gamma(k-j)\ +\
\sum\limits_{j=k}^{q}b_{j}\,\gamma_{M}(j-k),\\\\[11.62495pt] \gamma_{M}(k)\
=&\sum\limits_{i=1}^{\min{\\{k,p\\}}}a_{i}\,\gamma_{M}(|k-i|)\ +\
\sum\limits_{i=k+1}^{p}a_{i}\,\gamma(i-k)\ +\
\sum\limits_{j=1}^{q}b_{j}\,\gamma_{M}(|k-j|),\\\ \end{array}$ (A.2)

for $k\geq 1$, where $\gamma(h):=Cov[X_{t},X_{t-h}]$ and
$\gamma_{M}(h):=Cov[M_{t},M_{t-h}]$. Because of the conditional Poisson
distribution, the conditional dispersion ratio
$V(X_{t}|\mathcal{F}_{t-1})\big{/}E(X_{t}|\mathcal{F}_{t-1})=1$ constantly
(the unconditional distribution exhibits overdispersion, i. e.,
$\sigma^{2}/\mu>1$). To allow for a more flexible dispersion structure,
several extensions of the basic Poisson INGARCH model have been developed in
the literature, see Weiß (2018) for a survey. For this article, especially the
extension by Xu et al. (2012) is relevant, who defined the DINGARCH model by
introducing an additional dispersion parameter $\eta\geq 1$ into the
conditional variance:

$V(X_{t}|\mathcal{F}_{t-1})\ =\ \eta\,E(X_{t}|\mathcal{F}_{t-1}).$ (A.3)

So the Poisson INGARCH model corresponds to the choice $\eta=1$, while
specific examples for $\eta>1$ are given by the negative-binomial INGARCH
model of Xu et al. (2012) and the generalized-Poisson INGARCH model of Zhu
(2012). As a result, the mean and ACF of the DINGARCH model are the same as
for the Poisson INGARCH model, but the variance is inflated by the factor
$\eta$.

###  A.2 Proof of Proposition 1

To derive (3.9), we set $f(x)=\mathds{1}_{\mathbb{N}}(x)$ in (3.8) and obtain

$\begin{array}[]{rl}E\big{(}X^{*}\,\mathds{1}_{\mathbb{N}}(X^{*})\big{)}\
=&\lambda_{1}\,E\big{(}\mathds{1}_{\mathbb{N}}(X^{*}+1)\big{)}\ -\
\lambda_{2}\,E\big{(}\mathds{1}_{\mathbb{N}}(X^{*}-1)\big{)}\\\\[3.87498pt]
=&\lambda_{1}\,P(X^{*}\geq 0)\ -\ \lambda_{2}\,P(X^{*}\geq 2)\\\\[3.87498pt]
=&\mu\,P(X^{*}\geq 0)\ +\
\lambda_{2}\,\big{(}P(X^{*}=0)+P(X^{*}=1)\big{)}.\end{array}$

Next, let us turn to the second-order partial moment (3.10). In a first step,
we apply (3.8) with $f(x)=x\,\mathds{1}_{\mathbb{N}}(x)$, leading to

$\begin{array}[]{@{}l}E\big{(}(X^{*})^{2}\,\mathds{1}_{\mathbb{N}}(X^{*})\big{)}\
=\ \lambda_{1}\,E\big{(}(X^{*}+1)\,\mathds{1}_{\mathbb{N}}(X^{*}+1)\big{)}\ -\
\lambda_{2}\,E\big{(}(X^{*}-1)\,\mathds{1}_{\mathbb{N}}(X^{*}-1)\big{)}\\\\[3.87498pt]
\quad=\ \lambda_{1}\,E\big{(}X^{*}\,\mathds{1}_{\mathbb{N}}(X^{*}+1)\big{)}\
-\ \lambda_{2}\,E\big{(}X^{*}\,\mathds{1}_{\mathbb{N}}(X^{*}-1)\big{)}\ +\
\lambda_{1}\,P(X^{*}\geq 0)\ +\ \lambda_{2}\,P(X^{*}\geq 2).\end{array}$

In a second step, we apply (3.8) with $f(x)=\mathds{1}_{\mathbb{N}}(x+1)$ or
$f(x)=\mathds{1}_{\mathbb{N}}(x-1)$, respectively, leading to

$\begin{array}[]{rl}E\big{(}X^{*}\,\mathds{1}_{\mathbb{N}}(X^{*}+1)\big{)}\
=&\lambda_{1}\,P(X^{*}\geq-1)\ -\ \lambda_{2}\,P(X^{*}\geq 1)\\\\[1.93748pt] \
=&\lambda_{1}\,P(X^{*}=-1)+\lambda_{1}\,P(X^{*}=0)\ +\
(\lambda_{1}-\lambda_{2})\,P(X^{*}\geq 1)\\\\[1.93748pt] \
\overset{\eqref{skellam2}}{=}&\lambda_{2}\,P(X^{*}=1)+\lambda_{1}\,P(X^{*}=0)\
+\ (\lambda_{1}-\lambda_{2})\,P(X^{*}\geq 1),\\\\[3.87498pt]
E\big{(}X^{*}\,\mathds{1}_{\mathbb{N}}(X^{*}-1)\big{)}\
=&\lambda_{1}\,P(X^{*}\geq 1)\ -\ \lambda_{2}\,P(X^{*}\geq 3)\\\\[1.93748pt]
=&\lambda_{2}\,P(X^{*}=1)\ +\ \lambda_{2}\,P(X^{*}=2)\ +\
(\lambda_{1}-\lambda_{2})\,P(X^{*}\geq 1)\\\\[1.93748pt] \
\overset{\eqref{skellam3}}{=}&(\lambda_{2}-1)\,P(X^{*}=1)\ +\
\lambda_{1}\,P(X^{*}=0)\ +\ (\lambda_{1}-\lambda_{2})\,P(X^{*}\geq
1).\end{array}$

So altogether, we get

$\begin{array}[]{@{}l}E\big{(}(X^{*})^{2}\,\mathds{1}_{\mathbb{N}}(X^{*})\big{)}\
=\ \lambda_{1}\,P(X^{*}\geq 0)\ +\ \lambda_{2}\,P(X^{*}\geq 2)\ +\
(\lambda_{1}-\lambda_{2})^{2}\,P(X^{*}\geq 1)\\\\[1.93748pt] \qquad+\
\lambda_{1}\lambda_{2}\,P(X^{*}=1)+\lambda_{1}^{2}\,P(X^{*}=0)\ -\
\lambda_{2}\,(\lambda_{2}-1)\,P(X^{*}=1)-\lambda_{1}\lambda_{2}\,P(X^{*}=0)\\\\[3.87498pt]
\quad=\ \lambda_{1}\,P(X^{*}\geq 0)\ +\ \lambda_{2}\,P(X^{*}\geq 1)\ +\
(\lambda_{1}-\lambda_{2})^{2}\,P(X^{*}\geq 1)\\\\[1.93748pt] \qquad+\
\lambda_{2}\,(\lambda_{1}-\lambda_{2})\,P(X^{*}=1)\ +\
\lambda_{1}\,(\lambda_{1}-\lambda_{2})\,P(X^{*}=0)\\\\[3.87498pt] \quad=\
(\sigma^{2}+\mu^{2})\,P(X^{*}\geq 1)\ +\ \lambda_{2}\,\mu\,P(X^{*}=1)\ +\
\lambda_{1}\,(1+\mu)\,P(X^{*}=0).\end{array}$

This completes the proof of Proposition 1.

###  A.3 Proof of Theorem 1

We just need to show the existence and stationarity of the process
$\\{X_{t}^{\ast}\\}$. Carallo et al. (2024) considered a process similar to
$\\{X_{t}^{\ast}\\}$, see their equation (5), where they assume that
$\alpha_{1},\ldots,\alpha_{p},\beta_{1},\ldots,\beta_{q}$ are non-negative but
$\alpha_{0}\in\mathbb{R}$. The non-negativity restriction on the parameters is
redundant, as it does not affect the proof of Carallo et al. (2024).

For $M_{t}$ defined in our (2.3), we have

$\displaystyle\textstyle E(X_{t}^{\ast}|\mathcal{F}_{t-1})=M_{t}\ =\
\alpha_{0}+\sum\limits_{i=1}^{p}\alpha_{i}X_{t-i}^{\ast}\mathds{1}_{\mathbb{N}}(X_{t-i}^{\ast})+\sum\limits_{j=1}^{q}\beta_{j}M_{t-j}.$

Based on the sign of $X_{t-i}^{\ast}$, $i=1,\ldots,p$, $M_{t}$ can have
$2^{p}$ different cases. We consider the special case
$M_{t}=\alpha_{0}+\sum_{i=1}^{p}\alpha_{i}X_{t-i}^{\ast}+\sum_{j=1}^{q}\beta_{j}M_{t-j}$.
Using arguments similar to those in Section 3.1 in Carallo et al. (2024), we
know that if condition (3.13) is satisfied, then the existence and
stationarity of the process $\\{X_{t}^{\ast}\\}$ holds. For remaining cases,
where at least one of the summands in $\sum_{i=1}^{p}\alpha_{i}X_{t-i}^{\ast}$
vanishes, we can obtain similar conditions. But condition (3.13) is the
strictest, so it is a sufficient condition for all cases.

###  A.4 Derivatives of Log-Likelihood Function (4.1)

In what follows, the derivatives of the log-likelihood function
$L(\mbox{\boldmath$\theta$})=\sum_{t=1}^{n}l_{t}(\mbox{\boldmath$\theta$})$
according to (4.1) are discussed. If $X_{t}>0$ and $M_{t}>0$, then

$\displaystyle\frac{\partial
l_{t}(\mbox{\boldmath$\theta$})}{\partial\mbox{\boldmath$\theta$}^{*}}$
$\displaystyle\ =\
\left(-1+\frac{X_{t}}{a_{t}}+\frac{\sqrt{\delta}}{\sqrt{a_{t}}}\,b_{t}\right)\frac{\partial
M_{t}}{\partial\mbox{\boldmath$\theta$}^{*}},$ $\displaystyle\frac{\partial
l_{t}(\mbox{\boldmath$\theta$})}{\partial\delta}$ $\displaystyle\ =\
-1+\frac{X_{t}}{2a_{t}}-\frac{X_{t}}{2\delta}+\frac{M_{t}+\delta}{\sqrt{\delta
a_{t}}}\,b_{t},$
$\displaystyle\frac{\partial^{2}l_{t}(\mbox{\boldmath$\theta$})}{\partial\mbox{\boldmath$\theta$}^{*}\partial\mbox{\boldmath$\theta$}^{*\top}}$
$\displaystyle\ =\
\frac{1}{a_{t}}\left(-\frac{2X_{t}}{a_{t}}-\frac{\sqrt{\delta}}{\sqrt{a}_{t}}\,b_{t}+\delta
c_{t}\right)\frac{\partial
M_{t}}{\partial\mbox{\boldmath$\theta$}^{*}}\,\frac{\partial
M_{t}}{\partial\mbox{\boldmath$\theta$}^{*\top}}+\left(-1+\frac{X_{t}}{a_{t}}+\frac{\sqrt{\delta}}{\sqrt{a_{t}}}\,b_{t}\right)\frac{\partial^{2}M_{t}}{\partial\mbox{\boldmath$\theta$}^{*}\partial\mbox{\boldmath$\theta$}^{*\top}},$
$\displaystyle\frac{\partial^{2}l_{t}(\mbox{\boldmath$\theta$})}{\partial^{2}\delta^{2}}$
$\displaystyle\ =\
-\frac{X_{t}}{2a_{t}^{2}}+\frac{X_{t}}{2\delta^{2}}-\frac{M_{t}^{2}}{(\delta
a_{t})^{3/2}}\,b_{t}+\frac{(M_{t}+\delta)^{2}}{\delta a_{t}}\,c_{t},$
$\displaystyle\frac{\partial^{2}l_{t}(\mbox{\boldmath$\theta$})}{\partial\mbox{\boldmath$\theta$}^{*}\partial\delta}$
$\displaystyle\ =\
\frac{1}{a_{t}}\left(-\frac{X_{t}}{a_{t}}+\frac{M_{t}}{\sqrt{\delta
a_{t}}}\,b_{t}+(M_{t}+\delta)\,c_{t}\right)\frac{\partial
M_{t}}{\partial\mbox{\boldmath$\theta$}^{*}},$

where

$a_{t}\ =\ 2M_{t}+\delta,\qquad b_{t}\ =\
\frac{I_{X_{t}}^{\prime}}{I_{X_{t}}},\qquad c_{t}\ =\
\frac{I_{X_{t}}^{\prime\prime}I_{X_{t}}-(I_{X_{t}}^{\prime})^{2}}{I_{X_{t}}^{2}}.$

Here, the functions $I_{X_{t}},I_{X_{t}}^{\prime}$ and
$I_{X_{t}}^{\prime\prime}$ omit their argument $\sqrt{\delta\,a_{t}}$ for
brevity. From (3.4), we can calculate the derivatives required for $b_{t}$ and
$c_{t}$ as follows:

$\displaystyle I_{X_{t}}^{\prime}$ $\displaystyle\ =\
\frac{X_{t}}{\sqrt{\delta a_{t}}}\,I_{X_{t}}+I_{X_{t}+1},$ $\displaystyle
I_{X_{t}}^{\prime\prime}$ $\displaystyle\ =\
\frac{1}{2}(I_{X_{t}-1}^{\prime}+I_{X_{t}+1}^{\prime})\ =\
I_{X_{t}}+\frac{1}{2\sqrt{\delta
a_{t}}}\,\Big{[}(X_{t}-1)\,I_{X_{t}-1}-(X_{t}+1)\,I_{X_{t}+1}\Big{]}.$

For $i,i^{\prime}=1,\ldots,p,j,j^{\prime}=1,\ldots,q$ and $i\neq
i^{\prime},j\neq j^{\prime}$, we have

$\displaystyle\frac{\partial M_{t}}{\partial\alpha_{0}}$ $\displaystyle\ =\
1+\sum_{l=1}^{q}\beta_{l}\,\frac{\partial
M_{t-l}}{\partial\alpha_{0}},\quad\frac{\partial M_{t}}{\partial\alpha_{i}}\
=\ X_{t-i}+\sum_{l=1}^{q}\beta_{l}\,\frac{\partial
M_{t-l}}{\partial\alpha_{i}},\quad\frac{\partial M_{t}}{\partial\beta_{j}}\ =\
M_{t-j}+\sum_{l=1}^{q}\beta_{l}\,\frac{\partial M_{t-l}}{\partial\beta_{j}},$
$\displaystyle\frac{\partial^{2}M_{t}}{\partial\alpha_{0}^{2}}$
$\displaystyle\ =\
\sum_{l=1}^{q}\beta_{l}\,\frac{\partial^{2}M_{t-l}}{\partial\alpha_{0}^{2}},\qquad\frac{\partial^{2}M_{t}}{\partial\alpha_{i}^{2}}\
=\
\sum_{l=1}^{q}\beta_{l}\,\frac{\partial^{2}M_{t-l}}{\partial\alpha_{i}^{2}},\qquad\frac{\partial^{2}M_{t}}{\partial\alpha_{0}\partial\alpha_{i}}\
=\
\sum_{l=1}^{q}\beta_{l}\,\frac{\partial^{2}M_{t-l}}{\partial\alpha_{0}\partial\alpha_{i}},$
$\displaystyle\frac{\partial^{2}M_{t}}{\partial\alpha_{i}\partial\alpha_{i^{\prime}}}$
$\displaystyle\ =\
\sum_{l=1}^{q}\beta_{l}\,\frac{\partial^{2}M_{t-l}}{\partial\alpha_{i}\partial\alpha_{i^{\prime}}},\qquad\frac{\partial^{2}M_{t}}{\partial\alpha_{0}\partial\beta_{j}}\
=\ \frac{\partial
M_{t-j}}{\partial\alpha_{0}}+\sum_{l=1}^{q}\beta_{l}\,\frac{\partial^{2}M_{t-l}}{\partial\alpha_{0}\partial\beta_{j}},$
$\displaystyle\frac{\partial^{2}M_{t}}{\partial\alpha_{i}\partial\beta_{j}}$
$\displaystyle\ =\ \frac{\partial
M_{t-j}}{\partial\alpha_{i}}+\sum_{l=1}^{q}\beta_{l}\,\frac{\partial^{2}M_{t-l}}{\partial\alpha_{i}\partial\beta_{j}},\qquad\frac{\partial^{2}M_{t}}{\partial\beta_{j}^{2}}\
=\ 2\,\frac{\partial
M_{t-j}}{\partial\beta_{j}}+\sum_{l=1}^{q}\beta_{l}\,\frac{\partial^{2}M_{t-l}}{\partial\beta_{j}^{2}},$
$\displaystyle\frac{\partial^{2}M_{t}}{\partial\beta_{j}\partial\beta_{j^{\prime}}}$
$\displaystyle\ =\ \frac{\partial
M_{t-j}}{\partial\beta_{j^{\prime}}}+\frac{\partial
M_{t-j^{\prime}}}{\partial\beta_{j}}+\sum_{l=1}^{q}\beta_{l}\,\frac{\partial^{2}M_{t-l}}{\partial\beta_{j}\partial\beta_{j^{\prime}}}.$

Similarly, with an abuse of notations, if $X_{t}>0$ and $M_{t}<0$, then

$\displaystyle\frac{\partial
l_{t}(\mbox{\boldmath$\theta$})}{\partial\mbox{\boldmath$\theta$}^{*}}$
$\displaystyle\ =\
\left(1+\frac{X_{t}}{a_{t}}-\frac{\sqrt{\delta}}{\sqrt{a_{t}}}\,b_{t}\right)\frac{\partial
M_{t}}{\partial\mbox{\boldmath$\theta$}^{*}},$ $\displaystyle\frac{\partial
l_{t}(\mbox{\boldmath$\theta$})}{\partial\delta}$ $\displaystyle\ =\
-1-\frac{X_{t}}{2a_{t}}+\frac{X_{t}}{2\delta}-\frac{M_{t}-\delta}{\sqrt{\delta
a_{t}}}\,b_{t},$
$\displaystyle\frac{\partial^{2}l_{t}(\mbox{\boldmath$\theta$})}{\partial\mbox{\boldmath$\theta$}^{*}\partial\mbox{\boldmath$\theta$}^{*\top}}$
$\displaystyle\ =\
\frac{1}{a_{t}}\left(\frac{2X_{t}}{a_{t}}+\frac{\sqrt{\delta}}{\sqrt{a}_{t}}\,b_{t}+\delta
c_{t}\right)\frac{\partial
M_{t}}{\partial\mbox{\boldmath$\theta$}^{*}}\frac{\partial
M_{t}}{\partial\mbox{\boldmath$\theta$}^{*\top}}\ +\
\left(1+\frac{X_{t}}{a_{t}}-\frac{\sqrt{\delta}}{\sqrt{a_{t}}}\,b_{t}\right)\frac{\partial^{2}M_{t}}{\partial\mbox{\boldmath$\theta$}^{*}\partial\mbox{\boldmath$\theta$}^{*\top}},$
$\displaystyle\frac{\partial^{2}l_{t}(\mbox{\boldmath$\theta$})}{\partial^{2}\delta^{2}}$
$\displaystyle\ =\
\frac{X_{t}}{2a_{t}^{2}}-\frac{X_{t}}{2\delta^{2}}-\frac{M_{t}^{2}}{(\delta
a_{t})^{3/2}}\,b_{t}+\frac{(M_{t}-\delta)^{2}}{\delta a_{t}}\,c_{t},$
$\displaystyle\frac{\partial^{2}l_{t}(\mbox{\boldmath$\theta$})}{\partial\mbox{\boldmath$\theta$}^{*}\partial\delta}$
$\displaystyle\ =\
\frac{1}{a_{t}}\left(-\frac{X_{t}}{a_{t}}+\frac{M_{t}}{\sqrt{\delta
a_{t}}}\,b_{t}+(M_{t}-\delta)\,c_{t}\right)\frac{\partial
M_{t}}{\partial\mbox{\boldmath$\theta$}^{*}}.$

Now, $a_{t}=-2M_{t}+\delta$, and $b_{t}$ and $c_{t}$ are defined similarly as
before. Here, the omitted variable in the functions
$I_{X_{t}},I_{X_{t}}^{\prime}$ and $I_{X_{t}}^{\prime\prime}$ in $b_{t}$ and
$c_{t}$ is $\sqrt{\delta\,a_{t}}$.

If $X_{t}=0$, then the expression for $l_{t}(\mbox{\boldmath$\theta$})$
involves the sum of the conditional probabilities, but each term has been
discussed above.

##  Appendix B Tables

Table B.1: Moment properties of STINARCH$(1)$ model, see Section 3.3: true
values (“tob”) vs. linear approximation (“lin”), mean close to 5.

|  |  | $\mu$ | $\sigma^{2}/\mu$ |  PACF$(1)$ |  PACF$(2)$ |  PACF$(3)$
---|---|---|---|---|---|---|---
$\alpha_{0}$ | $\alpha_{1}$ | $\delta$ | tob | lin | tob | lin | tob | lin | tob | lin | tob | lin
8.75 | -0.75 | 1 | 5.044 | 5.000 | 2.303 | 2.710 | -0.698 | -0.750 | 0.024 | 0.000 | 0.007 | 0.000
|  | 0.5 | 5.032 | 5.000 | 2.195 | 2.501 | -0.707 | -0.750 | 0.023 | 0.000 | 0.007 | 0.000
|  | 0.25 | 5.027 | 5.000 | 2.139 | 2.394 | -0.712 | -0.750 | 0.022 | 0.000 | 0.007 | 0.000
|  | 0 | 5.021 | 5.000 | 2.082 | 2.286 | -0.717 | -0.750 | 0.021 | 0.000 | 0.007 | 0.000
7.5 | -0.50 | 1 | 5.009 | 5.000 | 1.557 | 1.581 | -0.493 | -0.500 | 0.001 | 0.000 | 0.000 | 0.000
|  | 0.5 | 5.004 | 5.000 | 1.448 | 1.459 | -0.496 | -0.500 | 0.001 | 0.000 | 0.000 | 0.000
|  | 0.25 | 5.002 | 5.000 | 1.391 | 1.397 | -0.498 | -0.500 | 0.000 | 0.000 | 0.000 | 0.000
|  | 0 | 5.000 | 5.000 | 1.332 | 1.333 | -0.499 | -0.500 | 0.000 | 0.000 | 0.000 | 0.000
6.25 | -0.25 | 1 | 5.005 | 5.000 | 1.262 | 1.265 | -0.248 | -0.250 | 0.000 | 0.000 | 0.000 | 0.000
|  | 0.5 | 5.002 | 5.000 | 1.166 | 1.167 | -0.249 | -0.250 | 0.000 | 0.000 | 0.000 | 0.000
|  | 0.25 | 5.001 | 5.000 | 1.117 | 1.117 | -0.250 | -0.250 | 0.000 | 0.000 | 0.000 | 0.000
|  | 0 | 5.000 | 5.000 | 1.067 | 1.067 | -0.250 | -0.250 | 0.000 | 0.000 | 0.000 | 0.000
3.75 | 0.25 | 1 | 5.009 | 5.000 | 1.263 | 1.265 | 0.249 | 0.250 | 0.000 | 0.000 | 0.000 | 0.000
|  | 0.5 | 5.003 | 5.000 | 1.166 | 1.167 | 0.249 | 0.250 | 0.000 | 0.000 | 0.000 | 0.000
|  | 0.25 | 5.002 | 5.000 | 1.117 | 1.117 | 0.250 | 0.250 | 0.000 | 0.000 | 0.000 | 0.000
|  | 0 | 5.000 | 5.000 | 1.067 | 1.067 | 0.250 | 0.250 | 0.000 | 0.000 | 0.000 | 0.000
2.5 | 0.50 | 1 | 5.019 | 5.000 | 1.567 | 1.581 | 0.497 | 0.500 | 0.000 | 0.000 | 0.000 | 0.000
|  | 0.5 | 5.008 | 5.000 | 1.453 | 1.459 | 0.499 | 0.500 | 0.000 | 0.000 | 0.000 | 0.000
|  | 0.25 | 5.004 | 5.000 | 1.394 | 1.397 | 0.499 | 0.500 | 0.000 | 0.000 | 0.000 | 0.000
|  | 0 | 5.000 | 5.000 | 1.333 | 1.333 | 0.500 | 0.500 | 0.000 | 0.000 | 0.000 | 0.000
1.25 | 0.75 | 1 | 5.091 | 5.000 | 2.606 | 2.710 | 0.744 | 0.750 | 0.000 | 0.000 | 0.000 | 0.000
|  | 0.5 | 5.042 | 5.000 | 2.453 | 2.501 | 0.747 | 0.750 | 0.000 | 0.000 | 0.000 | 0.000
|  | 0.25 | 5.020 | 5.000 | 2.372 | 2.394 | 0.748 | 0.750 | 0.000 | 0.000 | 0.000 | 0.000
|  | 0 | 5.000 | 5.000 | 2.286 | 2.286 | 0.750 | 0.750 | 0.000 | 0.000 | 0.000 | 0.000

Table B.2: Moment properties of STINARCH$(1)$ model, see Section 3.3: true
values (“tob”) vs. linear approximation (“lin”), mean close to 10.

|  |  | $\mu$ | $\sigma^{2}/\mu$ |  PACF$(1)$ |  PACF$(2)$ |  PACF$(3)$
---|---|---|---|---|---|---|---
$\alpha_{0}$ | $\alpha_{1}$ | $\delta$ | tob | lin | tob | lin | tob | lin | tob | lin | tob | lin
17.5 | -0.75 | 1 | 10.008 | 10.000 | 2.438 | 2.514 | -0.741 | -0.750 | 0.005 | 0.000 | 0.002 | 0.000
|  | 0.5 | 10.006 | 10.000 | 2.344 | 2.400 | -0.743 | -0.750 | 0.004 | 0.000 | 0.002 | 0.000
|  | 0.25 | 10.005 | 10.000 | 2.297 | 2.343 | -0.744 | -0.750 | 0.004 | 0.000 | 0.002 | 0.000
|  | 0 | 10.004 | 10.000 | 2.249 | 2.286 | -0.745 | -0.750 | 0.003 | 0.000 | 0.001 | 0.000
15 | -0.50 | 1 | 10.000 | 10.000 | 1.465 | 1.466 | -0.500 | -0.500 | 0.000 | 0.000 | 0.000 | 0.000
|  | 0.5 | 10.000 | 10.000 | 1.400 | 1.400 | -0.500 | -0.500 | 0.000 | 0.000 | 0.000 | 0.000
|  | 0.25 | 10.000 | 10.000 | 1.366 | 1.367 | -0.500 | -0.500 | 0.000 | 0.000 | 0.000 | 0.000
|  | 0 | 10.000 | 10.000 | 1.333 | 1.333 | -0.500 | -0.500 | 0.000 | 0.000 | 0.000 | 0.000
12.5 | -0.25 | 1 | 10.000 | 10.000 | 1.173 | 1.173 | -0.250 | -0.250 | 0.000 | 0.000 | 0.000 | 0.000
|  | 0.5 | 10.000 | 10.000 | 1.120 | 1.120 | -0.250 | -0.250 | 0.000 | 0.000 | 0.000 | 0.000
|  | 0.25 | 10.000 | 10.000 | 1.093 | 1.093 | -0.250 | -0.250 | 0.000 | 0.000 | 0.000 | 0.000
|  | 0 | 10.000 | 10.000 | 1.067 | 1.067 | -0.250 | -0.250 | 0.000 | 0.000 | 0.000 | 0.000
7.5 | 0.25 | 1 | 10.000 | 10.000 | 1.173 | 1.173 | 0.250 | 0.250 | 0.000 | 0.000 | 0.000 | 0.000
|  | 0.5 | 10.000 | 10.000 | 1.120 | 1.120 | 0.250 | 0.250 | 0.000 | 0.000 | 0.000 | 0.000
|  | 0.25 | 10.000 | 10.000 | 1.093 | 1.093 | 0.250 | 0.250 | 0.000 | 0.000 | 0.000 | 0.000
|  | 0 | 10.000 | 10.000 | 1.067 | 1.067 | 0.250 | 0.250 | 0.000 | 0.000 | 0.000 | 0.000
5 | 0.50 | 1 | 10.000 | 10.000 | 1.466 | 1.466 | 0.500 | 0.500 | 0.000 | 0.000 | 0.000 | 0.000
|  | 0.5 | 10.000 | 10.000 | 1.400 | 1.400 | 0.500 | 0.500 | 0.000 | 0.000 | 0.000 | 0.000
|  | 0.25 | 10.000 | 10.000 | 1.367 | 1.367 | 0.500 | 0.500 | 0.000 | 0.000 | 0.000 | 0.000
|  | 0 | 10.000 | 10.000 | 1.333 | 1.333 | 0.500 | 0.500 | 0.000 | 0.000 | 0.000 | 0.000
2.5 | 0.75 | 1 | 10.006 | 10.000 | 2.506 | 2.514 | 0.750 | 0.750 | 0.000 | 0.000 | 0.000 | 0.000
|  | 0.5 | 10.002 | 10.000 | 2.397 | 2.400 | 0.750 | 0.750 | 0.000 | 0.000 | 0.000 | 0.000
|  | 0.25 | 10.001 | 10.000 | 2.342 | 2.343 | 0.750 | 0.750 | 0.000 | 0.000 | 0.000 | 0.000
|  | 0 | 10.000 | 10.000 | 2.286 | 2.286 | 0.750 | 0.750 | 0.000 | 0.000 | 0.000 | 0.000

Table B.3: Moment properties of STINGARCH$(1,1)$ model, see Section 3.3: true
values (“tob”) vs. linear approximation (“lin”), mean close to 5.

|  |  |  | $\mu$ | $\sigma^{2}/\mu$ |  ACF$(1)$ |  ACF$(2)$ |  ACF$(3)$
---|---|---|---|---|---|---|---|---
$\alpha_{0}$ | $\alpha_{1}$ | $\beta_{1}$ | $\delta$ | tob | lin | tob | lin | tob | lin | tob | lin | tob | lin
8.5 | -0.45 | -0.25 | 1 | 5.013 | 5.000 | 1.623 | 1.656 | -0.512 | -0.521 | 0.356 | 0.365 | -0.247 | -0.255
|  |  | 0.5 | 5.006 | 5.000 | 1.511 | 1.529 | -0.516 | -0.521 | 0.360 | 0.365 | -0.251 | -0.255
|  |  | 0.25 | 5.003 | 5.000 | 1.456 | 1.464 | -0.518 | -0.521 | 0.363 | 0.365 | -0.254 | -0.255
|  |  | 0 | 5.000 | 5.000 | 1.394 | 1.397 | -0.519 | -0.521 | 0.363 | 0.365 | -0.253 | -0.255
6 | -0.45 | 0.25 | 1 | 5.006 | 5.000 | 1.428 | 1.436 | -0.404 | -0.406 | 0.080 | 0.081 | -0.015 | -0.016
|  |  | 0.5 | 5.004 | 5.000 | 1.322 | 1.325 | -0.405 | -0.406 | 0.080 | 0.081 | -0.015 | -0.016
|  |  | 0.25 | 5.001 | 5.000 | 1.262 | 1.269 | -0.405 | -0.406 | 0.081 | 0.081 | -0.016 | -0.016
|  |  | 0 | 5.000 | 5.000 | 1.212 | 1.211 | -0.407 | -0.406 | 0.080 | 0.081 | -0.014 | -0.016
4 | 0.45 | -0.25 | 1 | 5.008 | 5.000 | 1.427 | 1.436 | 0.405 | 0.406 | 0.080 | 0.081 | 0.015 | 0.016
|  |  | 0.5 | 5.003 | 5.000 | 1.322 | 1.325 | 0.407 | 0.406 | 0.082 | 0.081 | 0.016 | 0.016
|  |  | 0.25 | 4.995 | 5.000 | 1.266 | 1.269 | 0.406 | 0.406 | 0.081 | 0.081 | 0.016 | 0.016
|  |  | 0 | 4.993 | 5.000 | 1.213 | 1.211 | 0.407 | 0.406 | 0.082 | 0.081 | 0.017 | 0.016
1.5 | 0.45 | 0.25 | 1 | 5.025 | 5.000 | 1.641 | 1.656 | 0.517 | 0.521 | 0.360 | 0.365 | 0.251 | 0.255
|  |  | 0.5 | 4.999 | 5.000 | 1.522 | 1.529 | 0.520 | 0.521 | 0.364 | 0.365 | 0.255 | 0.255
|  |  | 0.25 | 4.996 | 5.000 | 1.462 | 1.464 | 0.521 | 0.521 | 0.363 | 0.365 | 0.256 | 0.255
|  |  | 0 | 5.003 | 5.000 | 1.396 | 1.397 | 0.520 | 0.521 | 0.364 | 0.365 | 0.254 | 0.255
8.5 | -0.25 | -0.45 | 1 | 5.006 | 5.000 | 1.324 | 1.331 | -0.296 | -0.299 | 0.206 | 0.209 | -0.144 | -0.147
|  |  | 0.5 | 5.004 | 5.000 | 1.228 | 1.228 | -0.300 | -0.299 | 0.211 | 0.209 | -0.147 | -0.147
|  |  | 0.25 | 5.001 | 5.000 | 1.173 | 1.176 | -0.298 | -0.299 | 0.208 | 0.209 | -0.146 | -0.147
|  |  | 0 | 4.999 | 5.000 | 1.124 | 1.123 | -0.298 | -0.299 | 0.210 | 0.209 | -0.146 | -0.147
4 | -0.25 | 0.45 | 1 | 5.005 | 5.000 | 1.261 | 1.263 | -0.223 | -0.222 | -0.043 | -0.044 | -0.010 | -0.009
|  |  | 0.5 | 5.001 | 5.000 | 1.168 | 1.166 | -0.222 | -0.222 | -0.046 | -0.044 | -0.008 | -0.009
|  |  | 0.25 | 4.999 | 5.000 | 1.112 | 1.116 | -0.221 | -0.222 | -0.045 | -0.044 | -0.010 | -0.009
|  |  | 0 | 4.999 | 5.000 | 1.064 | 1.065 | -0.222 | -0.222 | -0.044 | -0.044 | -0.008 | -0.009
6 | 0.25 | -0.45 | 1 | 5.005 | 5.000 | 1.262 | 1.263 | 0.220 | 0.222 | -0.046 | -0.044 | 0.009 | 0.009
|  |  | 0.5 | 5.008 | 5.000 | 1.165 | 1.166 | 0.222 | 0.222 | -0.045 | -0.044 | 0.009 | 0.009
|  |  | 0.25 | 5.006 | 5.000 | 1.117 | 1.116 | 0.221 | 0.222 | -0.044 | -0.044 | 0.011 | 0.009
|  |  | 0 | 4.996 | 5.000 | 1.064 | 1.065 | 0.223 | 0.222 | -0.043 | -0.044 | 0.012 | 0.009
1.5 | 0.25 | 0.45 | 1 | 5.014 | 5.000 | 1.330 | 1.331 | 0.299 | 0.299 | 0.210 | 0.209 | 0.147 | 0.147
|  |  | 0.5 | 5.008 | 5.000 | 1.226 | 1.228 | 0.298 | 0.299 | 0.208 | 0.209 | 0.146 | 0.147
|  |  | 0.25 | 5.011 | 5.000 | 1.175 | 1.176 | 0.298 | 0.299 | 0.208 | 0.209 | 0.147 | 0.147
|  |  | 0 | 4.996 | 5.000 | 1.125 | 1.123 | 0.300 | 0.299 | 0.210 | 0.209 | 0.146 | 0.147

Table B.4: Moment properties of STINGARCH$(1,1)$ model, see Section 3.3: true
values (“tob”) vs. linear approximation (“lin”), mean close to 10.

|  |  |  | $\mu$ | $\sigma^{2}/\mu$ |  ACF$(1)$ |  ACF$(2)$ |  ACF$(3)$
---|---|---|---|---|---|---|---|---
$\alpha_{0}$ | $\alpha_{1}$ | $\beta_{1}$ | $\delta$ | tob | lin | tob | lin | tob | lin | tob | lin | tob | lin
17 | -0.45 | -0.25 | 1 | 10.00 | 10.00 | 1.535 | 1.537 | -0.520 | -0.521 | 0.364 | 0.365 | -0.255 | -0.255
|  |  | 0.5 | 9.997 | 10.00 | 1.460 | 1.467 | -0.518 | -0.521 | 0.362 | 0.365 | -0.253 | -0.255
|  |  | 0.25 | 9.998 | 10.00 | 1.430 | 1.432 | -0.520 | -0.521 | 0.364 | 0.365 | -0.255 | -0.255
|  |  | 0 | 9.996 | 10.00 | 1.388 | 1.397 | -0.519 | -0.521 | 0.364 | 0.365 | -0.254 | -0.255
12 | -0.45 | 0.25 | 1 | 9.997 | 10.00 | 1.333 | 1.332 | -0.406 | -0.406 | 0.080 | 0.081 | -0.016 | -0.016
|  |  | 0.5 | 9.999 | 10.00 | 1.270 | 1.271 | -0.405 | -0.406 | 0.081 | 0.081 | -0.016 | -0.016
|  |  | 0.25 | 10.00 | 10.00 | 1.245 | 1.241 | -0.408 | -0.406 | 0.083 | 0.081 | -0.017 | -0.016
|  |  | 0 | 10.00 | 10.00 | 1.209 | 1.211 | -0.406 | -0.406 | 0.081 | 0.081 | -0.017 | -0.016
8 | 0.45 | -0.25 | 1 | 9.991 | 10.00 | 1.331 | 1.332 | 0.406 | 0.406 | 0.080 | 0.081 | 0.018 | 0.016
|  |  | 0.5 | 9.989 | 10.00 | 1.270 | 1.271 | 0.407 | 0.406 | 0.084 | 0.081 | 0.018 | 0.016
|  |  | 0.25 | 10.00 | 10.00 | 1.243 | 1.241 | 0.406 | 0.406 | 0.081 | 0.081 | 0.016 | 0.016
|  |  | 0 | 10.00 | 10.00 | 1.213 | 1.211 | 0.407 | 0.406 | 0.081 | 0.081 | 0.015 | 0.016
3 | 0.45 | 0.25 | 1 | 10.01 | 10.00 | 1.540 | 1.537 | 0.521 | 0.521 | 0.366 | 0.365 | 0.256 | 0.255
|  |  | 0.5 | 10.01 | 10.00 | 1.464 | 1.467 | 0.520 | 0.521 | 0.365 | 0.365 | 0.256 | 0.255
|  |  | 0.25 | 9.987 | 10.00 | 1.437 | 1.432 | 0.522 | 0.521 | 0.365 | 0.365 | 0.255 | 0.255
|  |  | 0 | 10.01 | 10.00 | 1.402 | 1.397 | 0.521 | 0.521 | 0.365 | 0.365 | 0.255 | 0.255
17 | -0.25 | -0.45 | 1 | 9.999 | 10.00 | 1.231 | 1.235 | -0.298 | -0.299 | 0.208 | 0.209 | -0.146 | -0.147
|  |  | 0.5 | 9.999 | 10.00 | 1.178 | 1.179 | -0.299 | -0.299 | 0.211 | 0.209 | -0.147 | -0.147
|  |  | 0.25 | 9.996 | 10.00 | 1.150 | 1.151 | -0.299 | -0.299 | 0.210 | 0.209 | -0.147 | -0.147
|  |  | 0 | 10.00 | 10.00 | 1.124 | 1.123 | -0.302 | -0.299 | 0.210 | 0.209 | -0.148 | -0.147
8 | -0.25 | 0.45 | 1 | 9.999 | 10.00 | 1.171 | 1.171 | -0.221 | -0.222 | -0.045 | -0.044 | -0.009 | -0.009
|  |  | 0.5 | 10.00 | 10.00 | 1.119 | 1.118 | -0.222 | -0.222 | -0.046 | -0.044 | -0.008 | -0.009
|  |  | 0.25 | 9.997 | 10.00 | 1.092 | 1.092 | -0.221 | -0.222 | -0.047 | -0.044 | -0.008 | -0.009
|  |  | 0 | 9.998 | 10.00 | 1.065 | 1.065 | -0.222 | -0.222 | -0.046 | -0.044 | -0.010 | -0.009
12 | 0.25 | -0.45 | 1 | 10.00 | 10.00 | 1.172 | 1.171 | 0.221 | 0.222 | -0.045 | -0.044 | 0.009 | 0.009
|  |  | 0.5 | 10.01 | 10.00 | 1.115 | 1.118 | 0.221 | 0.222 | -0.045 | -0.044 | 0.010 | 0.009
|  |  | 0.25 | 9.996 | 10.00 | 1.092 | 1.092 | 0.222 | 0.222 | -0.045 | -0.044 | 0.010 | 0.009
|  |  | 0 | 10.00 | 10.00 | 1.066 | 1.065 | 0.222 | 0.222 | -0.046 | -0.044 | 0.010 | 0.009
3 | 0.25 | 0.45 | 1 | 9.997 | 10.00 | 1.236 | 1.235 | 0.300 | 0.299 | 0.208 | 0.209 | 0.146 | 0.147
|  |  | 0.5 | 10.01 | 10.00 | 1.179 | 1.179 | 0.299 | 0.299 | 0.210 | 0.209 | 0.147 | 0.147
|  |  | 0.25 | 10.00 | 10.00 | 1.150 | 1.151 | 0.300 | 0.299 | 0.211 | 0.209 | 0.147 | 0.147
|  |  | 0 | 10.00 | 10.00 | 1.122 | 1.123 | 0.300 | 0.299 | 0.210 | 0.209 | 0.147 | 0.147

Table B.5: Estimation of STINARCH$(1)$ model, see Section 4.3: mean and
standard error (SE) of simulated estimates, where MLE use $\delta:=0.25$
throughout.

| Mean: | Standard error:
---|---|---
| MLE | CLADE | CLS | MLE | CLADE | CLS
$n$ | $\hat{\alpha}_{0}$ | $\hat{\alpha}_{1}$ | $\check{\alpha}_{0}$ | $\check{\alpha}_{1}$ | $\breve{\alpha}_{0}$ | $\breve{\alpha}_{1}$ | $\hat{\alpha}_{0}$ | $\hat{\alpha}_{1}$ | $\check{\alpha}_{0}$ | $\check{\alpha}_{1}$ | $\breve{\alpha}_{0}$ | $\breve{\alpha}_{1}$
$(\alpha_{0},\alpha_{1})=(7.5,-0.5)$ with true $\delta=0$
100 | 7.470 | -0.497 | 7.336 | -0.500 | 7.468 | -0.496 | 0.498 | 0.077 | 0.488 | 0.073 | 0.518 | 0.083
250 | 7.505 | -0.503 | 7.354 | -0.499 | 7.497 | -0.500 | 0.303 | 0.046 | 0.294 | 0.037 | 0.326 | 0.052
500 | 7.511 | -0.503 | 7.392 | -0.500 | 7.493 | -0.498 | 0.217 | 0.033 | 0.228 | 0.023 | 0.234 | 0.037
1000 | 7.510 | -0.503 | 7.436 | -0.500 | 7.495 | -0.499 | 0.149 | 0.022 | 0.172 | 0.013 | 0.162 | 0.026
$(\alpha_{0},\alpha_{1})=(7.5,-0.5)$ with true $\delta=0.25$
100 | 7.483 | -0.496 | 7.343 | -0.497 | 7.474 | -0.494 | 0.487 | 0.077 | 0.513 | 0.076 | 0.517 | 0.084
250 | 7.496 | -0.500 | 7.372 | -0.499 | 7.490 | -0.498 | 0.301 | 0.047 | 0.282 | 0.037 | 0.320 | 0.052
500 | 7.502 | -0.500 | 7.394 | -0.500 | 7.493 | -0.498 | 0.220 | 0.032 | 0.227 | 0.023 | 0.237 | 0.037
1000 | 7.499 | -0.500 | 7.427 | -0.499 | 7.493 | -0.498 | 0.156 | 0.023 | 0.183 | 0.015 | 0.169 | 0.027
$(\alpha_{0},\alpha_{1})=(7.5,-0.5)$ with true $\delta=1$
100 | 7.469 | -0.491 | 7.360 | -0.498 | 7.482 | -0.494 | 0.523 | 0.077 | 0.563 | 0.085 | 0.559 | 0.088
250 | 7.455 | -0.489 | 7.389 | -0.499 | 7.465 | -0.492 | 0.317 | 0.048 | 0.302 | 0.040 | 0.343 | 0.055
500 | 7.472 | -0.491 | 7.406 | -0.499 | 7.480 | -0.493 | 0.226 | 0.034 | 0.216 | 0.024 | 0.243 | 0.038
1000 | 7.463 | -0.490 | 7.451 | -0.500 | 7.479 | -0.494 | 0.152 | 0.023 | 0.155 | 0.014 | 0.163 | 0.025
$(\alpha_{0},\alpha_{1})=(2.5,0.5)$ with true $\delta=0$
100 | 2.601 | 0.478 | 2.460 | 0.477 | 2.627 | 0.473 | 0.479 | 0.095 | 0.459 | 0.091 | 0.485 | 0.097
250 | 2.539 | 0.490 | 2.433 | 0.484 | 2.555 | 0.488 | 0.293 | 0.058 | 0.288 | 0.056 | 0.302 | 0.060
500 | 2.509 | 0.498 | 2.450 | 0.491 | 2.524 | 0.495 | 0.194 | 0.039 | 0.182 | 0.034 | 0.201 | 0.040
1000 | 2.500 | 0.499 | 2.467 | 0.494 | 2.511 | 0.498 | 0.145 | 0.029 | 0.135 | 0.022 | 0.150 | 0.030
$(\alpha_{0},\alpha_{1})=(2.5,0.5)$ with true $\delta=0.25$
100 | 2.613 | 0.474 | 2.467 | 0.474 | 2.636 | 0.470 | 0.476 | 0.093 | 0.467 | 0.096 | 0.480 | 0.095
250 | 2.546 | 0.490 | 2.444 | 0.487 | 2.551 | 0.490 | 0.292 | 0.057 | 0.260 | 0.051 | 0.306 | 0.061
500 | 2.523 | 0.495 | 2.454 | 0.490 | 2.533 | 0.493 | 0.205 | 0.040 | 0.178 | 0.034 | 0.211 | 0.042
1000 | 2.517 | 0.497 | 2.471 | 0.494 | 2.523 | 0.496 | 0.145 | 0.028 | 0.139 | 0.023 | 0.152 | 0.030
$(\alpha_{0},\alpha_{1})=(2.5,0.5)$ with true $\delta=1$
100 | 2.644 | 0.471 | 2.499 | 0.472 | 2.660 | 0.467 | 0.483 | 0.094 | 0.489 | 0.098 | 0.495 | 0.097
250 | 2.568 | 0.489 | 2.454 | 0.487 | 2.573 | 0.488 | 0.298 | 0.057 | 0.294 | 0.055 | 0.311 | 0.060
500 | 2.549 | 0.493 | 2.465 | 0.491 | 2.550 | 0.492 | 0.205 | 0.041 | 0.184 | 0.036 | 0.215 | 0.043
1000 | 2.546 | 0.494 | 2.475 | 0.494 | 2.544 | 0.494 | 0.144 | 0.028 | 0.129 | 0.023 | 0.149 | 0.030

Table B.6: Estimation of STINARCH$(1)$ model, see Section 4.3: standard error
(SE) of simulated MLE and mean of approximate SE, where MLE use $\delta:=0.25$
throughout.

|  True $\delta=0$: |  True $\delta=0.25$: |  True $\delta=1$:
---|---|---|---
| simul. SE | approx. SE | simul. SE | approx. SE | simul. SE | approx. SE
$n$ | $\hat{\alpha}_{0}$ | $\hat{\alpha}_{1}$ | $\hat{\alpha}_{0}$ | $\hat{\alpha}_{1}$ | $\hat{\alpha}_{0}$ | $\hat{\alpha}_{1}$ | $\hat{\alpha}_{0}$ | $\hat{\alpha}_{1}$ | $\hat{\alpha}_{0}$ | $\hat{\alpha}_{1}$ | $\hat{\alpha}_{0}$ | $\hat{\alpha}_{1}$
$(\alpha_{0},\alpha_{1})=(7.5,-0.5)$
100 | 0.498 | 0.077 | 0.499 | 0.078 | 0.487 | 0.077 | 0.522 | 0.081 | 0.523 | 0.077 | 0.530 | 0.082
250 | 0.303 | 0.046 | 0.306 | 0.047 | 0.301 | 0.047 | 0.316 | 0.048 | 0.317 | 0.048 | 0.324 | 0.048
500 | 0.217 | 0.033 | 0.215 | 0.032 | 0.220 | 0.032 | 0.211 | 0.032 | 0.226 | 0.034 | 0.220 | 0.032
1000 | 0.149 | 0.022 | 0.147 | 0.022 | 0.156 | 0.023 | 0.155 | 0.023 | 0.152 | 0.023 | 0.156 | 0.023
$(\alpha_{0},\alpha_{1})=(2.5,0.5)$
100 | 0.479 | 0.095 | 0.473 | 0.094 | 0.476 | 0.093 | 0.474 | 0.095 | 0.483 | 0.094 | 0.488 | 0.094
250 | 0.293 | 0.058 | 0.288 | 0.058 | 0.292 | 0.057 | 0.284 | 0.056 | 0.298 | 0.057 | 0.289 | 0.057
500 | 0.194 | 0.039 | 0.202 | 0.041 | 0.205 | 0.040 | 0.198 | 0.039 | 0.205 | 0.041 | 0.207 | 0.041
1000 | 0.145 | 0.029 | 0.141 | 0.027 | 0.145 | 0.028 | 0.140 | 0.028 | 0.144 | 0.028 | 0.141 | 0.028

Table B.7: Estimation of STINARCH$(1)$ model, see Section 4.3, where also
$\delta$ is estimated: mean and standard error (SE) of simulated MLE as well
as mean of approximate SE.

|  |  |  | Mean: | Standard error: | Mean approx. SE:
---|---|---|---|---|---|---
$\alpha_{0}$ | $\alpha_{1}$ | $\delta$ | $n$ | $\hat{\alpha}_{0}$ | $\hat{\alpha}_{1}$ | $\hat{\delta}$ | $\hat{\alpha}_{0}$ | $\hat{\alpha}_{1}$ | $\hat{\delta}$ | $\hat{\alpha}_{0}$ | $\hat{\alpha}_{1}$ | $\hat{\delta}$
7.5 | -0.5 | 0.25 | 100 | 7.495 | -0.499 | 0.343 | 0.498 | 0.076 | 0.449 | 0.497 | 0.077 | 0.445
|  |  | 250 | 7.499 | -0.499 | 0.273 | 0.302 | 0.046 | 0.295 | 0.295 | 0.045 | 0.297
|  |  | 500 | 7.498 | -0.500 | 0.258 | 0.229 | 0.034 | 0.226 | 0.223 | 0.034 | 0.202
|  |  | 1000 | 7.508 | -0.501 | 0.254 | 0.153 | 0.023 | 0.172 | 0.156 | 0.023 | 0.164
7.5 | -0.5 | 1 | 100 | 7.495 | -0.496 | 0.945 | 0.517 | 0.080 | 0.745 | 0.534 | 0.083 | 0.677
|  |  | 250 | 7.477 | -0.496 | 0.947 | 0.329 | 0.049 | 0.469 | 0.323 | 0.050 | 0.460
|  |  | 500 | 7.501 | -0.501 | 0.982 | 0.224 | 0.035 | 0.335 | 0.224 | 0.034 | 0.340
|  |  | 1000 | 7.494 | -0.499 | 0.987 | 0.162 | 0.024 | 0.242 | 0.159 | 0.025 | 0.251
2.5 | 0.5 | 0.25 | 100 | 2.630 | 0.473 | 0.312 | 0.474 | 0.091 | 0.426 | 0.446 | 0.090 | 0.431
|  |  | 250 | 2.554 | 0.488 | 0.284 | 0.288 | 0.058 | 0.286 | 0.282 | 0.057 | 0.283
|  |  | 500 | 2.526 | 0.493 | 0.264 | 0.203 | 0.040 | 0.221 | 0.199 | 0.039 | 0.217
|  |  | 1000 | 2.509 | 0.498 | 0.261 | 0.141 | 0.028 | 0.166 | 0.136 | 0.028 | 0.159
2.5 | 0.5 | 1 | 100 | 2.614 | 0.475 | 0.960 | 0.473 | 0.089 | 0.736 | 0.480 | 0.093 | 0.690
|  |  | 250 | 2.530 | 0.492 | 0.967 | 0.312 | 0.060 | 0.486 | 0.298 | 0.057 | 0.459
|  |  | 500 | 2.521 | 0.496 | 0.982 | 0.209 | 0.042 | 0.329 | 0.207 | 0.041 | 0.333
|  |  | 1000 | 2.509 | 0.497 | 0.995 | 0.144 | 0.029 | 0.240 | 0.147 | 0.028 | 0.241

Table B.8: Estimation of STINGARCH$(1,1)$ model, see Section 4.3: mean and
standard error (SE) of simulated estimates, where MLE use $\delta:=0.25$
throughout.

| Mean: | Standard error:
---|---|---
| MLE | CLADE | MLE | CLADE
$n$ | $\hat{\alpha}_{0}$ | $\hat{\alpha}_{1}$ | $\hat{\beta}_{1}$ | $\check{\alpha}_{0}$ | $\check{\alpha}_{1}$ | $\check{\beta}_{1}$ | $\hat{\alpha}_{0}$ | $\hat{\alpha}_{1}$ | $\hat{\beta}_{1}$ | $\check{\alpha}_{0}$ | $\check{\alpha}_{1}$ | $\check{\beta}_{1}$
$(\alpha_{0},\alpha_{1},\beta_{1})=(7.25,-0.45,0.25)$ with true $\delta=0$
250 | 7.152 | -0.451 | 0.267 | 7.041 | -0.450 | 0.264 | 0.937 | 0.059 | 0.153 | 1.032 | 0.076 | 0.173
500 | 7.192 | -0.451 | 0.259 | 7.091 | -0.451 | 0.256 | 0.638 | 0.040 | 0.105 | 0.788 | 0.051 | 0.134
1000 | 7.214 | -0.452 | 0.258 | 7.104 | -0.451 | 0.255 | 0.431 | 0.031 | 0.071 | 0.565 | 0.038 | 0.091
$(\alpha_{0},\alpha_{1},\beta_{1})=(7.25,-0.45,0.25)$ with true $\delta=0.25$
250 | 7.129 | -0.452 | 0.271 | 7.048 | -0.453 | 0.267 | 0.911 | 0.060 | 0.153 | 0.976 | 0.076 | 0.166
500 | 7.175 | -0.453 | 0.266 | 7.098 | -0.453 | 0.259 | 0.610 | 0.041 | 0.100 | 0.768 | 0.052 | 0.128
1000 | 7.243 | -0.452 | 0.253 | 7.136 | -0.452 | 0.251 | 0.439 | 0.028 | 0.073 | 0.563 | 0.036 | 0.093
$(\alpha_{0},\alpha_{1},\beta_{1})=(7.25,-0.45,0.25)$ with true $\delta=1$
250 | 7.113 | -0.449 | 0.272 | 7.064 | -0.451 | 0.264 | 0.930 | 0.060 | 0.154 | 1.000 | 0.074 | 0.170
500 | 7.158 | -0.451 | 0.266 | 7.100 | -0.452 | 0.259 | 0.625 | 0.041 | 0.102 | 0.743 | 0.054 | 0.125
1000 | 7.253 | -0.448 | 0.248 | 7.170 | -0.450 | 0.245 | 0.456 | 0.028 | 0.076 | 0.573 | 0.037 | 0.097
$(\alpha_{0},\alpha_{1},\beta_{1})=(2.75,0.45,-0.25)$ with true $\delta=0$
250 | 2.764 | 0.450 | -0.257 | 2.549 | 0.446 | -0.248 | 0.531 | 0.068 | 0.151 | 0.637 | 0.089 | 0.183
500 | 2.758 | 0.451 | -0.258 | 2.510 | 0.448 | -0.241 | 0.368 | 0.045 | 0.101 | 0.514 | 0.056 | 0.145
1000 | 2.762 | 0.450 | -0.256 | 2.532 | 0.449 | -0.248 | 0.246 | 0.033 | 0.069 | 0.372 | 0.044 | 0.105
$(\alpha_{0},\alpha_{1},\beta_{1})=(2.75,0.45,-0.25)$ with true $\delta=0.25$
250 | 2.801 | 0.448 | -0.264 | 2.598 | 0.445 | -0.257 | 0.525 | 0.066 | 0.146 | 0.626 | 0.080 | 0.178
500 | 2.806 | 0.445 | -0.262 | 2.568 | 0.444 | -0.250 | 0.371 | 0.046 | 0.100 | 0.491 | 0.059 | 0.135
1000 | 2.765 | 0.448 | -0.252 | 2.537 | 0.446 | -0.243 | 0.264 | 0.031 | 0.072 | 0.388 | 0.040 | 0.107
$(\alpha_{0},\alpha_{1},\beta_{1})=(2.75,0.45,-0.25)$ with true $\delta=1$
250 | 2.854 | 0.443 | -0.265 | 2.639 | 0.443 | -0.258 | 0.538 | 0.065 | 0.148 | 0.656 | 0.080 | 0.186
500 | 2.840 | 0.438 | -0.257 | 2.599 | 0.441 | -0.248 | 0.378 | 0.046 | 0.103 | 0.478 | 0.060 | 0.137
1000 | 2.803 | 0.442 | -0.248 | 2.577 | 0.442 | -0.241 | 0.271 | 0.031 | 0.072 | 0.339 | 0.039 | 0.098

Table B.9: Estimation of STINGARCH$(1,1)$ model, see Section 4.3: mean and
standard error (SE) of simulated estimates, where MLE use $\delta:=0.25$
throughout.

| Means: | Standard errors:
---|---|---
| MLE | CLADE | MLE | CLADE
$n$ | $\hat{\alpha}_{0}$ | $\hat{\alpha}_{1}$ | $\hat{\beta}_{1}$ | $\check{\alpha}_{0}$ | $\check{\alpha}_{1}$ | $\check{\beta}_{1}$ | $\hat{\alpha}_{0}$ | $\hat{\alpha}_{1}$ | $\hat{\beta}_{1}$ | $\check{\alpha}_{0}$ | $\check{\alpha}_{1}$ | $\check{\beta}_{1}$
$(\alpha_{0},\alpha_{1},\beta_{1})=(6.25,-0.25,0.45)$ with true $\delta=0$
250 | 6.460 | -0.260 | 0.434 | 6.358 | -0.260 | 0.435 | 1.506 | 0.060 | 0.185 | 1.220 | 0.075 | 0.159
500 | 6.335 | -0.253 | 0.442 | 6.285 | -0.257 | 0.440 | 1.161 | 0.042 | 0.141 | 1.087 | 0.050 | 0.138
1000 | 6.278 | -0.254 | 0.450 | 6.214 | -0.253 | 0.446 | 0.838 | 0.031 | 0.098 | 0.901 | 0.039 | 0.108
$(\alpha_{0},\alpha_{1},\beta_{1})=(6.25,-0.25,0.45)$ with true $\delta=0.25$
250 | 6.368 | -0.263 | 0.448 | 6.278 | -0.262 | 0.447 | 1.422 | 0.062 | 0.173 | 1.149 | 0.077 | 0.147
500 | 6.366 | -0.256 | 0.442 | 6.306 | -0.256 | 0.437 | 1.147 | 0.044 | 0.137 | 1.084 | 0.054 | 0.135
1000 | 6.310 | -0.252 | 0.444 | 6.207 | -0.254 | 0.447 | 0.848 | 0.031 | 0.102 | 0.942 | 0.038 | 0.116
$(\alpha_{0},\alpha_{1},\beta_{1})=(6.25,-0.25,0.45)$ with true $\delta=1$
250 | 6.417 | -0.261 | 0.440 | 6.298 | -0.258 | 0.442 | 1.451 | 0.062 | 0.181 | 1.177 | 0.078 | 0.151
500 | 6.391 | -0.256 | 0.439 | 6.354 | -0.254 | 0.429 | 1.154 | 0.044 | 0.142 | 1.149 | 0.056 | 0.144
1000 | 6.326 | -0.252 | 0.443 | 6.261 | -0.253 | 0.441 | 0.872 | 0.030 | 0.103 | 0.950 | 0.038 | 0.115
$(\alpha_{0},\alpha_{1},\beta_{1})=(3.75,0.25,-0.45)$ with true $\delta=0$
250 | 3.683 | 0.255 | -0.436 | 3.543 | 0.209 | -0.414 | 0.712 | 0.065 | 0.210 | 0.669 | 0.138 | 0.200
500 | 3.725 | 0.253 | -0.451 | 3.555 | 0.209 | -0.419 | 0.516 | 0.043 | 0.151 | 0.590 | 0.115 | 0.181
1000 | 3.744 | 0.251 | -0.453 | 3.608 | 0.204 | -0.427 | 0.347 | 0.032 | 0.100 | 0.442 | 0.109 | 0.141
$(\alpha_{0},\alpha_{1},\beta_{1})=(3.75,0.25,-0.45)$ with true $\delta=0.25$
250 | 3.664 | 0.256 | -0.428 | 3.504 | 0.217 | -0.404 | 0.708 | 0.065 | 0.209 | 0.702 | 0.131 | 0.211
500 | 3.771 | 0.248 | -0.455 | 3.594 | 0.213 | -0.431 | 0.501 | 0.044 | 0.144 | 0.576 | 0.108 | 0.172
1000 | 3.758 | 0.249 | -0.451 | 3.617 | 0.214 | -0.439 | 0.369 | 0.031 | 0.104 | 0.457 | 0.095 | 0.137
$(\alpha_{0},\alpha_{1},\beta_{1})=(3.75,0.25,-0.45)$ with true $\delta=1$
250 | 3.733 | 0.252 | -0.434 | 3.537 | 0.229 | -0.417 | 0.705 | 0.062 | 0.204 | 0.731 | 0.119 | 0.222
500 | 3.821 | 0.244 | -0.454 | 3.600 | 0.221 | -0.430 | 0.531 | 0.045 | 0.151 | 0.615 | 0.102 | 0.186
1000 | 3.773 | 0.246 | -0.441 | 3.576 | 0.232 | -0.433 | 0.385 | 0.033 | 0.109 | 0.496 | 0.079 | 0.144

Table B.10: Estimation of STINGARCH$(1,1)$ model, see Section 4.3: standard
error of simulated MLE and mean of approximate SE, where MLE use
$\delta:=0.25$ throughout.

| Simulated SE: | Approximate SE: | Simulated SE: | Approximate SE:
---|---|---|---|---
$n$ | $\hat{\alpha}_{0}$ | $\hat{\alpha}_{1}$ | $\hat{\beta}_{1}$ | $\hat{\alpha}_{0}$ | $\hat{\alpha}_{1}$ | $\hat{\beta}_{1}$ | $\hat{\alpha}_{0}$ | $\hat{\alpha}_{1}$ | $\hat{\beta}_{1}$ | $\hat{\alpha}_{0}$ | $\hat{\alpha}_{1}$ | $\hat{\beta}_{1}$
$(\alpha_{0},\alpha_{1},\beta_{1})=(7.25,-0.45,0.25)$ with true $\delta=0$ |  $(6.25,-0.25,0.45)$ with true $\delta=0$
250 | 0.937 | 0.059 | 0.153 | 0.942 | 0.061 | 0.152 | 1.506 | 0.060 | 0.185 | 1.425 | 0.061 | 0.171
500 | 0.638 | 0.040 | 0.105 | 0.633 | 0.039 | 0.103 | 1.161 | 0.042 | 0.141 | 1.154 | 0.042 | 0.140
1000 | 0.431 | 0.031 | 0.071 | 0.438 | 0.028 | 0.074 | 0.838 | 0.031 | 0.098 | 0.849 | 0.031 | 0.103
$(\alpha_{0},\alpha_{1},\beta_{1})=(7.25,-0.45,0.25)$ with true $\delta=0.25$ |  $(6.25,-0.25,0.45)$ with true $\delta=0.25$
250 | 0.911 | 0.060 | 0.153 | 0.905 | 0.059 | 0.149 | 1.422 | 0.062 | 0.173 | 1.445 | 0.061 | 0.177
500 | 0.610 | 0.041 | 0.100 | 0.629 | 0.040 | 0.104 | 1.147 | 0.044 | 0.137 | 1.132 | 0.043 | 0.135
1000 | 0.439 | 0.028 | 0.073 | 0.446 | 0.029 | 0.073 | 0.848 | 0.031 | 0.102 | 0.880 | 0.031 | 0.104
$(\alpha_{0},\alpha_{1},\beta_{1})=(7.25,-0.45,0.25)$ with true $\delta=1$ |  $(6.25,-0.25,0.45)$ with true $\delta=1$
250 | 0.930 | 0.060 | 0.154 | 0.931 | 0.062 | 0.151 | 1.451 | 0.062 | 0.181 | 1.455 | 0.060 | 0.178
500 | 0.625 | 0.041 | 0.102 | 0.607 | 0.042 | 0.101 | 1.154 | 0.044 | 0.142 | 1.152 | 0.043 | 0.137
1000 | 0.456 | 0.028 | 0.076 | 0.435 | 0.030 | 0.072 | 0.872 | 0.030 | 0.103 | 0.869 | 0.032 | 0.102
$(\alpha_{0},\alpha_{1},\beta_{1})=(2.75,0.45,-0.25)$ with true $\delta=0$ |  $(3.75,0.25,-0.45)$ with true $\delta=0$
250 | 0.531 | 0.068 | 0.151 | 0.529 | 0.066 | 0.149 | 0.712 | 0.065 | 0.210 | 0.715 | 0.065 | 0.213
500 | 0.368 | 0.045 | 0.101 | 0.356 | 0.045 | 0.100 | 0.516 | 0.043 | 0.151 | 0.500 | 0.044 | 0.146
1000 | 0.246 | 0.033 | 0.069 | 0.256 | 0.031 | 0.070 | 0.347 | 0.032 | 0.100 | 0.368 | 0.031 | 0.105
$(\alpha_{0},\alpha_{1},\beta_{1})=(2.75,0.45,-0.25)$ with true $\delta=0.25$ |  $(3.75,0.25,-0.45)$ with true $\delta=0.25$
250 | 0.525 | 0.066 | 0.146 | 0.523 | 0.066 | 0.148 | 0.708 | 0.065 | 0.209 | 0.708 | 0.063 | 0.212
500 | 0.371 | 0.046 | 0.100 | 0.362 | 0.047 | 0.102 | 0.501 | 0.044 | 0.144 | 0.513 | 0.046 | 0.148
1000 | 0.264 | 0.031 | 0.072 | 0.253 | 0.032 | 0.069 | 0.369 | 0.031 | 0.104 | 0.373 | 0.032 | 0.107
$(\alpha_{0},\alpha_{1},\beta_{1})=(2.75,0.45,-0.25)$ with true $\delta=1$ |  $(3.75,0.25,-0.45)$ with true $\delta=1$
250 | 0.538 | 0.065 | 0.148 | 0.550 | 0.065 | 0.148 | 0.705 | 0.062 | 0.204 | 0.775 | 0.065 | 0.234
500 | 0.378 | 0.046 | 0.103 | 0.375 | 0.045 | 0.102 | 0.531 | 0.045 | 0.151 | 0.525 | 0.046 | 0.150
1000 | 0.271 | 0.031 | 0.072 | 0.261 | 0.033 | 0.070 | 0.385 | 0.033 | 0.109 | 0.374 | 0.032 | 0.106