Laplacian Split-BREAK Process with Application in Dynamic Analysis of the World Oil and Gas Market

Abstract

This manuscript deals with a novel, nonlinear, and non-stationary stochastic model with symmetric, Laplacian distributed innovations. The obtained model, named Laplacian Split-BREAK (LSB) process, is intended for dynamic analysis of time series with pronounced and permanent fluctuations. By using the method of characteristic functions (CFs), the basic stochastic properties of the LSB process are proven, with a special emphasis on its asymptotic behaviour. The different procedures for estimating its parameters are also given, along with numerical simulations of the obtained estimators. Finally, it has been shown that the LSB process, as an adequate stochastic model, can be applied in the analysis of dynamics in the world market of crude oil and natural gas.

1. Introduction

One of the important topics in modern research is the stochastic modelling of time series with emphasized and persistent fluctuations. To this end, in recent decades, various stochastic models have been proposed, especially in econometrics and financial engineering (see, c.f. [1,2,3,4,5,6,7,8]). A particular problem arises when the observed time series have a non-stationary dynamic, which usually affects the greater complexity of their stochastic structure (for more recent studies, see, c.f. [9,10,11,12,13,14]). To solve these and similar problems, Engle and Smith [15] proposed the stochastic permanent breaking (STOPBREAK) process, which was later examined by many authors, especially in the domain of structural and permanent changes in the fluctuations of some real-world data [16,17,18,19,20]. Stojanović et al. [21] introduced a more general form for the STOPBREAK-based stochastic process, named the Split-BREAK process, which was also applied in modelling various time series with permanent and pronounced fluctuations. It is worth pointing out that some more general forms of the Split-BREAK process, the so-called General (or Gaussian) Split-BREAK (GSB) process, were introduced later and discussed in Stojanović et al. [22,23,24] as well as Jovanović et al. [25].

The main characteristics of the GSB process are based on its Gaussian innovations, which somewhat simplifies the examination of the properties of this model. Nevertheless, there are some specific time series that do not have the Gaussian property. To that end, as an alternative approach, some other stochastic distributions can be considered. They are usually assumed to have a symmetric distribution (as in the case of a zero-mean Gaussian distribution), but with some additional properties that certain non-Gaussian distributions have. For these reasons, a new, general form of the Split-BREAK process with symmetric Laplacian innovations, called the LSB process, is proposed here.

In the next section, the definition and main properties of the LSB process are outlined. Some additional stochastic characteristics of the LSB process, with special emphasis on examining its asymptotic properties, are discussed in Section 3. As a basic tool, the characteristic functions (CFs) method was used here. The CF calculation procedure is known for simplifying the finding of probability distributions, rather than calculating them directly, and is particularly useful in determining the cumulative distribution functions (CDFs). As will be seen below, it will be quite suitable for investigating time series distributions as constituents of the LSB process. Various procedures for parameter estimation of the LSB process, as well as examining the asymptotic properties of the obtained estimators, are given in Section 4. Section 5 describes the Monte Carlo simulations of the obtained estimators, while Section 6 presents the application of the LSB process in the dynamic analysis of prices and trading volumes of crude oil and natural gas in the world market. Lastly, some concluding remarks are highlighted in Section 7.

2. Definition and Structure of the LSB Process

Akin to previous works of GSB process [21,22,23,24,25], the underlying LSB time series is given by the following additive decomposition:

(1)$y_t=m_t+{\epsilon }_t,$

where $t=0,1,...,T$ is the finite time values set. In Equation (1), series $\left(m_t\right)$ represents the martingale means and $\left({\epsilon }_t\right)$ are the innovations, that is, the series of independent identical distributed (IID) random variables (RVs). It is usually assumed that innovations $\left({\epsilon }_t\right)$ are given on the probability space $\left(\mathsf{\Omega },\mathcal{F},P\right),$ which is expanded by some filtration $F=\left({\mathcal{F}}_t\right).$ Practically interpreted, filtration $\left({\mathcal{F}}_t\right)$ collects the ‘information’ about some financial index at a certain point in time $t$. Therefore, the RVs $\left({\epsilon }_t\right)$ will be ${\mathcal{F}}_t$-adaptive for each $t=0,1,...,T$. In the following, we assume that innovations $\left({\epsilon }_t\right)$ have the centred Laplace distribution, given with probability density function (PDF):

(2)$f_{\epsilon }\left(x\right)=\frac{1}{2\lambda }\exp \left(-\frac{\left|x\right|}{\lambda }\right),x\in \mathbb{R},$

where $\lambda >0$ is the scale parameter. Let us emphasize that the proposal of the Laplace distribution is motivated by the fact that this distribution, compared to the Gaussian, can more adequately fit the distributions of many financial time series with pronounced “peaks”. Such a situation is shown in Figure 1, which shows the histograms of the empirical distributions of the log-returns $h_t=\ln \left(P_t/P_{t-1}\right),t=1,\dots ,T.$ Here, the price of crude oil and natural gas on the world market is taken as the basic financial index $\left(P_t\right)$, which will be discussed in more detail below. According to the well-known properties of Laplace distribution (see, c.f. [26], pp. 60), the conditional mean and variance of innovations $\left({\epsilon }_t\right)$ are, respectively,

(3)$E\left({\epsilon }_t|{\mathcal{F}}_{t-1}\right)=0,\hspace{1em}\hspace{1em}\mathrm{Var}\left({\epsilon }_t|{\mathcal{F}}_{t-1}\right)=E\left({\epsilon }_t^2|{\mathcal{F}}_{t-1}\right)=2{\lambda }^2.$

After integrating the function $f_{\epsilon }\left(x\right)$, where, due to the absolute value, two symmetric cases are distinguished, its CDF is obtained as follows:

