key: cord-0581012-sb8zq1o7 authors: Bhootna, Niharika; Kumar, Arun title: Tempered Stable Autoregressive Models date: 2021-06-06 journal: nan DOI: nan sha: 3b258ca35707a036ebe4f2f24cf3d5eae1abe23e doc_id: 581012 cord_uid: sb8zq1o7 In this article, we introduce and study a one sided tempered stable first order autoregressive model called TAR(1). Under the assumption of stationarity of the model, the marginal probability density function of the error term is found. It is shown that the distribution of the error term is infinitely divisible. Parameter estimation of the introduced TAR(1) process is done by adopting the conditional least square and method of moments based approach and the performance of the proposed methods are evaluated on simulated data. Also we study an autoregressive model of order one with tempered stable innovations. Using appropriate test statistic it is shown that the model fit very well on real and simulated data. Our models generalize the inverse Gaussian and one-sided stable autoregressive models existing in the literature. Autoregressive models serve as eminent approach of modelling among different time-series methods. Many real life time-series data exhibits autoregressive behaviour. Classical autoregressive time series models are based on the assumption of normality of the innovation term. However, this assumption may not be true for all the real life scenarios. Many real life time series probability distributions display heavy-tailed behaviour or semi-heavy tailed behaviour. Semi-heavy tailed pdf are those where tails are heavier than the Gaussian and lighter than the power law, see e.g. Omeya et al. [15] . It is well known that series of counts, proportions, binary outcomes or non-negative observations are some examples of non-normal real life time-series data see e.g. Grunwald et al. [8] . The autoregressive model of order 1 denoted by AR(1) is a simple, useful and interpretable model in a wide range of real life applications. Several AR(1) models with different marginal distributions are considered in literature. For example, Gaver and Lewis [7] considered linear additive first-order autoregressive scheme with gamma distributed marginals. The authors show that error distribution is same as the distribution of a non-negative random variable which has a point mass at 0 and which is exponential if positive (see Gaver and Lewis [7] ). A more general exponential ARMA(p, q) model called EARMA(p, q) is considered in Lawrance and Lewis [13] . A second order autoregressive model is considered in Dewald and Lewis [4] , where Laplace distributed marginals are considered. Abraham and Balakrishna [1] considered inverse Gaussian autoregressive model where the marginals are inverse Gaussian distributed. The authors shows that after fixing a parameter of the inverse Gaussian to 0 the error term is also inverse Gaussian distributed. In this paper, we consider tempered stable autoregressive model of order 1 and denote it by TAR (1) . This model generalizes the work discussed in Abraham and Balakrishna [1] . Since one parameter inverse Gaussian (also called Lévy distribution) is a particular case of the one sided tempered stable distribution. We derive the explicit form of the error density. Further, we show that if the AR(1) series marginals are stable distributed the error is also stable distributed with some scaling. A step-by-step procedure of the estimation of the parameters of the proposed model is given. Moreover, we consider an AR(1) model with one sided tempered stable innovations. This model is different from the TAR(1) model where marginals are one sided tempered stable and here the innovations are positive tempered stable distributed. The application of this model is shown on real life power consumption data. The rest of the paper is organized as follows: Section 2 defines the TAR(1) model where the explicit form of the density of the error term , infinite divisibility, moments and parameter estimation methods of TAR(1) model are discussed. In Section 3, an AR(1) model with positive tempered stable innovations is introduced. The estimation procedure of the parameters of the introduced model is given based on method of moments and conditional least square. The simulations study, where efficacy of the estimation procedure is discussed based on simulated data and real life application are discussed in this section. The last section concludes. In this section, we introduce the tempered stable autoregressive model and discuss the main properties. Consider an autoregressive process of order 1, defined as: where |ρ| < 1 and {ǫ n , n ≥ 1} is sequence of i.i.d. random variables. Assuming marginals X n 's to be stationary positive tempered stable distributed. Then there exist a distribution of ǫ n 's. The class of stable distributions is denoted by Stable(β, ν, µ, σ), with parameter β ∈ (0, 2] is the stability index, ν ∈ [−1, 1] is the skewness parameter, µ ∈ R is the location parameter and σ > 0 is the shape parameter. The stable class probability density functions do not possess closed form except for three cases (Gaussian, Cauchy, and Lévy). Generally stable distributions are represented in terms of their characteristic functions or Laplace transforms. Stable distributions are infinitely divisible. The onesided stable random variable S has following Laplace transform (see Samorodnitsky and Taqqu [19] ) E(e −sS ) = e −s β , s > 0, β ∈ (0, 1). The right tail of the S behaves (see e.g. Samorodnitsky and Taqqu [19] ) , as x → ∞. Next, we introduce the positive tempered stable distribution. The positive tempered stable random variable T with tempering parameter λ > 0 and stability index β ∈ (0, 1), has the Laplace transform (LT) (see e.g. Meerschaert et al. [14] , Kumar and Vellaisamy [12] ) Note that tempered stable distributions are obtained by exponential tempering in the distributions of α-stable distributions (see Rosinski [18] ). The advantage of tempered stable distribution over an α-stable distribution is that it has moments of all order and its density is also infinitely divisible. The probability density function for T is given by where f β (x) is the PDF of an α-stable distribution (Uchaikin and Zolotarev [22] ). The tail probability of tempered stable distribution has the following asymptotic behavior where c β,λ = 1 βπ Γ(1 + β) sin(πβ)e λ β . The first two moments of tempered stable distribution are given by (2.7) A pdf u(x) is called a semi-heavy tailed pdf if where v is a regularly varying function (see e.g. Omeya et al. [15] ). Recall that a positive function v is regularly varying with index ν if Using (2.6), it is straight forward to show that and hence the tempered stable density function is semi-heavy tailed. In this subsection, we discuss about the distributional properties related to the introduced TAR(1) process. If each X n in (2.1) is positive tempered stable then the Laplace transform of X n is given by Φ Xn (s) = E(exp(−sX n )) = E(exp(−s(ρX n−1 + ǫ n ))) = Φ X n−1 (s)Φ ǫn (s). (2.8) Assuming X n 's to be stationary distributed, it follows Φ X (s) = exp{−((λ + s) β − λ β )}, λ > 0, β ∈ (0, 1), (2.9) and Φ X (ρs) = exp{−((λ + ρs) β − λ β )}, λ > 0, β ∈ (0, 1). (2.10) Putting values from (2.9) and (2.10), yields Proof. By Feller [5] , a pdf f (x), x ≥ 0 is infinitely divisible iff its Laplace transform is of the form F (s) = e −φ(s) , s > 0 where φ(s) has a completely monotone derivative. Further, a function ψ(·) is completely monotone if it possesses derivatives ψ (n) (·) of all orders and (−1) n ψ (n) (s) ≥ 0, s > 0, n = 0, 1, · · · (see e.g. Feller [5] , p. 439). For Φ ǫ (s), we have φ(s) = (s + λ) β − (ρs + λ) β which gives φ ′ (s) = β(s + λ) β−1 − βρ(ρs+λ) β−1 . Consider ψ(s) = (s+λ) β−1 , which yields to ψ (n) (s) = (β−1)(β−2) · · · (β−n)(s+λ) β−n−1 . Thus (−1) n ψ (n) (s) ≥ 0, n = 0, 1, 2, · · · . Similarly, considering ζ(s) = (ρs + λ) β−1 we can show that ζ(s) is completely monotone. Moreover, (−1) n φ n (s) = (−1) n ψ n (s) − (−1) n ζ n (s) and ψ n (s) is always greater than ζ n (s) as 0 < ρ < 1, this implies (−1) n φ n (s) ≥ 0. So φ(s) has completely monotone derivative. Next, we provide the explicit form of the pdf of the error term using complex inversion formula, which generalize the corresponding results in Abraham and Balakrishna [1] . the pdf of the innovation term has following integral form where λ > 0 and β ∈ (0, 1). Proof. Inverse Laplace transform corresponding to Φ ǫ (s) provides the density of innovation term. The density for innovation can be computed using complex inversion formula for Laplace transform given by where point x 0 > a for some a is taken in such a way that the integrand is analytic for Re(s) > a. Note that the integrand e sx Φ ǫ (s) is an exponential function which is analytic in whole complex plane. However, due to fractional power in the exponent the integrand e sx Φ ǫ (s) has branch points at s = −λ and s = −λ/ρ. Thus we take a branch cut along the non-positive real axis and consider a (singlevalued) analytic branch of the integrand. To calculate (2.16), consider a closed double-key-hole contour C: ABCDEFGHIJA ( Fig. 1) with branch points at O = (−λ, 0) and O ′ = (−λ/ρ, 0). The contour consists of following segments: AB and IJ are arcs of radius R, BC, DE, FG and HI are line segment parallel to x-axis, CD, EF and GH are arcs of circles with radius r and JA is the line segment from x 0 − iy to x 0 + iy (see Fig. 1 ). The integrand is analytic within and on the contour C so that by Cauchy's residue theorem (see e.g. Schiff [21] ) Using (Lemma 4.1 Schiff [21] , p. 154), the integral on circular arcs AB and IJ tend to 0 as R → ∞. The integral over CD, EF and GH are also zero as r goes to infinity. Thus, as r → 0 and R → ∞, we have 1 2πi and is called Lévy probability density function see Applebaum [2] , which is a special case of the inverse Gaussian density. Abraham and Balakrishna [1] show that if the inverse Gaussian autoregressive model is stationary the error term is also inverse Gaussian distributed. Our results complement their findings. Remark 2.3. For λ = 0, the error density function in (2.15) reduces to Proposition 2.2. For λ = 0 the innovation term has the following fractional qth order moment Proof. For a positive random variable X with a Laplace transform φ(s), and for p > 0, it follows Thus for q ∈ (n − 1, n), for an integer n, we have Using the relationship for q ∈ (0, 1), it follows which gives the qth order fractional moment of the innovation term. For λ = 0, it can be written in the explicit form as follows Proof. For the given TAR(1) process X n = ρX n−1 + ǫ n , multiplying both the side by X n−r , we get X n X n−r = ρX n−1 X n−r + ǫ n X n−r . (2.27) Here X n is independent of ǫ j for j > n. Taking expectation on both the sides, yields E(X n X n−r ) = ρE(X n−1 X n−r ) + E(ǫ n X n−r ), which gives Thus, . Dividing both side by Var(X n )Var(X n−r ), yields Cov(X n X n−r ) Var(X n−1 )Var(X n−r ) , Using, Var(X n )=Var(X n−1 )=Var(X n−r ), leads to since ρ 0 = 1 and ρ 1 = ρ. Here, we discuss the estimation procedure of the parameters of the tempered stable autoregressive model. We consider conditional least square and moments based estimation methods. Parameter estimation by conditional least square: The conditional least square method provides a straightforward procedure to estimate the parameters for dependent observations by minimizing the sum of square of deviations about conditional expectation, see e.g. Klimko and Nelson [10] . To estimate the parameter ρ of process (2.1) the conditional least square (CLS) method is used. The CLS estimator of parameter vector θ = (ρ, λ, β) is obtained by minimizing with respect to θ. For an AR(1) sequence due to Markovian property Minimising (2.32) w.r.t ρ, β, and λ gives the following estimate for ρ and an estimating equation for β and λ.ρ and The estimateρ can be easily calculated from the observed sample using (2.33). However, some numerical method which minimizes the squared difference between left and right hand sides is required to estimate β and λ from (2.34). The parameters of positive tempered stable distribution cannot be estimated by classical maximum likelihood estimation method due to non availability of closed from probability density function. Also the CLS estimator provides an estimate for ρ only. Here we use method of moments to estimate the parameters β and λ. The general idea used by Pearson [16] is, the moments of the PDF about the origin are equal to the corresponding moments of the sample data. In case of TAR(1) process moments for innovation distribution can be extracted using the corresponding Laplace transform as the density (2.15) is in complex form. The idea is to apply the method of moments on the error sample which can be easily obtained from the original series by first estimating the ρ using CLS method. The kth order moment can be extracted from taking kth derivative of Laplace transform as follows for the innovation term. Using the Laplace transform defined in (2.13), we have The sample estimates for E(ǫ i ) and E(ǫ 2 i ) are given by After some manipulation, the above two equations leads to following nonlinear equations We use sequential least square optimization technique (SLSQP) to solve these equations which is is an iterative method to optimize a nonlinear model using the nonlinear least square model and sequential quadratic programming model (see Fu et al. [6] ). We use scipy.optimize which use SLSQP as an inbuilt method in python to solve nonlinear equations that gives the estimates for β and λ. Parameter estimation by tail index: The parameter estimation for the case λ = 0 (corresponding density given by (2.24)) is done by using the right tail behaviour of S for the TAR(1) process. For the TAR(1) process we have a positive stable random variable (1 − ρ β ) 1/β S with corresponding Laplace transform e −(1−ρ β )s β . Using (2.3), it follows Taking log on both sides of (2.40) and replacing P(S > (1 − ρ β ) −1/β x) = G ǫ (x), where G ǫ (x) represents the empirical tail for ǫ, we get Considering Y = log(G ǫ (x)) and X = − log(x), which yields β as slope of the regressing Y on X. To check the efficacy of the estimation procedure we perform simulation. Here we have considered an AR(1) model with marginals to be distributed as one sided tempered stable distribution and the corresponding Laplace transform for error term ǫ n is given in (2.13) . In order to simulate a tempered stable AR(1) series, first we simulate data for innovation using the Laplace transform, see e.g. Ridout (2008) . Ridout [17] discussed the simulation of a continuous distribution function using its Laplace transform if the density of the distribution is not in closed form. We consider three different simulated datasets of length n = 10, 000 each with three different parameter combinations. For the first case consider, ρ = 0.9, λ = 1 and β = 0. The sample trajectory for an AR(1) process for positive tempered stable distribution with parameters ρ = 0.9, λ = 2 and for β = 0.5 is given in Fig. 2 . Next we use resampling method to generate several samples of the simulated data. To ensure the stationarity here we assume AR(1) process with parameter ρ < |1|, as the AR(1) process is stationary if |ρ| < 1. Consider ρ = 0.9 and take tempered stable distribution parameters as β = 0.5 and λ = 2, by iterative resampling of simulated data we generate 1000 simulated samples for the defined parameters. The parameters are estimated using method of moments for λ, β and based on conditional least square for ρ considering different samples. The corresponding box-plots for estimated parameter values are shown in Fig. 3 . For the particular case, λ = 0 which leads to one-sided positive stable autoregressive model, simulation is done by the method given in Ridout [17] and parameter estimation is done by the tail index method using (2.41). We simulate a sample of 10, 000 observations for ρ = 0.5 and β = 0.5. By taking values greater than a particular threshold such that we have some percentage of the sorted dataset greater than that threshold value, we get regression based estimator for β. Taking a cutoff value of 0.009 gives the 30% of observation greater than the cutoff and gives the corresponding estimate of 0.36 for β. For cutoff value 0.019, we get almost 20% of the observation greater than threshold value and corresponding estimate is 0.41 for β which is more close to true parameter value. Further, for k = 0.06 the estimated β is 0.49 which is even more closer to the actual value. This clearly indicates that the accuracy of the estimated parameter depend on the threshold selected. The results are summarized in Table 2 , where k represents the threshold value. Threshold Actual β Estimated β Consider the stationary AR(1) process X n = ρX n−1 + ǫ n , |ρ| < 1, (3.42) where ǫ n 's are i.i.d positive tempered stable random variables defined in (2.4) and independent of X n−1 . Let φ Xn (s) and φ ǫn (s) are Laplace transform of X n and ǫ n respectively. Now assuming X 0 d = ǫ 0 and recursively writing equation (3.42), we get a moving average representation of AR(1) process as follows The Laplace transform of X n is Proposition 3.1. For λ = 0, the qth order moment of X n is given by Proof. Using the proposition (2.2), for a positive random variable X with corresponding Laplace transform φ(s), and for p > 0, it follows Thus for q ∈ (m − 1, m), for an integer m, we have For q ∈ (0, 1), it follows Remark 3.1. Using same argument as in Prop. 2.1, we show that the marginal distributions of X n are infinitely divisible. It is easy to show that for e nλ β − n i=0 (λ+ρ i s) β , the exponent φ(s) = n i=0 (λ + ρ i s) β − nλ β has completely monotone derivative, that is (−1) n φ n (s) ≥ 0, n = 0, 1, 2, · · · . Remark 3.2. The index of dispersion for X n is ID(X n ) = Var(X n ) E(X n ) . Using the mean and variance of the process X n , given by E( where 0 < ρ < 1, β ∈ (0, 1) and λ > 0. For λ > 1 the ID(X n ) < 1, which implies that the marginal distributions of X n are under-dispersed for λ > 1. To estimate the parameter ρ of the process (3.42), the conditional least square method discussed in previous section is used, which minimize with respect to θ. Using the equations (2.28), (2.30) and mean of ǫ n , we can write Minimizing (3.47) w.r.t ρ gives the following estimatê (3.48) The parameter estimation for β and λ is done again using method of moments for ǫ n , which leads to where m 1 and m 2 are first and second order sample moments for ǫ n respectively. After some manipulation, we get the following nonlinear equations . Using the optimization technique SLSQP, discussed in Section 3, these equations are solved in python to obtain the estimatesλ andβ. The estimate of λ can also be obtained using the estimateβ, denoted byλ To check the performance of the estimation method, we use simulated dataset. Using the method discussed in Ridout [17] , we simulate an independent sample of size 10,000 for error term using the Laplace transform of positive tempered stable random variable defined in (2.4). Then we recursively simulate AR(1) series by using fixed value of ρ and letting ǫ 0 = X 0 . The performance of the estimation method for three different sets of parameters is summarized in Table 3 . The sample trajectory for AR(1) process with tempered stable innovations is given in Fig. 4 . Further, to assess the performance of the estimation method, box plots of the estimated parametersρ,β andλ are plotted by taking multiple simulated random samples of innovations and generating multiple AR (1) By simply looking at the data we can visualize the stationarity with some sharp spikes suggesting non-Gaussian behaviour. We have applied Augmented Dickey-Fuller (ADF) test to check the stationarity of the data and the p-value comes out to be less than 0.05 which indicates that the null hypothesis of ADF test is rejected at 5% significance level and implies that the data is stationary. To determine the appropriate time series model, ACF and PACF plots are given in Fig. 7 and Fig. 8 respectively which determines the appropriate time series model with significant lag. Clearly PACF plot for the dataset is significant till lag 1 and the ACF plots tails off, this means that AR(1) model would be a good fit for this data. We assume that the innovation term ǫ n follows one sided tempered stable distribution. We fit the proposed AR(1) model with tempered stable distribution on the dataset. The estimated parameters using the above described methods areβ = 0.91,λ = 2.9 and ρ = 0.60. Using these estimated values, we simulate a dataset to assess the accuracy of the results. We simulate a data sample forλ = 2.8 andβ = 0.92. Considering the relation ǫ n = X n − ρX n−1 generates a dataset for the innovation term. To check whether both the simulated data using estimated parameters values and the actual data set follow same distribution or not we use two sample Kolmogorov-Smirnov (K-S) test and Mann-Whitney U test. The two sample K-S and Mann-Whitney U test both are non parametric tests which compare the distribution of two datasets. The two sample K-S test gives a p-value 0.056 and Mann-Whitney U test gives the p-value 0.377. Thus the null hypothesis for both the tests are accepted at 5% level of significance. This indicates that both the samples are from same distributions. Hence it is appropriate to model the considered power consumption time-series using AR(1) model with positive tempered stable innovations. We introduced two kinds of tempered stable autoregressive models. The first model denoted by TAR(1) has one-sided tempered stable marginals. The second model discussed in Section 3 has one-sided tempered stable innovations. Autoregressive models with positively distributed innovations are quite common in literature for example inverse Gaussian, Gamma and one-sided stable are used. The autoregressive model with one-sided tempered stable innovations generalizes the inverse Gaussian and positive stable innovations based autoregressive models. We provide the parameter estimation method of the introduced model and assess the efficacy using simulated data. Further, we show the applicability of the introduced model on real life power consumption data. We hope that introduced models will be helpful in modeling of several real life series emerging from different areas like finance, economics, natural hazards, power consumption and others. Inverse Gaussian Autoregressive Models Lévy Processes and Stochastic Calculus Financial Modeling with Jump Processes A new Laplace second-order autoregressive time series model-NLAR(2) Introduction to Probability Theory and its Applications Sequential Quadratic Programming Method for Nonlinear Least Squares Estimation and Its Application First order autoregressive gamma sequences and point processes A unified view of linear AR(1) models, Clayton Applications of the empirical characteristic function to estimation and detection problems On Conditional Least Squares Estimation for Stochastic Processes Modified Yule-Walker method for α-stable time series models Inverse tempered stable subordinator The exponential autoregressive-moving average process EARMA(p, q) Transient anamolous subdiffusions on bounded domains Semi-heavy tails Contributions to the mathematical theory of evolution Generating random numbers from a distribution specified by its Laplace transform Tempering stable processes Stable Non-Gaussian Random Processes Lévy processes and infinitely divisible distributions The Laplace Transform: Theory and Applications