$F_{\epsilon }\left(x\right)≔P\left\{{\epsilon }_t<x\right\}={\int }_{-\mathrm{\infty }}^x{f}_{\epsilon }\left(z\right)dz=\left\{\begin{array}{ll}\frac{1}{2}\exp \left(x/\lambda \right),& x<0\\ 1-\frac{1}{2}\exp \left(-x/\lambda \right),& x\ge 0\end{array}\right.$

(4)$=\frac{1}{2}+\frac{1}{2}\mathrm{s}\mathrm{g}\mathrm{n}\left(x\right)\left(1-\exp \left(-\frac{\left|x\right|}{\lambda }\right)\right).$

On the other hand, for a series of martingale means $\left(m_t\right),$ it is assumed that they are given recurrently by the relation:

(5)$m_t=m_{t-1}+q_{t-1}{\epsilon }_{t-1}=m_0+\sum _{j=0}^{t-1}{q}_j{\epsilon }_j,$

where it is almost certainly (as) $m_0\stackrel{as}{=}\mu \left(const.\right)$ and ${\epsilon }_{-1}={\epsilon }_0\stackrel{as}{=}0$. Moreover, ($q_t$) is a series of noise indicators, that is, the RVs that depend on innovation series $\left({\epsilon }_t\right)$ as follows:

(6)$q_t=I\left({\epsilon }_{t-1}^2>c\right)=\left\{\begin{array}{ll}1,& {\epsilon }_{t-1}^2>c\\ 0,& {\epsilon }_{t-1}^2\le c.\end{array}\right.$

Here, $c>0$ is the so-called critical value of the reaction, i.e., the parameter that indicates the significance of previous realizations of innovations $\left({\epsilon }_t\right)$, which would allow that their current values to be included in Equation (6). More precisely, in the case $q_{t-1}=0$, the value of martingale mean $m_t$ is equal to its previous value $m_{t-1}$. Therefore, the value of the basic LSB series ($y_t$), given by Equation (1), is obtained with a “small” fluctuation, dependent on (${\epsilon }_t$) only. On the contrary, when $q_t=1$, a pronounced fluctuation of the value $y_t$ is registered. In this way, the previous realizations of innovations $\left({\epsilon }_t\right)$ impact the level of variations of the basic LSB time series ($y_t$), that is, the intensity of fluctuations of the LSB process. Using some simple calculations, in the same way as in Jovanović et al. [25], the following properties of the mentioned time series can be easily proven:

Theorem 1. 

Let $\left(y_t\right)$ and $\left(m_t\right)$ be the series defined by Equations (1) and (5), respectively. Then, both series have the constant and equal mean:

$E\left(y_t\right)=E\left(m_t\right)=\mu ,$

as well as the variances:

(7) $\mathrm{V}\mathrm{a}\mathrm{r}\left(y_t\right)=2\left(t{a}_c+1\right){\lambda }^2, \hspace{1em}\mathrm{V}\mathrm{a}\mathrm{r}\left(m_t\right)=2t{a}_c{\lambda }^2, \hspace{1em}t\ge 0,$

where $a_c=E\left(q_t\right)=E\left(q_t^2\right)=P\left\{{\epsilon }_t^2>c\right\}.$ In addition, the correlation functions of the series $\left( y_t\right)$ and $\left(m_t\right)$ are, respectively,

${\rho }_y\left(s,t\right)=\frac{a_c\min \left(s,t\right)+1}{\sqrt{\left(a_{c}s+1\right)·\left(a_{c}t+1\right)}},\hspace{1em}\hspace{1em}{\rho }_m\left(s,t\right)=\frac{\min \left(s,t\right)}{\sqrt{s·t}}.\mathrm{∎}$

It is noticeable that previous theorem gives some additional information about the stochastic properties of the LSB process. Since $\left(m_t\right)$ is measurable in relation to the field ${\mathcal{F}}_{t-1}$, this series represents the component that describes the predictability and stability of the LSB process (see Figure 2 below). On the other hand, the innovation series $\left({\epsilon }_t\right)$ is the noise (deviation) component of the main LSB series $\left(y_t\right)$, compared to the martingale means $\left(m_t\right)$. Furthermore, according to Equation (7), it follows that the variances of the series $\left(y_t\right)$ and $\left(m_t\right)$ are non-constant and their order of dependence is equal to time moment $\left(t\right)$ when they are observed. Thereby, the correlation functions of $\left(y_t\right)$ and $\left(m_t\right)$ depend on both time variables $t,s$, thus indicating that these series are non-stationary. On the other hand, when $s>t\ge 0$, it follows:

$\underset{s\to t}{\lim }{\rho }_m\left(s,t\right)=\frac{t}{\sqrt{t^2}}=1,\hspace{1em}\hspace{1em}\underset{s\to t}{\lim }{\rho }_y\left(s,t\right)=\frac{a_{c}t+1}{\sqrt{{\left(a_{c}t+1\right)}^2}}=1,$

that is, both correlation functions satisfy the condition of L²-continuity.

Finally, we define another important LSB series, the so-called increments of the LSB process:

(8)$X_t=y_t-y_{t-1},\hspace{1em}t=1,\dots ,T.$

This series, according to Equations (1), (5), and (6), can be represented in the following way:

(9)$X_t={\epsilon }_t-{\theta }_{t-1}{\epsilon }_{t-1},$

where ${\theta }_t=1-q_t=I\left({\epsilon }_{t-1}^2\le c\right)$. Obviously, the series $\left(X_t\right)$ represents a bilinear, as well as stationary stochastic process with a random coefficient ${\theta }_t$, which is usually called Splitting moving average process (of order 1), that is, a Split-MA(1) process. Note that the term ‘split’, as in the case of the main LSB series $\left(y_t\right)$, is justified by the fact that $\left(X_t\right)$ operates in two modes:

• (a). If the fluctuations of innovations (${\epsilon }_t$) were emphasized in the previous moment in time, it follows ${\theta }_{t-1}=0$. Thus, equality (9) becomes $X_t={\epsilon }_t$.

• (b). Fluctuations ${\epsilon }_{t-1}$ whose square do not exceed the critical value $c$ imply ${\theta }_{t-1}=1$. Thus, the value of $X_t$ is given as a linear, integrated MA(1) process $X_t={\epsilon }_t-{\epsilon }_{t-1}.$.

Obviously, the Split-MA(1) series $\left(X_t\right)$ has a similar structure to MA(1) processes, which can be applied in their examination. Based on earlier assumptions, the basic properties of this series, obtained by some simple computations, can be expressed as follows:

Theorem 2. 

Let$\left(X_t\right)$be the time series defined by Equations (8) and (9). Then, the mean and the variance of this series are, respectively,

$E\left(X_t\right)=0, \mathrm{V}\mathrm{a}\mathrm{r}\left(X_t\right)=E\left(X_t^2\right)=2{\lambda }^2\left(b_c+1\right),$

where$b_c=1-a_c=P\left({\epsilon }_{t-1}^2\le c\right).$In addition, the correlation of this series is as follows:

${\rho }_X\left(h\right)≔\frac{\mathrm{C}\mathrm{o}\mathrm{v}\left(X_t,X_{t+h}\right)}{\mathrm{V}\mathrm{a}\mathrm{r}\left(X_t\right)}=\left\{\begin{array}{cc}1,\hfill & \hfill h=0\\ -\frac{b_c}{b_c+1},\hfill & \hfill h=\pm 1\\ 0,\hfill & \hfill \mathrm{otherwise},\end{array}\right.$

and it has the finite moments (of the order$n$):

$E\left(X_t^n\right)=\left\{\begin{array}{cc}0,\hfill & \hfill n=2k+1,\\ \frac{n!{\lambda }^n}{2}\left(2+n{b}_c\right),\hfill & \hfill n=2k,\end{array} k=0,1,2,\dots \right.∎$

It is worth noticing that, according to Equations (8) and (9), it follows that:

$y_t-y_{t-1}={\epsilon }_t-{\theta }_{t-1}{\epsilon }_{t-1},\hspace{1em}t=1,\dots ,T,$

which represents a nonlinear integrated autoregressive moving average (ARIMA) model with “temporary” components $\left({\theta }_{t-1}{\epsilon }_{t-1}\right)$. As will be seen later, this implies the specific structure of the (stationary) series $\left(X_t\right)$, as well as the other (non-stationary) components of the LSB process. Nevertheless, as will be seen below, due to its stationary quality, the Split-MA(1) process $\left(X_t\right)$ plays an important role in examining the stochastic properties of the LSB process, as well as in estimating its parameters. Figure 2 shows the realizations of all LSB time series obtained by Monte Carlo simulations. Note that non-stationary time series $\left(y_t\right)$ and $\left(m_t\right)$ (graphs on the left) can have very different trajectories, in contrast to realizations of the stationary series $\left(X_t\right)$ and $\left({\epsilon }_t\right)$ (graphs on the right).

3. Distributional Features of the LSB Process

Some important stochastic properties of LSB series, regarding their distributions and asymptotic behaviour, are discussed here. To this end, as mentioned earlier, the series $\left(X_t\right)$ plays a significant role due to its stationary quality. Using the method of characteristic functions (CFs), the basic stochastic properties of this series can be expressed in the following way.

Theorem 3. 

Let $\left(X_t\right)$ be the time series defined by Equations (8) and (9). Then, for any $x\in \mathbb{R}$ and $t=0, 1,\dots ,T$ , the CDF of the RVs $\left(X_t\right)$ is given by:

(10) $F_X\left(x\right)≔P\left\{X_t<x\right\}=\left(1-\frac{b_c}{2}\right)F_{\epsilon }\left(x\right)+\frac{b_c}{4}G\left(x\right),$

where $F_{\epsilon }\left(x\right)$ is the CDF of the Laplacian distributed RVs $\left({\epsilon }_t\right)$ , given by Equation (4), and

(11) $G\left(x\right)≔1+\mathrm{s}\mathrm{g}\mathrm{n}\left(x\right)\left(1-\left(\frac{\left|x\right|}{\lambda }+1\right)\exp \left(-\frac{\left|x\right|}{\lambda }\right)\right),\hspace{1em}x\in \mathbb{R}.$

Proof. 

To begin with, let us define the series of RVs ${\xi }_t={\theta }_t{\epsilon }_t$, $t=0,1,\dots ,T.$ Thereby, the RVs (${\theta }_t$) and (${\epsilon }_t$) are mutually independent, so it follows:

$E\left({\xi }_t\right)=E\left({\theta }_t\right)E\left({\epsilon }_t\right)=0,\hspace{1em}\hspace{1em}Var\left({\xi }_t\right)=E\left({\theta }_t^2\right)E\left({\epsilon }_t^2\right)=2b_c{\lambda }^2.$

It can easily be shown that it is $\mathrm{C}\mathrm{o}\mathrm{v}\left({\xi }_t,{\xi }_{t+h}\right)=0$, i.e., the RVs $\left({\xi }_t\right)$ are uncorrelated for any $h\ne 0$. Using the conditional probabilities, for the CDF of the RVs $\left({\xi }_t\right)$ one obtains:

$\begin{array}{ll}{F}_{\xi }\left(x\right)& ≔P\left\{{\xi }_t<x\right\}\\ & =P\left\{{\xi }_t<x|{\theta }_t=1\right\}·P\left\{{\theta }_t=1\right\}+P\left\{{\xi }_t<x|{\theta }_t=0\right\}·P\left\{{\theta }_t=0\right\}\\ & =P\left\{{\epsilon }_t<x\right\}·P\left\{{\theta }_t=1\right\}+P\left\{x>0\right\}·P\left\{{\theta }_t=0\right\}\\ & =b_c{F}_{\epsilon }\left(x\right)+\left(1-b_c\right)F_0\left(x\right),\end{array}$

where $F_0\left(x\right)=I\left(x>0\right)$ is the CDF of the RV $I_0\stackrel{as}{=}0$. According to this, the CF of the RVs $\left({\xi }_t\right)$ is obtained as follows:

$\begin{array}{ll}{\phi }_{\xi }\left(u\right)& ≔\underset{-\mathrm{\infty }}{\overset{+\mathrm{\infty }}{\int }}{e}^{iux}{F}_{\xi }\left(dx\right)=\underset{-\mathrm{\infty }}{\overset{+\mathrm{\infty }}{\int }}{e}^{iux}\left[b_c{F}_{\epsilon }+\left(1-b_c\right)F_0\right]\left(dx\right)\\ & =b_c{\phi }_{\epsilon }\left(u\right)+\left(1-b_c\right){\phi }_0\left(u\right).\end{array}$

Here, ${\phi }_{\epsilon }\left(u\right)={\left(1+{\lambda }^2{u}^2\right)}^{-1}$ and ${\phi }_0\left(u\right)\equiv 1$ are, respectively, CFs of the RVs $\left({\epsilon }_t\right)$ and $I_0\stackrel{as}{=}0$. Substituting these CFs into the previous equality gives:

${\phi }_{\xi }\left(u\right)=1+b_c\left(\frac{1}{1+{\lambda }^2{u}^2}-1\right).$

According to this equality and Equation (9), as well as the result of using the partial decomposition of rational functions, for the CF of the RVs $\left(X_t\right)$ one obtains:

(12)$\begin{array}{ll}{\phi }_X\left(u\right)& ={\phi }_{\epsilon }\left(u\right)·{\phi }_{\xi }\left(u\right)=\frac{1-b_c}{1+{\lambda }^2{u}^2}+\frac{b_c}{{\left(1+{\lambda }^2{u}^2\right)}^2}\\ & =\frac{1-b_c}{1+{\lambda }^2{u}^2}+\frac{b_c}{4}\left(\frac{1}{1+i\lambda u}+\frac{1}{1-i\lambda u}+\frac{1}{{\left(1+i\lambda u\right)}^2}+\frac{1}{{\left(1-i\lambda u\right)}^2}\right)\\ & =\frac{1-b_c}{1+{\lambda }^2{u}^2}+\frac{b_c}{2\left(1+{\lambda }^2{u}^2\right)}+\frac{b_c}{4}\left[{\left(1+2i\frac{\lambda u}{2}\right)}^{-2}+{\left(1-2i\frac{\lambda u}{2}\right)}^{-2}\right]\\ & =\frac{1-\frac{b_c}{2}}{1+{\lambda }^2{u}^2}+\frac{b_c}{4}\left[{\psi }_k\left(-\frac{\lambda u}{2}\right)+{\psi }_k\left(\frac{\lambda u}{2}\right)\right].\end{array}$

Here, ${\psi }_k\left(u\right)={\left(1-2iu\right)}^{-k/2}$ are the CFs of RVs ${\chi }_k^2$ with chi-square distribution and $k=4$ degrees of freedom (DF). According to Lévy’s correspondence theorem (see, c.f., [27]), the CDFs which correspond to CFs ${\psi }_k\left(-\lambda u/2\right)$ and ${\psi }_k\left(\lambda u/2\right)$ are chi-square-based RVs $-\lambda {\chi }_k^2/2$ and $\lambda {\chi }_k^2/2$, respectively. By applying the well-known facts about ${\chi }_k^2$ distribution, these CDFs are, respectively,

$\begin{array}{c}{G}_1\left(x\right)≔P\left\{-\frac{\lambda }{2}{\chi }_k^2<x\right\}=P\left\{{\chi }_k^2>-\frac{2x}{\lambda }\right\}=\left\{\begin{array}{ll}\left(\frac{x}{\lambda }+1\right)\exp \left(-\frac{x}{\lambda }\right),& x\le 0\\ 1,& x>0;\end{array}\right.\\ G_2\left(x\right)≔P\left\{\frac{\lambda }{2}{\chi }_k^2<x\right\}=P\left\{{\chi }_k^2<\frac{2x}{\lambda }\right\}=\left\{\begin{array}{ll}1-\left(\frac{x}{\lambda }+1\right)\exp \left(-\frac{x}{\lambda }\right),& x\ge 0\\ 0,& x<0.\end{array}\right.\end{array}$

Note that the functions $G_1\left(x\right)$ and $G_2\left(x\right)$ satisfy the equality $G\left(x\right)=G_1\left(x\right)+G_2\left(x\right)$, where $G\left(x\right)$ is the function given by Equation (11). Hence, by again applying the Lévy’s correspondence theorem to the last expression in (12), Equation (10) immediately follows. □

Remark 1. 

Let us emphasize that the function $G\left(x\right)$ is not an ‘ordinary’ CDF, but the sum of two chi-square-based CDFs $G_1\left(x\right)$ and $G_2\left(x\right)$ , which can also be seen in the left plot in Figure 3. Furthermore, by differentiating Equation (10) when $x\ne 0$ , and after some computations, the probability density function (PDF) of the series $\left(X_t\right)$ can be obtained as follows:

(13) $f_X\left(x\right)=\left(1-\frac{b_c}{2}\right)\frac{d{F}_{\epsilon }\left(x\right)}{dx}+\frac{b_c}{4}\frac{dG\left(x\right)}{dx}=\frac{1}{2\lambda }\exp \left(-\frac{\left|x\right|}{\lambda }\right)\left(1+\frac{b_c}{2}\left(\frac{\left|x\right|}{\lambda }-1\right)\right).$

Notice that PDF $f_X\left(x\right)$ differs from the PDF of the innovations $\left({\epsilon }_t\right)$ with respect to the last multiplicative term in Equation (13). Moreover, these PDFs become equal when $b_c=0$ , and both are symmetric functions, as shown in the right-hand plots of Figure 3 . ■

Using a similar procedure, the distributive properties of non-stationary LSB series $\left(m_t\right)$ and $\left(y_t\right)$ can be shown as follows.

Theorem 4. 

Let$\left(y_t\right)$and$\left(m_t\right)$be the time series defined by Equations (1) and (5), respectively, where$m_0\stackrel{as}{=}\mu \left(const\right)$. Then, for any$x\in \mathbb{R}$and$t=0, 1,\dots ,T$, the CDFs of the series$\left(m_t\right)$and$\left(y_t\right)$are, respectively,

(14)$F_m\left(x,t\right)≔P\left\{m_t<x\right\}=\underset{j=1}{\stackrel{t}{\otimes }}\left[\left(1-b_c\right)F_{\epsilon }^{\left(j\right)}\left(x\right)+b_c{F}_0\left(x\right)\right]\otimes F_{\mu }\left(x\right),$

(15)$F_y\left(x,t\right)≔P\left\{y_t<x\right\}=\underset{j=1}{\stackrel{t}{\otimes }}\left[\left(1-b_c\right)F_{\epsilon }^{\left(j\right)}\left(x\right)+b_c{F}_0\left(x\right)\right]\otimes F_{\mu }\left(x\right)\otimes F_{\epsilon }\left(x\right),$

where$“\otimes ”$is the convolution operator, and$F_0\left(x\right)$, $F_{\epsilon }^{\left(j\right)}\left(x\right)$are the CDFs of the RVs$I_0$, ${\epsilon }_t$, respectively. Additionally, the following convergences (in the distribution) hold when$T=+\infty$:

(16)$\frac{1}{\sqrt{t}}{m}_t\stackrel{d}{⟶}\mathcal{N}\left(0,2a_c{\lambda }^2\right),\hspace{1em}\hspace{1em}\frac{1}{\sqrt{t}}{y}_t\stackrel{d}{⟶}\mathcal{N}\left(0,2a_c{\lambda }^2\right),\hspace{1em}\hspace{1em}t\to +\infty .$

Proof. 

Similarly to the previous theorem, let us introduce RVs ${\eta }_t=q_t{\epsilon }_t$, $t=0,1,\dots ,T.$ It is easy to show that $\left({\eta }_t\right)$ is a series of mutually uncorrelated RVs, with $E\left({\eta }_t\right)=0,\mathrm{V}\mathrm{a}\mathrm{r}\left({\eta }_t\right)=2a_c{\lambda }^2,$ where $a_c=E\left(q_t\right)=P\left\{{\epsilon }_t^2>c\right\}=1-b_c.$ Applying the conditional probabilities again, for the CDF of $\left({\eta }_t\right)$ one obtains:

$\begin{array}{ll}{F}_{\eta }\left(x\right)& ≔P\left\{{\eta }_t<x\right\}\\ & =P\left\{{\eta }_t<x|q_t=1\right\}·P\left\{q_t=1\right\}+P\left\{{\eta }_t<x|q_t=0\right\}·P\left\{q_t=0\right\}\\ & =P\left\{{\epsilon }_t<x\right\}·P\left\{q_t=1\right\}+P\left\{x>0\right\}·P\left\{q_t=0\right\}\\ & =a_c{F}_{\epsilon }\left(x\right)+\left(1-a_c\right)F_0\left(x\right).\end{array}$

According to this CDF, the appropriate CF of the RVs $\left({\eta }_t\right)$ can be obtained as follows:

$\begin{array}{ll}{\phi }_{\eta }\left(u\right)& =\underset{-\mathrm{\infty }}{\overset{+\mathrm{\infty }}{\int }}{e}^{iux}{F}_{\eta }\left(dx\right)=\underset{-\mathrm{\infty }}{\overset{+\mathrm{\infty }}{\int }}{e}^{iux}\left[a_c{F}_{\epsilon }+\left(1-a_c\right)F_0\right]\left(dx\right)\\ & =a_c{\phi }_{\epsilon }\left(u\right)+\left(1-a_c\right){\phi }_0\left(u\right)=1+a_c\left(\frac{1}{1+{\lambda }^2{u}^2}-1\right)\\ & =\frac{1-b_c}{1+{\lambda }^2{u}^2}+b_c.\end{array}$

Now, by using Equation (5), for the CFs of the martingale means $\left(m_t\right)$ one obtains:

(17)${\phi }_m\left(u,t\right)={\phi }_{\mu }\left(u\right)\prod _{j=0}^{t-1}{\phi }_{\eta }\left(u\right)=e^{iu\mu }{\left(\frac{1-b_c}{1+{\lambda }^2{u}^2}+b_c\right)}^t,$

where ${\phi }_{\mu }\left(u\right)=e^{iu\mu }$ is CF of the RV $m_0\stackrel{as}{=}\mu .$ According to Equation (17) and the correspondence theorem of Lévy [27], Equation (14) follows. In a similar way, according to Equation (1), the CFs of the series $\left(y_t\right)$ can be written as follows:

(18)${\phi }_y\left(u,t\right)={\phi }_m\left(u\right)·{\phi }_{\epsilon }\left(u\right)=\frac{e^{iu\mu }}{1+{\lambda }^2{u}^2}·{\left(\frac{1-b_c}{1+{\lambda }^2{u}^2}+b_c\right)}^t.$

From here, Equation (15) is immediately followed by reapplying the correspondence theorem of Lévy.

Let us now prove the second part of the theorem, that is, Equation (16). First, using previous Equations (17) and (18), for the CFs of RVs $m_t/\sqrt{t}$ and $y_t/\sqrt{t}$, when $t=1,2,\dots ,$ one obtains:

$\begin{array}{l}{\phi }_m\left(\frac{u}{\sqrt{t}},t\right)=e^{iu\mu /\sqrt{t}}{\left[1+a_c\left(\frac{t}{t+{\lambda }^2{u}^2}-1\right)\right]}^t=e^{\frac{iu\mu }{\sqrt{t}}}{\left(1-a_c\frac{{\lambda }^2{u}^2}{t+{\lambda }^2{u}^2}\right)}^t,\\ {\phi }_y\left(\frac{u}{\sqrt{t}},t\right)=\frac{t{e}^{iu\mu /\sqrt{t}}}{t+{\lambda }^2{u}^2}{\left[1+a_c\left(\frac{t}{t+{\lambda }^2{u}^2}-1\right)\right]}^t=\frac{t{e}^{iu\mu /\sqrt{t}}}{t+{\lambda }^2{u}^2}{\left(1-a_c\frac{{\lambda }^2{u}^2}{t+{\lambda }^2{u}^2}\right)}^t.\end{array}$

From here, taking the limit values, when $t\to +\mathrm{\infty }$ and $u\in \mathbb{R}$ is a fixed (but arbitrary) value, the following is obtained:

$\underset{t\to +\mathrm{\infty }}{\lim }{\phi }_m\left(\frac{u}{\sqrt{t}},t\right)=\underset{t\to +\mathrm{\infty }}{\lim }{\phi }_y\left(\frac{u}{\sqrt{t}},t\right)=e^{-a_c{\lambda }^2{u}^2}.$

Thus, obtained limit is the obvious the CF of the Gaussian distribution $\mathcal{N}\left(0,2a_c{\lambda }^2\right)$, which confirms the convergences in (16). □

Remark 2. 

In the proofs of the previous two theorems, the uncorrelated series of RVs $\left({\xi }_t\right)$ and $\left({\eta }_t\right)$ can be interpreted as innovations with non-zero “optional” values. As the relation ${\xi }_t+{\eta }_t\stackrel{as}{=}{\epsilon }_t$ holds $,$ it is sufficient to take only one of these two series for consideration. Moreover, it is easy to see that these RVs have the following CDFs:

$\begin{array}{l}{F}_{\xi }\left(x\right)=b_c{F}_{\epsilon }\left(x\right)+\left(1-b_c\right)F_0\left(x\right),\\ F_{\eta }\left(x\right)=\left(1-b_c\right)F_{\epsilon }\left(x\right)+b_c{F}_0\left(x\right)=F_{\epsilon }\left(x\right)+F_0\left(x\right)-F_{\xi }\left(x\right).\end{array}$

Obviously, both of these functions are continuous almost everywhere, with the only breaking point$x=0$where ‘jumps’ of size$1-b_c$and$b_c$occur, respectively (for more details see Stojanović et al. [28,29]). For these reasons, the CDFs of the series$\left({\xi }_t\right)$and$\left({\eta }_t\right)$represent mixtures of the Laplacian and discrete distribution concentrated at zero, which we call the Contaminated Laplacian Distribution (CLD). This fact prevents the application of some of the standard procedures in the examination the properties of non-stationary series$\left(y_t\right)$and$\left(m_t\right)$.

On the other hand, asymptotic relations in Equation (16) show that even non-stationary time series$\left(m_t\right)$and$\left(y_t\right)$, obtained by non-Gaussian innovations$\left({\epsilon }_t\right)$, one can generate series$\left(m_t/\sqrt{t}\right)$and$\left(y_t/\sqrt{t}\right)$, which will converge to Gaussian distribution when$t\to +\infty$. These facts are of practical importance when applying the LSB process and can be easily observed by the convergence of the appropriate CFs${\phi }_m\left(u/\sqrt{t},t\right)$and${\phi }_y\left(u/\sqrt{t},t\right)$. As an illustration of this, convergences of the modulus of both of these CFs, for different time values, are shown in Figure 4. ∎

In the following, the asymptotic properties of some other time series, obtained by linear transformations of the non-stationary time series $\left(m_t\right)$ and $\left(y_t\right)$ are presented. It primarily refers to the possibility of finding their asymptotically normal (AN) distributions, and can be proven by the following proposition:

Theorem 5. 

Let us define, for an arbitrary $\alpha \ge 1$ , the $\alpha$ -mean time series:

${\stackrel{-}{M}}_{t;\alpha }=\frac{1}{t^{\alpha }}\sum _{j=1}^t{m}_j,\hspace{1em}{\stackrel{-}{Y}}_{t;\alpha }=\frac{1}{t^{\alpha }}\sum _{j=1}^t{y}_j,$

where $\left(y_t\right)$ and $\left(m_t\right)$ are the non-stationary time series given by (1) and (5), respectively. Then the following statements are valid:

• (i)

  When $1\le \alpha \le 3/2$ , both series ${\stackrel{-}{M}}_{t;\alpha }$ and ${\stackrel{-}{Y}}_{t;\alpha }$ have an asymptotically normal distribution, i.e., the following relations, when $t\to +\infty$ , are valid:

  (19) ${\stackrel{-}{M}}_{t;\alpha }~\mathcal{N}\left(\mu t^{1-\alpha },\frac{2a_c{\lambda }^2{t}^{3-2\alpha }}{3}\right),{\stackrel{-}{Y}}_{t;\alpha }~\mathcal{N}\left(\mu t^{1-\alpha },\frac{2a_c{\lambda }^2{t}^{3-2\alpha }}{3}\right).$

• (ii)

  When $\alpha >3/2$ , both series ${\stackrel{-}{M}}_{t;\alpha }$ and ${\stackrel{-}{Y}}_{t;\alpha }$ are asymptotically vanishing, i.e.,

  (20) ${\stackrel{-}{M}}_{t;\alpha }\stackrel{d}{⟶}{I}_0,{\stackrel{-}{Y}}_{t;\alpha }\stackrel{d}{⟶}{I}_0,t\to +\mathrm{\infty }.$

Proof. 

First, we prove the statement of the theorem in the case of series ${\stackrel{-}{M}}_{t;\alpha }$. According to the definition of series $\left(m_t\right)$, given by Equation (5), the following is obtained:

$\begin{array}{ll}{\stackrel{-}{M}}_{t;\alpha }& =\frac{1}{t^{\alpha }}{\displaystyle \sum _{j=1}^t{m}_j}=\frac{1}{t^{\alpha }}{\displaystyle \sum _{j=1}^t\left(m_0+{\displaystyle \sum _{k=0}^{j-1}{q}_k{\epsilon }_k}\right)}=\frac{1}{t^{\alpha }}\left(t{m}_0+{\displaystyle \sum _{j=0}^{t-1}\left(t-j\right)q_j{\epsilon }_j}\right)\\ & =t^{1-\alpha }{m}_0+{\displaystyle \sum _{k=1}^t\frac{k}{t^{\alpha }}{\eta }_{t-k}}.\end{array}$

Therefore, the time series ${\stackrel{-}{M}}_{t;\alpha }$ is represented as the sum of uncorrelated RVs ${\eta }_{t-k}$, when $k=1,\dots ,t$. Applying the well-known features of the CFs, as well as the CF of the series $\left({\eta }_t\right)$, for the CFs of RVs ${\stackrel{-}{M}}_{t;\alpha }$ one obtains:

$\begin{array}{lll}{\phi }_{\stackrel{-}{M};\alpha }\left(u,t\right)={\phi }_m& \left(\frac{u}{t^{\alpha -1}},0\right){\displaystyle \prod _{k=1}^t{\phi }_{\eta }\left(\frac{ku}{t^{\alpha }}\right)}& =e^{iu\mu t^{1-\alpha }}{\displaystyle \prod _{k=1}^t\left[1+a_c\left(\frac{t^{2\alpha }}{t^{2\alpha }+k^2{\lambda }^2{u}^2}-1\right)\right]}.\end{array}$

From here, by taking the logarithm of the function ${\phi }_{\stackrel{-}{M};\alpha }\left(u,t\right)$, the following function is obtained:

${\psi }_M\left(u,t,\alpha \right)≔\ln {\phi }_{\stackrel{-}{M};\alpha }\left(u,t\right)=iu\mu t^{1-\alpha }+\sum _{k=1}^t{f}_k\left(u,t,\alpha \right),$

where:

$f_k\left(u,t,\alpha \right)≔\ln \left[1+a_c\left(\frac{t^{2\alpha }}{t^{2\alpha }+k^2{\lambda }^2{u}^2}-1\right)\right]=\ln \left(1-\frac{a_c{k}^2{\lambda }^2{u}^2}{t^{2\alpha }+k^2{\lambda }^2{u}^2}\right).$

After some computation, we find that:

$\begin{array}{c}\frac{\partial f_k\left(0,t,\alpha \right)}{\partial u}={\left.-\frac{2a_c{k}^2{\lambda }^2{t}^{2\alpha }u}{\left(t^{2\alpha }+k^2{\lambda }^2{u}^2\right)\left(t^{2\alpha }+\left(1-a_c\right)k^2{\lambda }^2{u}^2\right)}\right|}_{u=0}=0,\\ \begin{array}{ll}\frac{{\partial }^2{f}_k\left(0,t,\alpha \right)}{\partial u^2}& =-{\left.\frac{2{a_{c}k}^2{\lambda }^2{t}^{2\alpha }\left(t^{4\alpha }+\left(a_c-2\right)k^2{\lambda }^2{t}^{2\alpha }{u}^2+3\left(a_c-1\right)k^4{\lambda }^4{u}^4\right)}{{\left(t^{2\alpha }+k^2{\lambda }^2{u}^2\right)}^2{\left(t^{2\alpha }+\left(1-a_c\right)k^2{\lambda }^2{u}^2\right)}^2}\right|}_{u=0}\\ & =-\frac{2a_c{k}^2{\lambda }^2}{t^{2\alpha }}<0.\end{array}\end{array}$

Thus, the functions $f_k\left(u,t,\alpha \right)$ have a local maximum at the point $u=0.$ By applying the Laplace approximation of functions $f_k\left(u,t,\alpha \right)$ at $u=0$ (see, c.f. [25,28]), it follows:

$\begin{array}{ll}{\psi }_M\left(u,t,\alpha \right)& =iu\mu t^{1-\alpha }+{\displaystyle \sum _{k=1}^t\left[\frac{{\partial }^2{f}_k\left(0,t,\alpha \right)}{\partial u^2}·\frac{u^2}{2}+{\mathcal{o}}_k\left(u^2\right)\right]}\\ & =iu\mu t^{1-\alpha }-{\displaystyle \sum _{k=1}^t\left[\frac{a_c{k}^2{\lambda }^2}{t^{2\alpha }}{u}^2+{\mathcal{o}}_k\left(t^{-2\alpha }{u}^2\right)\right]}\\ & =iu\mu t^{1-\alpha }-\frac{a_c{\lambda }^2{u}^2}{t^{2\alpha }}·\frac{t\left(t+1\right)\left(2t+1\right)}{6}+\mathcal{o}\left(t^{2-2\alpha }{u}^2\right).\end{array}$

Here $\mathcal{o}\left(z\right)$ denotes an infinitesimally small value of a higher order in relation to $z$, when $z\to 0$. Taking the asymptotic value in the previous expression, when $t\to +\infty$, gives:

${\psi }_M\left(u,t,\alpha \right)~\left\{\begin{array}{cc}iu\mu t^{1-\alpha }-a_c{\lambda }^2{t}^{3-2\alpha }/3,& 1\le \alpha \le \frac{3}{2}\\ 0,& \alpha >\frac{3}{2}.\end{array}\right.$

Substituting this expression in CFs ${\phi }_{\stackrel{-}{M};\alpha }\left(u,t\right)$, it is easy to see that the first asymptotic relation in Equation (19) is valid.

The proof for the series ${\stackrel{-}{Y}}_{t;\mathsf{\alpha }}$ can be conducted in an analogous manner. Using the previously proven facts and Equation (1), we find that:

$\begin{array}{ll}{\stackrel{-}{Y}}_{t;\alpha }& =\frac{1}{t^{\alpha }}{\displaystyle \sum _{j=1}^t\left(m_j+{\epsilon }_j\right)}={\stackrel{-}{M}}_{t;\alpha }+{\displaystyle \sum _{j=1}^t\frac{{\epsilon }_j}{t^{\alpha }}}=t^{1-\alpha }{m}_0+{\displaystyle \sum _{k=1}^t\frac{k}{t^{\alpha }}{\eta }_{t-k}}+{\displaystyle \sum _{k=0}^{t-1}\frac{{\epsilon }_{t-k}}{t^{\alpha }}}\\ & =t^{1-\alpha }{m}_0+\frac{{\epsilon }_t}{t^{\alpha }}+{\displaystyle \sum _{k=1}^t\left(1+k{q}_{t-k}\right)\frac{{\epsilon }_{t-k}}{t^{\alpha }}}.\end{array}$

Since ${\epsilon }_{t-k}$, $k=0,1,\dots ,t$, are mutually independent RVs, the CFs of time series ${\stackrel{-}{Y}}_{t;\alpha }$ can be obtained, after some calculations, as follows:

$\begin{array}{ll}{\phi }_{\stackrel{-}{Y};\alpha }\left(u,t\right)& ={\phi }_m\left(\frac{u}{t^{\alpha -1}},0\right){\phi }_{\epsilon }\left(\frac{u}{t^{\alpha }}\right){\displaystyle \prod _{k=1}^t}\left[\left(1-a_c\right){\phi }_{\epsilon }\left(\frac{u}{t^{\alpha }}\right)+a_c{\phi }_{\epsilon }\left(\frac{\left(k+1\right)u}{t^{\alpha }}\right)\right]\\ & =\frac{t^{2\alpha }{e}^{iu\mu t^{1-\alpha }}}{t^{2\alpha }+{\lambda }^2{u}^2}{\displaystyle \prod _{k=1}^t}[\frac{t^{2\alpha }}{t^{2\alpha }+{\lambda }^2{u}^2} + a_c(\frac{t^{2\alpha }}{t^{2\alpha }+{\lambda }^2{\left(k+1\right)}^2{u}^2}-\frac{t^{2\alpha }}{t^{2\alpha }+{\lambda }^2{u}^2}\left)\right]\\ & =e^{iu\mu t^{1-\alpha }}{\left(\frac{t^{2\alpha }}{t^{2\alpha }+{\lambda }^2{u}^2}\right)}^{t+1}{\displaystyle \prod _{k=1}^t}\left[1+a_c\left(\frac{t^{2\alpha }+{\lambda }^2{u}^2}{t^{2\alpha }+{\lambda }^2{\left(k+1\right)}^2{u}^2}-1\right)\right].\end{array}$

Based on this, and applying the same procedure as in the previous part of the proof, that is, by taking the logarithm of the function ${\phi }_{\stackrel{-}{Y};\alpha }\left(u,t\right)$, and by expanding the function ${\psi }_Y\left(u,t,\mathsf{\alpha }\right)≔\ln {\phi }_{\stackrel{-}{Y};\alpha }\left(u,t\right)$ at $u=0$, we have:

$\begin{array}{ll}{\psi }_Y\left(u,t,\alpha \right)& =iu\mu t^{1-\alpha }+\left(\mathrm{t}+1\right)\ln \left(\frac{t^{2\alpha }}{t^{2\alpha }+{\lambda }^2{u}^2}\right)\\ & \hspace{1em}+{\displaystyle \sum _{k=1}^t\mathrm{l}\mathrm{n}\left[1+a_c\left(\frac{t^{2\alpha }+{\lambda }^2{u}^2}{t^{2\alpha }+{\lambda }^2{\left(k+1\right)}^2{u}^2}-1\right)\right]}\\ & =iu\mu t^{1-\alpha }+\left(\mathrm{t}+1\right)\ln \left(1-\frac{{\lambda }^2{u}^2}{t^{2\alpha }+{\lambda }^2{u}^2}\right)\\ & \hspace{1em}-{\displaystyle \sum _{k=1}^t\left[\frac{a_c\left(k^2+2k\right){\lambda }^2{u}^2}{t^{2\alpha }}+{\mathcal{o}}_k\left(t^{-2\alpha }{u}^2\right)\right]}\\ & =iu\mu t^{1-\alpha }-\frac{{\lambda }^2{u}^2\left(t+1\right)}{t^{2\alpha }+{\lambda }^2{u}^2}-\frac{a_c{\lambda }^2{u}^2}{t^{2\alpha }}·\frac{t\left(t+1\right)\left(2t+7\right)}{6}+\mathcal{o}\left(t^{2-2\alpha }{u}^2\right).\end{array}$

Taking the asymptotic values, when $t\to +\infty$, the following is obtained:

${\psi }_Y\left(u,t,\alpha \right)~\left\{\begin{array}{cc}iu\mu t^{1-\alpha }-{\lambda }^2{u}^2\left(t^{1-2\alpha }+a_c{t}^{-2\alpha }+\frac{a_c{t}^{3-2\alpha }}{3}\right),& 1\le \alpha \le \frac{3}{2}\\ 0,& \alpha >\frac{3}{2}.\end{array}\right.$

Replacing this expression in the CFs ${\phi }_{\stackrel{-}{Y};\alpha }\left(u,t\right)$, the theorem is completely proven. □

Remark 3. 

As with the GSB process, the case $\alpha =3/2$ is of special interest in the previous theorem. Asymptotic relations (19) then give:

(21) $\frac{1}{t^{3/2}}\sum _{j=1}^t{m}_j\stackrel{d}{⟶}\mathcal{N}\left(0,\frac{2a_c{\lambda }^2}{3}\right),\frac{1}{t^{3/2}}\sum _{j=1}^t{y}_j\stackrel{d}{⟶}\mathcal{N}\left(0,\frac{2a_c{\lambda }^2}{3}\right),t\to +\mathrm{\infty }.$

The convergences (21), as in the case of GSB process, we called the central limit theorems (CLTs) for the LSB process. ■

4. Procedures for Parameter Estimation

In this section, we consider the estimations of the (unknown) parameters of the LSB process, that is the critical value $\left(c\right),$ the scale parameter $\left(\lambda \right)$, and the mean value $\left(\mu \right)$. In order to estimate the first parameter $c$, it should be emphasized that a series of increments $\left(X_t\right)$, as the only observable and stationary series of the LSB process, will be used. As has already been pointed out, this series is close to standard, linear MA models, but with Laplacian distributed innovations $\left({\epsilon }_t\right)$. For these reasons, the estimation procedures like those used for the mentioned time series will be used here. However, the specificity of the series $\left(X_t\right)$ requires some additional estimation procedures, as well as examination of the quality of the obtained estimators.

4.1. Moments-Based Estimators

According to Theorem 2, the first correlation of the series $\left(X_t\right)$ is as follows:

${\rho }_X\left(1\right)=-\frac{b_c}{1+b_c},\hspace{1em}0<b_c<1,$

and solving on $b_c$, the estimated value of this parameter is obtained:

(22)${\tilde{b}}_c=-\frac{{\widehat{\rho }}_X\left(1\right)}{1+{\widehat{\rho }}_X\left(1\right)}.$

Note that it is here:

${\widehat{\rho }}_X\left(1\right)=\left(\sum _{t=1}^T{X}_t{X}_{t-1}\right){\left(\sum _{t=1}^T{X}_t^2\right)}^{-1}$

the sample autocorrelation of the Split-MA(1) process $\left(X_t\right)$. According to Equation (22), it is easy to see that ${\tilde{b}}_c$ is the appropriate estimate if and only if

$0<{\tilde{b}}_c<1⟺-0.5<{\widehat{\rho }}_X\left(1\right)<0.$

Using the estimator ${\tilde{b}}_c,$ the squared scale parameter ${\lambda }^2$ can be estimated by applying the equality:

$\mathrm{V}\mathrm{a}\mathrm{r}\left(X_t\right)=E\left(X_t^2\right)=2{\lambda }^2\left(b_c+1\right).$

Applying Theorem 2 again, as well as the sample variance for $\left(X_t\right)$, we get the estimator:

(23)${\tilde{\lambda }}^2=\frac{1}{2T\left({\tilde{b}}_c+1\right)}\sum _{t=1}^T{X}_t^2.$

It is worth pointing out that in Stojanović et al. [21] it is proven that the estimators ${\tilde{b}}_c$ and ${\tilde{\lambda }}^2$ are strictly consistent if innovations $\left({\epsilon }_t\right)$ have a continuous distribution. Moreover, these estimators are asymptotically normal (AN) in the case of symmetrically distributed RVs $\left({\mathsf{\epsilon }}_t\right)$. Note that both conditions are satisfied for Laplace innovations $\left({\epsilon }_t\right)$ and, hence, it follows:

Theorem 6. 

The estimators ${\tilde{b}}_c$ and ${\tilde{\lambda }}^2$ , given by Equations (22) and (23), respectively, are strictly consistent estimates for the parameters $b_c$ and ${\lambda }^2$ , i.e.,

(24) ${\tilde{b}}_c\stackrel{as}{⟶}{b}_c,{\tilde{\lambda }}^2\stackrel{as}{⟶}{\lambda }^2,T\to +\mathrm{\infty }.$

In addition, ${\tilde{b}}_c$ and ${\tilde{\lambda }}^2$ are asymptotically normal estimators for $b_c$ and ${\lambda }^2$ , respectively, i.e.,

(25) $\sqrt{T}\left({\tilde{b}}_c-b_c\right)\stackrel{d}{⟶}\mathcal{N}\left(0,{\tilde{V}}_1\right),\sqrt{T}\left({\tilde{\lambda }}^2-{\lambda }^2\right)\stackrel{d}{⟶}\mathcal{N}\left(0,{\tilde{V}}_2\right),T\to +\mathrm{\infty },$

where:

$\begin{array}{l}{\tilde{V}}_1={\left(1+b_c\right)}^2\left(1+7b_c-3b_c\left(L+b_c\right)\right),\\ {\tilde{V}}_2=\frac{4{\lambda }^4\left(5+20b_c+2b_{c}L-3b_c^2\right)}{{\tilde{V}}_1},\\ L=\frac{1}{2{\lambda }^2}E\left({\theta }_t{\epsilon }_{t-1}^2\right)=\frac{1}{2{\lambda }^2}E\left({\theta }_t^2{\epsilon }_{t-1}^2\right)=1-e^{-\frac{\sqrt{c}}{\lambda }}\left(\frac{c}{2{\lambda }^2}+\frac{\sqrt{c}}{\lambda }+1\right).\end{array}$

Proof. 

According to Theorem 2, $\left(X_t\right)$ is a stationary series of RVs and its spectral density is:

$s_X\left(\omega \right)=\frac{1}{2\pi }\left({\gamma }_X\left(-1\right)e^{-i\omega }+{\gamma }_X\left(0\right)+{\gamma }_X\left(1\right)e^{i\omega }\right)=\frac{{\lambda }^2}{\pi }\left(1-2b_c\cos \omega +b_c\right),$

where ${\gamma }_X\left(h\right)=\mathrm{C}\mathrm{o}\mathrm{v}\left(X_t,X_{t+h}\right)=E\left(X_t{X}_{t+h}\right)$ is the covariance function of $\left(X_t\right).$ Obviously, the function $s_X\left(\omega \right)$ is continuous at the point $\omega =0$, so $\left(X_t\right)$ is also an ergodic series of RVs. Thus, by applying the mean ergodic theorem (see, c.f. [30], pp.131), for any $h\ge 0$, one obtains:

(26)$\frac{1}{T}\sum _{t=1}^T{X}_t{X}_{t+h}\stackrel{as}{⟶}{\gamma }_X\left(h\right),T\to +\mathrm{\infty }.$

Hence, by applying the convergence (26), when $h=0,1$, to the estimators ${\tilde{b}}_c$ and ${\tilde{\lambda }}^2$, given by Equations (22) and (23), respectively, the almost certain convergence in (24) is easily obtained.

Now, let us define series $Z_t=X_t{X}_{t+1}$, $t=0,1,\dots ,T-1,$ which has 1-dependent property (see, c.f. Definition 6.3.1 in [31], pp.245). Obviously, the series $\left(Z_t\right)$ is also stationary and ergodic, with:

$E\left(Z_t\right)=\mathrm{C}\mathrm{o}\mathrm{v}\left(X_t,X_{t+1}\right)={\gamma }_X\left(1\right)=-2b_c{\lambda }^2.$

According to the Cauchy–Swartz and Minkowski inequalities, it follows:

$E{\left|Z_t-{\gamma }_X\left(1\right)\right|}^3\le {\left[{\left(E{\left|Z_t\right|}^3\right)}^{1/3}+2b_c{\lambda }^2\right]}^3<+\mathrm{\infty },$

and by applying the Hoeffding–Robbins theorem [32], we have:

(27)$\sqrt{T}\left[Z_t-{\gamma }_X\left(1\right)\right]=T^{-\frac{1}{2}}\sum _{t=1}^T\left(X_t{X}_{t+1}+2b_c{\lambda }^2\right)\stackrel{d}{⟶}\mathcal{N}\left(0,A\right),\hspace{1em}\hspace{1em}T\to +\mathrm{\infty }.$

Here, after some computations, one obtains:

$\begin{array}{ll}A& =\mathrm{V}\mathrm{a}\mathrm{r}\left(Z_t\right)+2\mathrm{C}\mathrm{o}\mathrm{v}\left(Z_t,Z_{t+1}\right)=E\left(X_t^2{X}_{t+1}^2\right)+2E\left(X_t{X}_{t+1}^2{X}_{t+2}\right)-3{\gamma }_X^2\left(1\right)\\ & =4{\lambda }^4\left(1+7b_c+b_{c}L\right)-16{\lambda }^4{b}_{c}L+12{\lambda }^4{b}_c^2\\ & =4{\lambda }^4\left(1+7b_c-3b_c\left(L+b_c\right)\right).\end{array}$

In addition, by again applying the mean ergodic theorem, i.e., Equation (26) when $h=0$, it follows:

(28)$\frac{1}{T}\sum _{t=1}^T{X}_t^2\stackrel{as}{⟶}E\left(X_t^2\right)=2{\lambda }^2\left(b_c+1\right),\hspace{1em}\hspace{1em}T\to +\mathrm{\infty }.$

Thus, the convergences (27) and (28) imply:

(29)$\sqrt{T}\left[{\widehat{\rho }}_X\left(1\right)-{\rho }_X\left(1\right)\right]\stackrel{d}{⟶}\mathcal{N}\left(0,V_0\right),\hspace{1em}\hspace{1em}T\to +\mathrm{\infty },$

wherein:

$V_0=\frac{A}{{\gamma }_X^2\left(0\right)}=\frac{1+7b_c-3b_c\left(L+b_c\right)}{{\left(1+b_c\right)}^2}.$

Finally, according to the continuity of convergence in distribution (see, c.f. [33], pp.118), it follows:

$\sqrt{T}\left({\tilde{b}}_c-b_c\right)\stackrel{d}{⟶}\mathcal{N}\left(0,{\tilde{V}}_1\right),\hspace{1em}\hspace{1em}T\to +\mathrm{\infty },$

wherein:

${\tilde{V}}_1={\left(\frac{d{\rho }_X\left(1\right)}{d{b}_c}\right)}^{-2}{V}_0={\left(1+b_c\right)}^4{V}_0.$

In that way, the first convergence in (25) is proved.

To prove the second convergence in (25), note that series $\left(X_t^2\right)$ is also 1-dependent. According to previous mentioned inequalities, as well as Theorem 2, one obtains:

$E{\left|X_t^2-{\gamma }_X\left(0\right)\right|}^3\le {\left[{\left(E{\left|X_t\right|}^6\right)}^{1/3}-2{\lambda }^2\left(b_c+1\right)\right]}^3<+\mathrm{\infty }.$

Thus, by applying again the Hoeffding–Robbins theorem [32], we get:

(30)$T^{-\frac{1}{2}}\sum _{t=1}^T\left(X_t^2-2{\lambda }^2\left(b_c+1\right)\right)\stackrel{d}{⟶}\mathcal{N}\left(0,D\right),\hspace{1em}\hspace{1em}T\to +\mathrm{\infty },$

where, after some computations, one obtains:

$D=\mathrm{V}\mathrm{a}\mathrm{r}\left(X_t^2\right)+2\mathrm{C}\mathrm{o}\mathrm{v}\left(X_t^2,X_{t+1}^2\right)=4{\lambda }^4\left(5+20b_c+2b_{c}L-3b_c^2\right).$

Finally, according to the first almost certain convergence in (24), as well as the previously obtained convergence (30), it follows:

$\sqrt{T}{\tilde{\lambda }}^2=\frac{1}{2\sqrt{T}\left({\tilde{b}}_c+1\right)}\sum _{t=1}^T{X}_t^2\stackrel{d}{⟶}\mathcal{N}\left({\lambda }^2,{\tilde{V}}_2\right),\hspace{1em}\hspace{1em}T\to +\mathrm{\infty },$

where ${\tilde{V}}_2=D/{\tilde{V}}_1$. Thus, the second convergence in (25) is obtained. ■

Remark 4.

Using the estimators${\tilde{b}}_c$and${\tilde{\lambda }}^2$, the appropriate estimator of the critical value$c=\tilde{c}$can be easily obtained as a solution to the equation:

$P\left\{{\epsilon }_t^2\le c\right\}={\tilde{b}}_c,$

or, equivalently,

(31)$F_{\epsilon }\left(\sqrt{c}\right)-F_{\epsilon }\left(-\sqrt{c}\right)=2F_{\epsilon }\left(\sqrt{c}\right)-1={\tilde{b}}_c⟺\tilde{c}={\left[F_{\epsilon }^{-1}\left(\frac{{\tilde{b}}_c+1}{2}\right)\right]}^2.$

Here,$F_{\epsilon }\left(x\right)$is the CDF of the Laplacian innovations$\left({\epsilon }_t\right)$, given by Equation (4), and we also note that in Equation (31) the property of symmetry of its distribution is used. It is also noticeable that, from the fact that$\tilde{c}$is a continuous function of the estimator${\tilde{b}}_c$, the strong consistency and AN property of this estimator follows (see, c.f. [33], pp. 24). ∎

4.2. Gauss–Newton and Maximum Likelihood Estimators

It is worth pointing out that ${\tilde{b}}_c$ and ${\tilde{\lambda }}^2$ are not the most efficient estimators for $b_c$ and ${\lambda }^2$. In order to obtain more efficient estimators for the observed parameters, similarly to the case of the GSB process (see, c.f. [21,25]), the so-called modified Gauss–Newton method of parameter estimation of nonlinear functions can be used. First, we can write Equation (9) as:

${\epsilon }_t=X_t+{\theta }_{t-1}{\epsilon }_{t-1},\hspace{1em}\hspace{1em}t=1,\dots ,T,$

or, in the functional form,

${\epsilon }_t\left(X,\theta \right)=X_t+{\theta }_{t-1}{\epsilon }_{t-1}\left(X,\theta \right).$

Let ${\tilde{b}}_c$ be the estimated value of the parameter $b_c$, obtained by the previously described procedure. According to this, for each $t=1,2,\dots ,T$, we can recursively compute:

(32)${\tilde{\theta }}_t≔I\left({\epsilon }_{t-1}^2\left(X,\tilde{\theta }\right)\le \tilde{c}\right),{\epsilon }_t\left(X,\tilde{\theta }\right)≔X_t+{\tilde{\theta }}_{t-1}{\epsilon }_{t-1}\left(X,\tilde{\theta }\right),$

with the previously mentioned initial values ${\tilde{\theta }}_0=1$, ${\epsilon }_0\left(X,\tilde{\theta }\right)={\epsilon }_{-1}\left(X,\tilde{\theta }\right)=0.$ Additionally, if we define a series of RVs:

(33)$W_t\left(X,\theta \right)={\theta }_t{W}_{t-1}\left(X,\theta \right)+{\epsilon }_{t-1}\left(X,\theta \right),$

where $t=1,2,\dots ,T,$ and $W_0\left(X,\tilde{\theta }\right)=0$, it is obvious that RVs $W_t\left(X,\theta \right)$ are ${\mathcal{F}}_{t-1}$ adapted, and thus independent of ${\epsilon }_t$ and ${\theta }_{t+1}$. Furthermore, according to mentioned properties of RVs $\left({\theta }_t\right)$ and $\left({\epsilon }_t\right)$, it follows that $\left(W_t\left(X,\theta \right)\right)$ is a stationary and ergodic series of RVs with:

$E\left(W_t\left(X,\theta \right)\right)=0,\mathrm{V}\mathrm{a}\mathrm{r}\left(W_t\left(X,\theta \right)\right)=E{\left(W_t\left(X,\theta \right)\right)}^2=\frac{2{\lambda }^2}{1-b_c},$

and correlation function ${\rho }_W\left(h\right)=b_c^{\left|h\right|},h=0,\pm 1,\dots$ Now, using a procedure similar to the one in Lawrence and Lewis [34], we define the so-called residual series:

(34)$R_t\left(X,\theta \right)=W_t\left(X,\theta \right)-b_c{W}_{t-1}\left(X,\theta \right).$

It can be shown easily that $R_t\left(X,\theta \right)$ are also ${\mathcal{F}}_{t-1}$ adapted, but mutually non-correlated RVs (see, c.f. Jovanović et al. [25]). Thus, Equation (34) represents a linear autoregressive (AR) process $\left(W_t\left(X,\theta \right)\right)$ with innovations $\left(R_t\left(X,\theta \right)\right)$. From here, by applying the usual regression procedure, another estimator of the unknown, threshold parameter $b_c\in \left(\mathrm{0,1}\right)$ can be obtained by the equality:

(35)${\widehat{b}}_c=\left(\sum _{t=0}^{T-1}{W}_t\left(X,\tilde{\theta }\right)W_{t+1}\left(X,\tilde{\theta }\right)\right){\left(\sum _{t=0}^{T-1}{W}_t^2\left(X,\tilde{\theta }\right)\right)}^{-1},$

where, according to Equations (33) and (34),

$\begin{array}{l}{W}_t\left(X,\tilde{\theta }\right)≔{\left.W_t\left(X,\theta \right)\right|}_{\theta =\tilde{\theta }}={\tilde{\theta }}_t{W}_{t-1}\left(X,\tilde{\theta }\right)+{\epsilon }_{t-1}\left(X,\tilde{\theta }\right),\\ R_t\left(X,\tilde{\theta }\right)≔{\left.R_t\left(X,\theta \right)\right|}_{\theta =\tilde{\theta }}=W_t\left(X,\tilde{\theta }\right)-b_c{W}_{t-1}\left(X,\tilde{\theta }\right).\end{array}$

As in the case of previous estimators, the estimator ${\widehat{b}}_c$ gives the estimator of critical value $c=\widehat{c}$ as a solution of the equation:

(36)$P\left\{{\epsilon }_t^2\le c\right\}={\widehat{b}}_c⟺\widehat{c}={\left[F_{\epsilon }^{-1}\left(\frac{{\widehat{b}}_c+1}{2}\right)\right]}^2.$

Notice that using the previously obtained estimators ${\tilde{b}}_c$ (as well as ${\widehat{b}}_c$), that is, the modelled values of innovations $\left({\epsilon }_t\right)$, defined by Equation (32), the scale parameter $\lambda$ can be estimated. It is well-known that for the Laplace innovations $\left({\epsilon }_t\right)$ the most efficient estimates are obtained by using the maximum likelihood (ML) estimators (see, c.f. [35,36]). In the case of our LSB process, the ML estimator can be obtained according to Equations (1) and (2), that is, based on the maximization of the log-likelihood function:

$\mathcal{l}\left(y_1,\dots ,y_T;\lambda \right)=\ln \prod _{t=1}^T\left[\frac{1}{2\lambda }\exp \left(-\frac{\left|y_t-m_t\right|}{\lambda }\right)\right]=-T\ln \left(2\lambda \right)-\frac{1}{\mathsf{\lambda }}\sum _{t=1}^T\left|{\epsilon }_t\left(X,\tilde{\theta }\right)\right|.$

By solving the equation $\mathcal{l}\left(y_1,\dots ,y_T;\lambda \right)/\partial \lambda =0$ and using Equation (32), the estimator of scale parameter $\lambda$ is obtained as the mean absolute deviation (MAD):

(37)$\widehat{\lambda }=\frac{1}{T}\sum _{t=1}^T\left|{\epsilon }_t\left(X,\tilde{\theta }\right)\right|.$

Strict consistency and AN of the estimates ${\widehat{b}}_c$ and $\widehat{\lambda }$ can be proven as follows:

Theorem 7. 

Estimators ${\widehat{b}}_c$ and $\widehat{\lambda }$ , defined by Equations (34) and (36), respectively, are strictly consistent estimators for the parameters $b_c$ and $\lambda$ , i.e.,

(38) ${\widehat{b}}_c\stackrel{as}{⟶}{b}_c,\widehat{\lambda }\stackrel{as}{⟶}\mathsf{\lambda },T\to +\mathrm{\infty }.$

Moreover, ${\widehat{b}}_c$ and $\widehat{\lambda }$ are asymptotically normal estimators for $b_c$ and $\lambda$ , respectively, i.e.,

(39) $\sqrt{T}\left({\widehat{b}}_c-b_c\right)\stackrel{d}{⟶}\mathcal{N}\left(0,1-b_c^2\right),\sqrt{T}\left(\widehat{\lambda }-\mathsf{\lambda }\right)\stackrel{d}{⟶}\mathcal{N}\left(0,{\lambda }^2\right),T\to +\mathrm{\infty }.$

Proof. 

According to definitions of the residuals ($R_t(X,\tilde{\theta })$) and the estimator ${\widehat{b}}_c$, given by Equations (34) and (35), respectively, we obtain:

${\widehat{b}}_c-b_c=\frac{\sum _{t=0}^{T-1}{W}_t\left(X,\tilde{\theta }\right)W_{t+1}\left(X,\tilde{\theta }\right)-b_c\sum _{t=0}^{T-1}{W}_t^2\left(X,\tilde{\theta }\right)}{\sum _{t=1}^T{W}_t^2\left(X,\tilde{\theta }\right)}$

(40)$=\frac{\sum _{t=0}^{T-1}{W}_t\left(X,\tilde{\theta }\right)R_{t+1}\left(X,\tilde{\theta }\right)}{\sum _{t=0}^{T-1}{W}_t^2\left(X,\tilde{\theta }\right)}.$

Similarly to the previous theorem, it can be shown that $\left(W_t\right(X,\tilde{\theta }\left)\right),\left(R_t\right(X,\tilde{\theta }\left)\right),$ and (${\epsilon }_t(X,\tilde{\theta })$) are stationary and ergodic time series. Thus, by applying the mean ergodic theorem, it follows:

(41)$\begin{array}{l}\frac{1}{T}{\displaystyle \sum _{t=0}^{T-1}{W}_t\left(X,\tilde{\theta }\right)R_{t+1}\left(X,\tilde{\theta }\right)}\stackrel{as}{⟶}E\left(W_t\left(X,\tilde{\theta }\right)R_{t+1}\left(X,\tilde{\theta }\right)\right)=0\hspace{1em}\hspace{1em},T\to +\mathrm{\infty },\\ \frac{1}{T}{\displaystyle \sum _{t=0}^{T-1}{W}_t^2\left(X,\tilde{\theta }\right)}\stackrel{as}{⟶}E\left(W_t^2\left(X,\tilde{\theta }\right)\right)=\frac{2{\lambda }^2}{1-b_c},\hspace{1em}\hspace{1em}T\to +\mathrm{\infty },\end{array}$

$\frac{1}{T}\sum _{t=0}^{T-1}\left|{\epsilon }_t\left(X,\tilde{\theta }\right)\right|\stackrel{as}{⟶}E\left|{\epsilon }_t\left(X,\tilde{\theta }\right)\right|=\lambda ,T\to +\mathrm{\infty }.$

These convergences obviously imply the almost certain convergences in Equation (38).

Furthermore, according to Equation (40), in the following decomposition one obtains:

$\sqrt{T}\left({\widehat{b}}_c-b_c\right)=\frac{T^{-\frac{1}{2}}{\mathbb{M}}_{T-1}}{T^{-1}{\mathbb{W}}_{T-1}},$

wherein:

${\mathbb{M}}_{T-1}≔\sum _{t=0}^{T-1}{W}_t\left(X,\tilde{\theta }\right)R_{t+1}\left(X,\tilde{\theta }\right),{\mathbb{W}}_{T-1}≔\sum _{t=0}^{T-1}{W}_t^2\left(X,\tilde{\theta }\right).$

As for each $t=1,2,\dots ,T$,

$E\left({\mathbb{M}}_t|{\mathcal{F}}_{t-1}\right)={\mathbb{M}}_{t-1}+W_t\left(X,\tilde{\theta }\right)E\left(R_{t+1}\left(X,\tilde{\theta }\right)\right)={\mathbb{M}}_{t-1},$

the series $\left({\mathbb{M}}_t\right)$ is a martingale. Applying Billingsley’s CLT for martingales [37], we get:

(42)$T^{-\frac{1}{2}}{\mathbb{M}}_{T-1}\stackrel{d}{⟶}\mathcal{N}\left(0,D_1\right),$

wherein, according to Equation (34),

$\begin{array}{ll}{D}_1& =E{\left(W_t\left(X,\tilde{\theta }\right)R_{t+1}\left(X,\tilde{\theta }\right)\right)}^2=E\left(W_t^2\left(X,\tilde{\theta }\right)\right)E{\left(W_{t+1}\left(X,\theta \right)-b_c{W}_t\left(X,\theta \right)\right)}^2\\ & =\frac{2{\lambda }^2}{1-b_c}E\left(W_{t+1}^2\left(X,\tilde{\theta }\right)-2b_c{W}_t\left(X,\theta \right)W_{t+1}\left(X,\theta \right)+b_c^2{W}_t^2\left(X,\tilde{\theta }\right)\right)\\ & =\frac{2{\lambda }^2}{1-b_c}\left(\frac{2{\lambda }^2\left(1+b_c^2\right)}{1-b_c}-2b_c\mathrm{C}\mathrm{o}\mathrm{v}\left(W_t\left(X,\theta \right),W_{t+1}\left(X,\theta \right)\right)\right)\\ & =\frac{2{\lambda }^2}{1-b_c}\left(\frac{2{\lambda }^2\left(1+b_c^2\right)}{1-b_c}-\frac{4{\lambda }^2{b}_c^2}{1-b_c}\right)\\ & =\frac{4{\lambda }^4\left(1+b_c\right)}{1-b_c}.\end{array}$

On the other hand, the almost certain convergence in Equation (41) is the same as:

(43)$T^{-1}{\mathbb{W}}_{T-1}\stackrel{as}{⟶}\frac{2{\lambda }^2}{1-b_c},T\to +\mathrm{\infty }.$

Thus, the convergences (42) and (43) give the first convergence in (39). Finally, according to Equation (33), and after some computations, one obtains:

$\begin{array}{ll}\mathrm{V}\mathrm{a}\mathrm{r}\left(\widehat{\lambda }\right)& =\frac{1}{T^2}{\displaystyle \sum _{t=1}^T}\left[E\left({\epsilon }_t^2\left(X,\tilde{\theta }\right)\right)-{\lambda }^2\right]=\frac{1}{T}\left(E{\left(W_{t+1}\left(X,\theta \right)-{\theta }_{t+1}{W}_t\left(X,\theta \right)\right)}^2-{\lambda }^2\right)\\ & =\frac{1}{T}\left(E\left(W_{t+1}^2\left(X,\theta \right)\right)-2E\left({\theta }_{t+1}{W}_{t+1}\left(X,\theta \right)W_t\left(X,\theta \right)\right)+E\left({\theta }_{t+1}{W}_t^2\left(X,\theta \right)\right)-{\lambda }^2\right)\\ & =\frac{1}{T}\left(\frac{2{\lambda }^2\left(1+b_c\right)}{1-b_c}-2E\left({\theta }_{t+1}\left({\theta }_{t+1}{W}_t\left(X,\theta \right)+{\epsilon }_t\left(X,\theta \right)\right)W_t\left(X,\theta \right)\right)-{\lambda }^2\right)\\ & =\frac{1}{T}\left(\frac{2{\lambda }^2\left(1+b_c\right)}{1-b_c}-2E\left({\theta }_{t+1}{W}_t^2\left(X,\theta \right)\right)-{\lambda }^2\right)\\ & =\frac{1}{T}\left(\frac{2{\lambda }^2\left(1+b_c\right)}{1-b_c}-\frac{4{\lambda }^2{b}_c}{1-b_c}-{\lambda }^2\right)=\frac{{\lambda }^2}{T},\end{array}$

and the second convergence in (39) immediately follows □

Remark 5. 

Similarly to Theorem 6, the estimator$\widehat{c},$defined by Equation (36), is a continuous function of${\widehat{b}}_c$. Thus, the strong consistency and AN properties of$\widehat{c}$are also valid. Additionally, based on the scale parameter estimator$\widehat{\lambda }$given by Equation (37), the estimator${\widehat{\lambda }}^2$can be defined. It is also strictly consistent and AN but, at the same time, the most efficient estimator for the parameter${\lambda }^2.$Namely, by applying the continuity of the convergence in distribution [33], (pp. 118), it follows that:

$\sqrt{T}\left({\widehat{\lambda }}^2-{\mathsf{\lambda }}^2\right)\stackrel{d}{⟶}\mathcal{N}\left(0,{\widehat{V}}_2\right),\hspace{1em}\hspace{1em}T\to +\mathrm{\infty },$

where${\widehat{V}}_2={\left({\left({\lambda }^2\right)}^{\prime }\right)}^2{\lambda }^2=4{\lambda }^4$. Thus, the variance of the estimator${\widehat{\lambda }}^2$is$4{\lambda }^4/T$and it is the minimum of the variance for the unbiased estimator of the parameter${\lambda }^2$, obtained from Rao–Cramér’s inequality (see, c.f. [35,36]). In that way, for the previously defined estimators, inequalities:

${\widehat{V}}_1\left(b_c\right)\le {\tilde{V}}_1\left(b_c\right), {\widehat{V}}_2\left({\lambda }^2\right)\le {\tilde{V}}_2\left(b_c,{\lambda }^2\right)$

hold, where we denoted${\widehat{V}}_1\left(b_c\right)=1-b_c^2$. This means that${\widehat{b}}_c$and${\widehat{\lambda }}^2$are more efficient estimators than${\tilde{b}}_c$and${\tilde{\lambda }}^2$, respectively, which can also be seen inFigure 5, where asymptotic variances of these estimators are shown. ∎

4.3. Estimators of the Mean Value

In the last part of this section, we consider estimators for the mean value of the LSB process $\mu =E\left(y_t\right)$. For this aim, the sample mean of series is usually used $\left(y_t\right)$:

(44)$\tilde{\mu }:={\stackrel{-}{y}}_T=\frac{1}{T}\sum _{t=1}^T{y}_t,$

which is an unbiased estimator, that is, $E\left(\tilde{\mu }\right)=E\left({\stackrel{-}{y}}_T\right)=\mu$, but its variance is unbounded. To prove this fact, we use the $\alpha$-mean series ${\stackrel{-}{Y}}_{T;\alpha }$, defined in Theorem 5, in the case where $\alpha =1$. Then, the estimator $\tilde{\mu }$ can be represented as a sum of uncorrelated RVs:

$\tilde{\mu }=m_0+\frac{1}{T}\left[\sum _{k=1}^T\left(1+k{q}_{T-k}\right){\epsilon }_{T-k}+{\epsilon }_T\right],$

and using the same procedure as in Theorem 5, the variance of $\tilde{\mu }$ is as follows:

$\tilde{V}≔V\left(\tilde{\mu }\right)=\frac{2a_c{\lambda }^{2}T}{3}+\mathcal{O}\left(T^{-1}\right)⟶+\mathrm{\infty },\hspace{1em}\hspace{1em}T\to +\mathrm{\infty }.$

Let us emphasize that the variance $\tilde{V}=V\left(\tilde{\mu }\right)$ is asymptotically equal to that in Equation (19), when $\alpha =1$, as expected.

In order to obtain a more efficient estimator for the parameter $\mu$, a sample mean of the mean series ${\stackrel{-}{y}}_t$, when $t=1,\dots ,T$, can be taken, i.e.,

(45)$\widehat{\mu }:=\frac{1}{T}\sum _{t=1}^T{\stackrel{-}{y}}_t=\frac{1}{T}\sum _{t=1}^T{\omega }_t{y}_t,$

where ${\omega }_t≔H\left(T\right)-H(t-1)$ and $H\left(t\right)≔{\sum }_{j=1}^t{j}^{-1},t=1,\dots ,T$ are the harmonic numbers, with $H\left(0\right)=0$. The estimator $\widehat{\mu }$ is also unbiased for the parameter $\mu$, but its weights are more pronounced at the ‘older’ time points of the realization of the series $\left(y_t\right)$. This is consistent with the fact that the covariances of RVs $y_t$ depend on these ‘older’ time indices. For these reasons, the estimator $\widehat{\mu }$ is more efficient than $\tilde{\mu }$, which can be shown using a procedure similar to that of Jovanović et al. [25]. First, we represent estimator $\widehat{\mu }$ as a sum of uncorrelated RVs:

(46)$\widehat{\mu }=\frac{1}{T}\left[m_0\sum _{t=1}^T{\omega }_t+\sum _{j=0}^{T-1}\left(q_j{\epsilon }_j\sum _{t=j+1}^T{\omega }_t\right)+\sum _{t=1}^T{\omega }_t{\epsilon }_t\right],$

and after some computations (see, for more detail [25]), the variance of $\widehat{\mu }$ can be obtained as follows:

$\widehat{V}≔V\left(\widehat{\mu }\right)=2a_c{\lambda }^2{H}^2\left(T\right)+\mathcal{o}\left(H^{-2}\left(T\right)\right)⟶+\mathrm{\infty },\hspace{1em}\hspace{1em}T\to +\mathrm{\infty }.$

In this way, it can be seen that the asymptotic variance $\widehat{V}:=V\left(\widehat{\mu }\right)$ is also unbounded, but is asymptotically smaller than $\tilde{V}=V\left(\tilde{\mu }\right)$, because:

$\underset{T\to +\infty }{\lim }\frac{V\left(\widehat{\mu }\right)}{V\left(\tilde{\mu }\right)}=\underset{T\to +\infty }{\lim }\frac{H^2\left(T\right)}{T}=0.$

Thus, the estimator $\widehat{\mu }$ is (asymptotically) more efficient than $\tilde{\mu }$.

5. Numerical Simulations of the LSB Estimators

In this section, the parameter estimation procedures of the LSB model are discussed, where $N=1500$ Monte Carlo simulations of the basic LSB series are generated, of the length $T=1000$. At the same time, the main goal is to examine the quality of previously proposed estimators, that is, their asymptotic properties, which were analysed and shown in the previous section. To this end, appropriate estimation errors and normality testing procedures are also used. The summarized values of the estimated parameters, i.e., the mean (Mean), minimum (Min.), and maximum (Max.), as well as the corresponding mean square estimated errors (MSEE), are shown in the left part of Table 1. In addition, the obtained estimates were tested for their AN properties using Anderson–Darling and Cramér–von Mises normality tests. The appropriate test statistics (denoted by AD and W, respectively), as well as their $p$-values, were computed using the R-package “nortest” [38] and are presented in the right part of Table 1.

Based on these obtained estimates, it is evident that almost all of them have the AN property, which is also confirmed by the previous theoretical results. It is worth pointing out that even the estimates of the mean values $\tilde{\mu }$ and $\widehat{\mu }$, obtained by the realizations of the non-stationary series $\left(y_t\right)$, have the AN properties. This is already explained by the theoretical findings provided by Theorems 4 and 5, which describe the AN properties of this series. Therefore, the effectiveness of both of these estimators due to their non-stationarity is not pronounced and, due to unlimited variance, there are a wide range of obtained estimated values. On the other hand, we note that the AN property is not particularly emphasized in the case of estimates of the critical value $\left(c\right)$. This is a consequence of the three-step procedure for estimating this parameter, because the estimates of the parameter $c$ are obtained after the estimates of the parameters $b_c$ and $\lambda$ have been computed. Nevertheless, it is clear that, in accordance with the previous theoretical findings, first of all in Theorem 7, the estimate $\widehat{c}$ is more efficient and has a more pronounced AN property than estimate $\tilde{c}$.

Finally, in the case of estimates of the scale parameter $\left(\lambda \right)$, their efficiency and AN properties are clearly visible. At the same time, it should be noticed that the moment-based estimate $\tilde{\lambda }$ has a somewhat slower efficiency and AN property compared to the ML estimate $\widehat{\lambda }$, obtained using Equation (37) and the modelled innovations $\left({\epsilon }_t\right)$. This is fully consistent with the previously explained theoretical findings, that is, the definitions of both of these two estimates, as well as Theorems 6 and 7. Some visual confirmation of these facts can be seen in Figure 6, where the histograms of the frequency distribution of the obtained estimates are shown. Thus, for instance, in all cases (with the exception $\tilde{c}$) the presence of AN properties clearly can be observed, as well as the efficiency of the obtained estimates.

6. Application in Dynamic Analysis of the World Oil and Gas Market

In the following, the application of the LSB process in the dynamic analysis of prices and trading volumes of crude oil and natural gas on the world market is considered. We emphasize that the temporal dynamics of these two energy sources is of particular importance, and they are greatly influenced by global external factors, such as, for instance, the recent COVID-19 pandemic and the war in Ukraine. Precisely for these reasons, it can be assumed that all of these factors cause permanent and pronounced fluctuations in the dynamics of the price and trading volume of these two energy products, which can be seen in the following Figure 7. As mentioned in the introduction, it will be shown here that the LSB process can be an appropriate stochastic model for describing dynamics of these kinds. To this end, based on official data from the National Association of Securities Dealers Automated Quotations (NASDAQ) Stock Market [39], we observed daily changes in crude oil prices and trading volumes (in US dollars per barrel) and natural gas (in US dollars per cubic meter) from 2 April 2018 to 31 March 2023.

In this way, time series of real-world data of length $T=1261$ were obtained, and their main statistical indicators can be seen in the following Table 2. Based on these, it can be easily concluded that in both cases there are pronounced and permanent fluctuations. For example, the average price of crude oil is (approximately) 62.57 US dollars per barrel. However, it varies from 9.06 US dollars (on 21 April 2020, just over a month after the official announcement of the COVID-19 pandemic), all the way to 123.7 US dollars (on 8 March 2022, a few weeks after the start of the war in Ukraine). Let us notice that the price of natural gas also has pronounced price ranges, although less pronounced than in the case of the price of crude oil.

Note that in addition to the basic financial data (price and volumes of trading), Table 2 also shows descriptive statistics of the so-called log-volumes. They represent an aggregate financial indicator, obtained as a natural logarithm of the total monetary value of trading volumes, i.e.,

$y_t^{\left(j\right)}:=\ln \left(P_t^{\left(j\right)}·V_t^{\left(j\right)}\right),\hspace{1em}\hspace{1em}t=0,1,\dots ,T_j,\hspace{1em}\hspace{1em}j=1,2.$

Here, $\left(P_t^{\left(1\right)}\right)$ and $\left(V_t^{\left(1\right)}\right)$ are, respectively, the price and trading volumes of crude oil, and $\left(P_t^{\left(2\right)}\right)$ and $\left(V_t^{\left(2\right)}\right)$ are the price and trading volumes of natural gas, observed at some point in time $t=1,2,\dots ,T$. As is stated in two studies [40,41], the usage of log-volumes changes the interpretation of activity shocks because unexpected values are not affected by the growth trend in their dynamics. In addition, the variance of log-volatility shocks is then more uniform across the sample (that is, over the timeline of the observed data). This can also be seen through the sample variance and standard deviation of both observed log-volume series, which are shown in Table 2. Additionally, the corresponding Split-MA(1) processes for these series are as follows:

$X_t^{\left(j\right)}≔y_t^{\left(j\right)}-y_{t-1}^{\left(j\right)}=\ln \frac{P_t^{\left(j\right)}}{P_{t-1}^{\left(j\right)}}+\ln \frac{V_t^{\left(j\right)}}{V_{t-1}^{\left(j\right)}},\hspace{1em}\hspace{1em}t=1,\dots ,T,\hspace{1em}\hspace{1em}j=1,2,$

i.e., they represent the sum of the log-returns of prices and trading volumes.

We further consider the possibility of using the LSB process as a suitable stochastic model of logarithmic volume dynamics. To that end, the basic LSB series, i.e., realizations of log-volumes of crude oil and natural gas, will be referred to as Series A and Series B, respectively. According to these, as well as the results of using Equations (1) and (5), the martingale means $\left(m_t^{\left(j\right)}\right)$ and innovations $\left({\epsilon }_t^{\left(j\right)}\right)$ can be obtained by iterative procedure:

(47)$\left\{\begin{array}{l}{\epsilon }_t^{\left(j\right)}=y_t^{\left(j\right)}-m_t^{\left(j\right)},\\ m_t^{\left(j\right)}=m_{t-1}^{\left(j\right)}+{\epsilon }_{t-1}^{\left(j\right)}I\left\{{\left({\epsilon }_{t-2}^{\left(j\right)}\right)}^2\ge \tilde{c}\right\},\end{array}\right.$

where $j=1,2$ and $\tilde{c}$ is the estimated critical value, obtained by using Equation (31). As initial values in (47), as before, we have taken ${\epsilon }_0^{\left(j\right)}={\epsilon }_{-1}^{\left(j\right)}=0,$ as well as $m_0^{\left(j\right)}=y_0^{\left(j\right)}=\widehat{\mu }$, $j=1,2$. The estimated values of basic statistical indicators of the increment series $\left(X_t^{\left(j\right)}\right)$, $j=1,2$, as well as two modelled series, martingale means $\left(m_t^{\left(j\right)}\right)$, $j=1,2,$ and innovation series $\left({\epsilon }_t^{\left(j\right)}\right)$, $j=1,2$, are presented in the following Table 3.

Based on the obtained estimated parameter values, certain observations can be made, which also derive from previously obtained theoretical results. Note first that the average values of log-volumes are “close” to the average values of martingale means, and that is consistent with equality $E\left(y_t\right)=E\left(m_t\right)$. Moreover, both series A and B have other similar statistical indicators (variance or standard deviation, for example), which indicates a certain similarity in their dynamics and other stochastic characteristics. This can be seen by comparing statistical indicators of the increments ($X_t^{\left(j\right)}$), $j=1,2$, and the innovation series (${\epsilon }_t^{\left(j\right)}$), $j=1,2$. It is noticeable that the sample means of both series are close to zero, that is, they have the property of symmetry of their empirical distributions. Finally, it is worth pointing out that both series (${\epsilon }_t^{\left(j\right)}$), $j=1,2$, have the estimated values of kurtosis $K_j\approx 6$, $j=1,2,$ which could indicate their suitability for stochastic modelling with the Laplace distribution.

In the following, using the estimation procedure described in the previous section, the parameters of both real-world time series can be estimated, assuming that their dynamics are subjected to the LSB model. Table 4 shows the estimates obtained by applying the above-mentioned procedures, that is, two kinds of parameter estimates of the LSB model. Additionally, some other estimates, such as the first-order sample correlation ${\widehat{\rho }}_X\left(1\right)$ and the estimates of the threshold parameter $b_c$, are shown. Let us notice that the condition $-0.5<{\widehat{\rho }}_X\left(1\right)<0$ is satisfied in both series cases, which enables the estimation of the parameter $b_c$.

As has already been pointed out, the values of the modelled series ($m_t^{\left(j\right)}$) and (${\epsilon }_t^{\left(j\right)}$) were computed by using the most robust estimators of the LSB process $\widehat{c},\widehat{\mu },{\widehat{\sigma }}^2$. The agreement between the modelled and actual data can be seen in Figure 8a, where in addition to the observed log-volume values ($y_t^{\left(j\right)}$), the modelled martingale mean values ($m_t^{\left(j\right)}$) are shown. At the same time, the agreement between the increment series ($X_t^{\left(j\right)}$) and innovations (${\epsilon }_t^{\left(j\right)}$) is shown in Figure 8b. It is worth noting that high correlations between the actual and modelled time series are clearly observable, which can also be explained by the theoretical findings presented in Section 2. Namely, the martingale means ($m_t^{\left(j\right)}$) are equal to the log-volumes in cases when there were no pronounced fluctuations of the series ($y_t^{\left(j\right)}$) in the previous time period. On the other hand, if emphatic fluctuations occur, the values of the series ($m_t^{\left(j\right)}$) and ($y_t^{\left(j\right)}$) become different, and the resulting deviations indicate the existence of significant fluctuations and potential risk in the market. Similarly, if at some point in time point’s innovation series (${\epsilon }_t^{\left(j\right)}$) has a pronounced fluctuation, the next value of (${\epsilon }_t^{\left(j\right)}$) will be equal to increments ($X_t^{\left(j\right)}$). It is obvious that the agreement of realizations between these time series is better if, in addition to permanent and emphatic fluctuations of (${\epsilon }_t^{\left(j\right)}$), the critical value c is relatively small.

Notice that the CDFs and PDFs of the series ($X_t^{\left(j\right)}$), $j=1,2$ can be obtained by using the results given in Theorem 3, that is, Equations (10) and (13), respectively. The above plots in Figure 9 show fitted PDFs of both empirical distributions of the series ($X_t^{\left(j\right)}$), $j=1,2$. Finally, using the results of Theorem 4, that is, Equation (15), the fitted CDFs $H_t\left(x\right)$ of the log-volumes ($y_t^{\left(j\right)}$) can be obtained by the iterative procedure:

(48)$H_t^{\left(j\right)}\left(x\right)≔\left(H_{t-1}^{\left(j\right)}\otimes F_{\eta }\right)\left(x\right)=\underset{-\infty }{\overset{+\infty }{\int }}{H}_{t-1}^{\left(j\right)}\left(x-u\right)F_{\eta }\left(du\right),t=1,2,\dots ,T,$

where:

(49)$\begin{array}{ll}{H}_0^{\left(j\right)}\left(x\right)& ≔\left(F_{\mu }\otimes F_{\epsilon }\right)\left(x\right)=\underset{-\infty }{\overset{+\infty }{\int }}{F}_{\mu }\left(x-u\right)F_{\epsilon }\left(du\right)=\underset{-\infty }{\overset{+\infty }{\int }}I\left(x-u>\mu \right)f_{\epsilon }\left(u\right)du\\ & =\underset{-\infty }{\overset{x-\mu }{\int }}{f}_{\epsilon }\left(u\right)du=F_{\epsilon }\left(x-\mu \right).\end{array}$

From here, by differentiating the CDFs obtained using Equation (48), the corresponding PDFs of the log-volumes ($y_t^{\left(j\right)}$) can also be easily computed. Therefore, due to the non-stationarity, these PDFs are dependent on the time argument $t$. The graphs below in Figure 9 show the theoretical PDFs of the series ($y_t^{\left(j\right)}$), obtained by using the numerical procedure in the R-package “distr” [42]. The PDFs of the length $t<T=1261$ are shown with dashed lines, and the PDFs of the length $T=1261$ are shown with a solid line.

7. Conclusions

The main stochastic properties of the Laplacian Split-BREAK (LSB) process are presented here, along with the investigation of the asymptotic properties of the corresponding LSB series, as well as the procedures for estimating their parameters. It is useful to point out that one of the advantages of this stochastic model, as in the case of the Gaussian Split-BREAK (GSB) process, is that it enables the usage of appropriate stationary and non-stationary components, which provide different procedures for estimating its unknown parameters. At the same time, of particular importance is the asymptotic behaviour of LSB series. This is considered as well as the obtained parameter estimators.

Let us point out that one of the important features of the LSB process, as well as the class of STOPBREAK-based processes in general, is the ability to “remove” the sharp boundary between stochastic processes with permanent shocks and those in which they remain transient. Therefore, these stochastic processes can vary between different well-known non-linear stochastic models (see, for more details, Stojanović et al. [23,24]). For instance, the LSB process can vary from an IID (white noise) series, for a larger critical value of reaction $c$, to a random walk process, as $c$ approaches zero.

In addition, some of the possibilities of applying the LSB process in modelling dynamics of the real-world series with emphasized and persistent fluctuations are also described. This provides opportunities for potential future research based on the various kinds of Split-BREAK processes. At the same time, it is worth pointing out that the dynamic analysis of log-volumes, as composite time series, may represent a certain limitation, due to a possibility of omitting some other characteristics of the oil and gas market.

Author Contributions

Conceptualization, V.S.S. and H.S.B.; methodology, V.S.S. and H.S.B.; software, V.S.S.; validation, E.L. and I.B.; formal analysis, V.S.S. and H.S.B.; data curation, V.S.S. and I.B.; writing---original draft preparation, V.S.S., H.S.B., and E.L.; writing---review and editing, E.L. and I.B.; visualization, V.S.S. and H.S.B.; supervision, V.S.S. and H.S.B.; project administration, E.L. and I.B. All authors have read and agreed to the published version of the manuscript.

Data Availability Statement

Official data on prices and trading volumes on the Oil and Gas World market can be found on the website https://nasdaq.com.

Conflicts of Interest

The authors declare no conflict of interest.

Abbreviations

The following abbreviations are used in this manuscript:

AD test	Anderson–Darling test (statistics)
AN	Asymptotic Normality
AR process	Autoregressive process
CDF	Cumulative Distribution Function
CF	Characteristic Function
CLD	Contaminated Laplacian Distribution
DF	Degrees of Freedom
GSB process	Gaussian Split-BREAK process
IID	Independent Identical Distributed
LSB process	Laplacian Split-BREAK process
MA process	Moving Average process
ML method	Maximum Likelihood method
MSEE	Mean Square Estimated Error
NASDAQ	National Association of Securities Dealers Automated Quotations
PDF	Probability Density Function
RV	Random Variable
Split-MA process	Splitting Moving Average process
STOPBREAK	Stochastic Permanent Breaking
W test	Cramér–von Mises test (statistics)

Footnotes

Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

Figures and Tables

Figure 1. Histograms of the empirical distributions of log-returns fitted with Gaussian and Laplace distributions.

Figure 2. Dynamics of the non-stationary and stationary series of the LSB process (parameter values are: [Forumla omitted. See PDF.] and [Forumla omitted. See PDF.]).

Figure 3. Plots of CDFs (a) and PDFs (b) of increments [Forumla omitted. See PDF.], regarded as a mixed distribution of Laplacian distributed innovations [Forumla omitted. See PDF.], and two chi-square-based distributions (parameter values are: [Forumla omitted. See PDF.] [Forumla omitted. See PDF.], [Forumla omitted. See PDF.]).

Figure 4. Convergences of modulus of the characteristic functions [Forumla omitted. See PDF.] and [Forumla omitted. See PDF.], when [Forumla omitted. See PDF.] (parameter values are the same as in Figure 3).

Figure 5. (a) Plots of the asymptotic variances of the estimators [Forumla omitted. See PDF.] (dashed line) and [Forumla omitted. See PDF.] (solid line). (b) Plot of the ratio [Forumla omitted. See PDF.] of the asymptotic variances of the estimators [Forumla omitted. See PDF.] and [Forumla omitted. See PDF.].

Figure 6. Histograms of empirical distributions of the estimated parameters (true parameter values are the same as in Table 1).

Figure 7. Dynamics of prices and trading volumes of energy products on the world market in the period from 2018 to 2023: (a) crude oil; (b) natural gas.

Figure 8. Dynamic diagrams of empirical and modelled data: (a) log-volumes (solid lines) and martingale means (dashed lines); (b) increments (solid lines) and innovation series (dashed lines). The diagrams above present the dynamics of Series A, and below are the dynamics of Series B.

Figure 9. Empirical distributions of actual data (given by histograms), along with the corresponding fitted PDFs (given by lines): (a) Series A; (b) Series B.

Table 1

Summary statistics of the estimated parameter values of the LSB process, along with the realized statistics of their normality tests (true parameter values are: $\mu =0$ and $c=\lambda =1$).

ParameterEstimators		Statistics	Estimated Values	AD(p-Value)	W(p-Value)
Critical Value	$$\tilde{c}$$	Min.	0.5769	0.9908 *(0.0129)	0.1791 **(0.0099)
		Mean	1.0257
		(MSEE)	(3.63 × 10 −2)
		Max.	2.0964
	$$\widehat{c}$$	Min.	0.7142	0.3202(0.5320)	0.0441(0.6059)
		Mean	0.9944
		(MSEE)	(7.19 × 10 −3)
		Max.	1.2685
Scale Parameter	$$\tilde{\lambda }$$	Min.	0.7829	0.5413(0.1641)	0.0702(0.2773)
		Mean	1.0026
		(MSEE)	(5.36 × 10 −3)
		Max.	1.2394
	$$\widehat{\lambda }$$	Min.	0.8592	0.4244(0.3173)	0.0588(0.3920)
		Mean	1.0028
		(MSEE)	(2.52 × 10 −3)
		Max.	1.1702
Mean Value	$$\tilde{\mu }$$	Min.	−56.420	0.1843(0.9088)	0.0282(0.8701)
		Mean	0.3491
		(MSEE)	(252.01)
		Max.	55.208
	$$\widehat{\mu }$$	Min.	−38.595	0.2508(0.7417)	0.0289(0.8621)
		Mean	0.1647
		(MSEE)	(96.62)
		Max.	33.0875

* $p<0.05$, ** $p<0.01$.

Table 2

Basic statistical and market indicators of crude oil and natural gas in the last five years.

Statistics	Crude Oil			Natural Gas
	Price	Volumes	Log-Volumes	Price	Volumes	Log-Volumes
Mean	65.42	4.22 × 105	16.832	3.591	1.22 × 105	12.672
Median	64.62	3.88 × 105	17.102	2.855	1.17 × 105	12.789
Mode	67.04	N/A	N/A	2.662	1.66 × 105	12.717
Sample Variance	397.1	5.41 × 1010	0.8172	3.274	3.56 × 109	0.8194
Stand. Deviation	19.94	2.33 × 105	0.9040	1.809	5.97 × 104	0.9052
Minimum	9.060	1.23 × 104	11.617	1.482	1.20 × 103	8.281
Maximum	123.7	1.69 × 106	18.273	9.680	4.35 × 105	14.56

Table 3

Statistical indicators of increments, martingale means, and innovations of LSB processes.

Statistics	Series A			Series B
	$\left({\mathit{X}}_{\mathit{t}}^{\left(1\right)}\right)$	$\left({\mathit{m}}_{\mathit{t}}^{\left(1\right)}\right)$	$\left({\mathit{\epsilon }}_{\mathit{t}}^{\left(1\right)}\right)$	$\left({\mathit{X}}_{\mathit{t}}^{\left(2\right)}\right)$	$\left({\mathit{m}}_{\mathit{t}}^{\left(2\right)}\right)$	$\left({\mathit{\epsilon }}_{\mathit{t}}^{\left(2\right)}\right)$
Mean	3.52 × 10−5	16.898	−0.0196	1.17 × 10−5	12.688	−0.0713
Median	−0.0577	17.127	−0.0537	−0.0660	12.830	−0.0702
Mode	N/A	16.588	N/A	N/A	12.709	N/A
Sample Variance	0.9688	0.7227	0.7156	0.8492	0.8975	0.6484
Stand. Deviation	0.8178	0.8501	0.8460	0.9215	0.9473	0.8053
Minimum	−3.5600	11.617	−3.5560	−4.0362	8.2808	−3.8895
Maximum	4.9542	18.185	4.9542	4.6061	14.560	4.2371
Range	9.5142	6.5680	8.5142	8.6423	6.2794	8.1266
Skewness	2.3966	−2.1680	1.0276	1.5001	−1.9954	0.6335
Kurtosis	5.1593	5.3253	6.1547	7.7933	5.2876	6.0679

Table 4

Estimated parameter values of the log-volumes series.

Parameter Estimates		Series A	Series B
Mean Value	$\tilde{\mu }$	16.832	12.672
	$\widehat{\mu }$	17.004	12.587
Sample Correlation	${\widehat{\rho }}_X\left(1\right)$	−0.1911	−0.3948
Threshold Parameter	${\tilde{b}}_c$	0.2362	0.6523
	${\widehat{b}}_c$	0.3766	0.5546
Critical Value	$\tilde{c}$	0.0196	0.2868
	$\widehat{c}$	0.0459	0.1487
Scale Parameter	$\tilde{\lambda }$	0.5201	0.5069
	$\widehat{\lambda }$	0.4520	0.5113

References

1. Hardouvelis, G.A.; Kim, D. Margin Requirements, Price Fluctuations, and Market Participation in Metal Futures. J. Money Credit Bank.; 1995; 27, pp. 659-671. [DOI: https://dx.doi.org/10.2307/2077741]

2. Fatás, A. Endogenous Growth and Stochastic Trends. J. Monet. Econ.; 2000; 45, pp. 107-128. [DOI: https://dx.doi.org/10.1016/S0304-3932(99)00043-4]

3. Chang, Y.; Schorfheide, F. Labor-Supply Shifts and Economic Fluctuations. J. Monet. Econ.; 2003; 50, pp. 1751-1768. [DOI: https://dx.doi.org/10.1016/j.jmoneco.2003.02.001]

4. Artis, M.J.; Hoffmann, M. Financial Globalization, International Business Cycles and Consumption Risk Sharing. Scand. J. Econ.; 2008; 110, pp. 447-471. [DOI: https://dx.doi.org/10.1111/j.1467-9442.2008.00546.x]

5. Lamey, L.; Deleersnyder, B.; Steenkamp, J.B.E.; Dekimpe, M.G. The Effect of Business-Cycle Fluctuations on Private-Label Share: What Has Marketing Conduct Got to Do with It?. J. Mark.; 2012; 76, pp. 1-19. [DOI: https://dx.doi.org/10.1509/jm.09.0320]

6. Blanchard, O.J.; L’Huillier, J.P.; Lorenzoni, G. News, Noise, and Fluctuations: An Empirical Exploration. Am. Econ. Rev.; 2013; 103, pp. 3045-3070. [DOI: https://dx.doi.org/10.1257/aer.103.7.3045]

7. Rodríguez, G.; Villanueva Vega, P.; Castillo Bardalez, P. Driving Economic Fluctuations in Peru: The Role of the Terms of Trade. Emp. Econ.; 2018; 55, pp. 1089-1119. [DOI: https://dx.doi.org/10.1007/s00181-017-1318-2]

8. Rebei, N.; Sbia, R. Transitory and Permanent Shocks in the Global Market for Crude Oil. J. Appl. Econom.; 2021; 36, pp. 1047-1064. [DOI: https://dx.doi.org/10.1002/jae.2863]

9. Xu, Z.; Wang, H.; Zhang, H.; Zhao, K.; Gao, H.; Zhu, Q. Non-Stationary Turbulent Wind Field Simulation of Long-Span Bridges Using the Updated Non-Negative Matrix Factorization-Based Spectral Representation Method. Appl. Sci.; 2019; 9, 5506. [DOI: https://dx.doi.org/10.3390/app9245506]

10. Granero-Belinchón, C.; Roux, S.G.; Garnier, N.B. Information Theory for Non-Stationary Processes with Stationary Increments. Entropy; 2019; 21, 1223. [DOI: https://dx.doi.org/10.3390/e21121223]

11. Zhao, D.; Gelman, L.; Chu, F.; Ball, A. Novel Method for Vibration Sensor-Based Instantaneous Defect Frequency Estimation for Rolling Bearings under Non-Stationary Conditions. Sensors; 2020; 20, 5201. [DOI: https://dx.doi.org/10.3390/s20185201] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/32933089]

12. Qu, C.; Li, J.; Yan, L.; Yan, P.; Cheng, F.; Lu, D. Non-Stationary Flood Frequency Analysis Using Cubic B-Spline-Based GAMLSS Model. Water; 2020; 12, 1867. [DOI: https://dx.doi.org/10.3390/w12071867]

13. Aguejdad, R. The Influence of the Calibration Interval on Simulating Non-Stationary Urban Growth Dynamic Using CA-Markov Model. Remote Sens.; 2021; 13, 468. [DOI: https://dx.doi.org/10.3390/rs13030468]

14. Narr, C.F.; Chernyavskiy, P.; Collins, S.M. Partitioning Macroscale and Microscale Ecological Processes Using Covariate-Driven Non-Stationary Spatial Models. Ecological Applications; 2022; 32, e02485. [DOI: https://dx.doi.org/10.1002/eap.2485] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/34676934]

15. Engle, R.F.; Smith, A.D. Stochastic Permanent Breaks. Rev. Econ. Stat.; 1999; 81, pp. 553-574. [DOI: https://dx.doi.org/10.1162/003465399558382]

16. Diebold, F.X.; Inoue, A. Long Memory and Regime Switching. J. Econom.; 2001; 105, pp. 131-159. [DOI: https://dx.doi.org/10.1016/S0304-4076(01)00073-2]

17. Gonzalo, J.; Martínez, O. Large Shocks vs. Small Shocks. (Or does size matter? May be so.). J. Econom.; 2006; 135, pp. 311-347. [DOI: https://dx.doi.org/10.1016/j.jeconom.2005.07.022]

18. Kapetanios, G.; Tzavalis, E. Modeling Structural Breaks in Economic Relationships Using Large Shocks. J. Econ. Dynam. Control; 2010; 34, pp. 417-436. [DOI: https://dx.doi.org/10.1016/j.jedc.2009.10.001]

19. Dendramis, Y.; Kapetanios, G.; Tzavalis, E. Level Shifts in Stock Returns Driven by Large Shocks. J. Empir. Financ.; 2014; 29, pp. 41-51. [DOI: https://dx.doi.org/10.1016/j.jempfin.2014.04.001]

20. Dendramis, Y.; Kapetanios, G.; Tzavalis, E. Shifts in Volatility Driven by Large Stock Market Shocks. J. Econ. Dynam. Control; 2015; 55, pp. 130-147. [DOI: https://dx.doi.org/10.1016/j.jedc.2015.03.006]

21. Stojanović, V.; Popović, B.Č.; Popović, P. The Split-BREAK Model. Braz. J. Probab. Stat.; 2011; 25, pp. 44-63. [DOI: https://dx.doi.org/10.1214/09-BJPS025]

22. Stojanović, V.; Popović, B.Č.; Popović, P. Stochastic Analysis of GSB Process. Publ. Inst. Math.; 2014; 95, pp. 149-159. [DOI: https://dx.doi.org/10.2298/PIM1409149S]

23. Stojanović, V.; Popović, B.Č.; Popović, P. Model of General Split-BREAK Process. REVSTAT– Stat. J.; 2015; 13, pp. 145-168.

24. Stojanović, V.; Milovanović, G.V.; Jelić, G. Distributional Properties and Parameters Estimation of GSB Process: An Approach Based on Characteristic Functions. ALEA –Lat. Am. J. Probab. Math. Stat.; 2016; 13, pp. 835-861. [DOI: https://dx.doi.org/10.30757/ALEA.v13-33]

25. Jovanović, M.; Stojanović, V.; Kuk, K.; Popović, B.; Čisar, P. Asymptotic Properties and Application of GSB Process: A Case Study of the COVID-19 Dynamics in Serbia. Mathematics; 2022; 10, 3849. [DOI: https://dx.doi.org/10.3390/math10203849]

26. Johnson, N.L.; Kotz, S.; Balakrishnan, N. Continuous Univariate Distributions; Wiley: Hoboken, NJ, USA, 1994; ISBN 0-471-58495-9

27. Williams, D. Probability with Martingales; Cambridge University Press: Cambridge, UK, 1991; Section 18.1.

28. Stojanović, V.; Popović, B.Č.; Milovanović, G.V. The Split-SV model. Comput. Statist. Data Anal.; 2016; 100, pp. 560-581. [DOI: https://dx.doi.org/10.1016/j.csda.2014.08.010]

29. Stojanović, V.; Kevkić, T.; Jelić, G. Application of the Homotopy Analysis Method in Approximation of Convolutions Stochastic Distributions. Univ. Politehnica Bucharest Sci. Bull.; 2017; 79, pp. 103-112.

30. Varadhan, S.R.S. Probability Theory, in: Courant Lecture Notes; American Mathematical Society: Providence, RI, USA, 2001; 131.

31. Fuller, W.A. Introduction to Statistical Time Series; John Wiley & Sons: New York, NY, USA, 1996.

32. Hoeffding, W.; Robbins, H. The central limit theorem for dependent random variables. Duke Math J.; 1948; 15, pp. 773-780.

33. Serfling, R.J. Approximation Theorems of Mathematical Statistics; 2nd ed. John Wiley & Sons: New York, NY, USA, 2002.

34. Lawrence, A.J.; Lewis, P.A.W. Reversed Residuals in Autoregressive Time Series Analysis. J. Time Ser. Anal.; 1992; 13, pp. 253-266. [DOI: https://dx.doi.org/10.1111/j.1467-9892.1992.tb00105.x]

35. Norton, R.M. The Double Exponential Distribution: Using Calculus to Find a Maximum Likelihood Estimator. Am. Stat.; 1984; 38, pp. 135-136. [DOI: https://dx.doi.org/10.2307/2683252]

36. Sharma, D. On Estimating the Variance of a Generalized Laplace Distribution. Metrika; 1984; 31, pp. 85-88. [DOI: https://dx.doi.org/10.1007/BF01915188]

37. Billingsley, P. The Lindeberg–Levy Theorem for Martingales. Proc. Am. Math. Soc.; 1961; 12, pp. 788-792.

38. Tests for Normality. R Package Version 1.0-2. 2013; Available online: http://CRAN.R-project.org/package=nortest (accessed on 31 March 2023).

39. Nasdaq. Available online: https://nasdaq.com (accessed on 3 April 2023).

40. Bessembinder, H.; Seguin, P.J. Futures-Trading Activity and Stock Price Volatility. J. Finan.; 1992; 47, pp. 2015-2034. [DOI: https://dx.doi.org/10.1111/j.1540-6261.1992.tb04695.x]

41. So, M.K.; Chen, C.W.; Chiang, T.C.; Lin, D.S. Modelling Financial Time Series with Threshold Nonlinearity in Returns and Trading Volume. Appl. Stoch. Model. Bus. Ind.; 2007; 23, pp. 319-338. [DOI: https://dx.doi.org/10.1002/asmb.674]

42. Ruckdeschel, P.; Kohl, M.; Stabla, T.; Camphausen, F. S4 Classes for Distributions. R News; 2006; 6, pp. 2-6. Available online: https://CRAN.R-project.org/doc/Rnews (accessed on 31 March 2023).