A Bayesian analysis of generalized threshold autoregressive models

STATISTIC~S& ELSEVIER

Statistics & Probability Letters 40 (1998) 15-22

A Bayesian analysis of generalized threshold autoregressive models Cathy W.S. Chen* Department of Statistics, Feng-Chia University, Taichung 407, Taiwan

Abstract

The threshold autoregressive (TAR) model is generalized which results in more flexibility in applications. We construct a Bayesian framework to show that Markov chain Monte Carlo method can be applied to estimating parameters with success. © 1998 Elsevier Science B.V. All rights reserved Keywords: M C M C method; Gibbs sampler; Metropolis algorithm; Generalized T A R models

1. Introduction

Numerous nonlinear time series models have been proposed during the last three decades, among which threshold autoregressive (TAR) model and bilinear model are perhaps the most popular in the literature. The TAR model was introduced by Tong (1978) and Tong and Lim (1980) as an alternative model for describing periodic time series, which has been further developed by Tong (1983, 1990) and Tsay (1989). Some interesting statistical and physical properties of this class of models include limit cycles, amplitude dependent frequencies and jump phenomena which linear models usually fail to capture. McCulloch and Tsay (1993) proposed a Bayesian procedure for detecting threshold values in a TAR model via some posterior probability plots. Chen and Lee (1995) introduced a Bayesian analysis of threshold autoregression models by transforming a TAR model into a changepoint problem in linear regression via arranged autoregression. Broemeling and Cook (1992) and Geweke and Terui (1993) also proposed Bayesian analysis of two-regime threshold autoregressions. The above-mentioned work on univariate series {Yt} is limited to the TAR model which is a piecewise linear model in the space of Yt-d where d is the delay parameter. However, the series is often affected by other variables in practice. Hence, it is appropriate to add certain relevant exogeneous variables and to allow the threshold variable to be a function of the exogenous variables. Therefore, the main objective of this paper is to consider the above generalized threshold autoregressive (GTAR) model. The framework presented here is motivated by Gibbs sampling with Metropolis algorithm scheme. Detailed discussions can be found elsewhere (Geman and Geman, 1984; Tanner and Wong, 1987; Gelfand and Smith, * Fax: +886 4451 7092; e-mail: [email protected]. 0167-7152/98/$ - see front matter (~) 1998 Elsevier Science B.V. All rights reserved PII S 0 1 6 7 - 7 1 5 2 ( 9 8 ) 0 0 0 7 7 - 7

C W.S. Chen I Statistics & Probability Letters 40 (1998) 15-22

16

1990; and Gelfand et al., 1990). The method is extremely useful for exploring intractable joint distributions that have convenient conditional structure. The paper is organized as follows. In Section 2, we give the general framework of the models considered in this paper and show that many time series models considered in the literature are special eases of the GTAR model. In Section 3, the necessary conditional distributions for implementing the Gibbs sampler are given. Finally, we illustrate the proposed approach by using nominal GNP and money supply data from Taiwan in Section 4.

2. The generalized threshold autoregressive model We now consider a two-regime generalized threshold autoregressive (GTAR) model. Note that we assume the threshold variable to be xl without loss of generality. Let

yt =

p~ k ~""x ~( 1 ) + Z °O) ¢~(1) ..[_2...,t q)i Yt-i Pi Xi, t + al l) i=1 i=1 P2 k a (2) ~(02)"[- Z ~12)yt-i+ Z fl}2)xi, t i=1 i=1

if Xl,t-d <~r,

(I) if XI,t-- d > r,

where r is a real number and d is a positive integer commonly referred to as the delay (or threshold lag) of the model. The innovational series {a~i)} is a sequence of independently and identically distributed normal random variates with mean zero and variance a 2, where {al l)} and {al z)} are independent. Finally, { y t - l . . . . . y t - p , } is a set of lagged variables in regime i, and {xl,t . . . . . xk, t} denotes a set of exogenous variables in regime i. Model (1) is sufficiently flexible to accommodate the following interesting cases of practical importance. 1. Model (1) reduces to a TAR (threshold autoregressive) model if the exogenous variables {xi, t, i = 1 . . . . . k } are deleted. 2. Model (1) becomes a nonhomogeneous linear regression model if only the noise variance a/2 is different for different regimes. 3. Model (1) reduces to a random-level shift model if only the constant terms ¢k~i) are distinct for different regimes. In addition to the constant terms, slopes fl~) and fl~.2) may also be different for different regime. 4. Model (1) reduces to a switching regression or regime model (see Deutsch et al., 1994; Quandt, 1958; Goldfeld and Quandt, 1972; Granger and Teriisvirta, 1993) if the lagged variables {yj, t, j = 1 . . . . . pi; i = 1,2} are deleted. The main idea of our modeling procedure is to transform a GTAR model into a regular changepoint problem in linear regression analysis. This is achieved by using the concept of arranged autoregression which is an autoregression with cases rearranged based on the values of a particular regressor. For the GTAR model in (1), arranged autoregression becomes useful if it is arranged according to the threshold variable XI. Let p = max{pl, P2}. We assume that the first p observations, y l , y z . . . . . y p , are fixed. Let rri be the time index of the ith smallest observation of {Xl,p+l-d . . . . . Xl,,,-d}. By conditioning on the first p observations, one can write the likelihood function as L(Ol,O2,a~,cxZ~,r,d[~)

(x ~71scr2(n-p-s) 2 xexp

-~al2 i=1

Yni+d-

Z.(1) ~ ,,(1) k P) xj:,+d I - pl qg) Y~,+a-j -- 2..,, j=l j=l

C. W.S. Chen I Statistics

& Probability

Letters 40 (1998)

15-22

k

P2

Yni+d - $r’ - C

4y)Yn,+&j

- C

.j=l

wheres

(1)

satisfies qn,Gr

(1)

17

,

fiy’~j,~,+d

(2)

j=l

(1)

(2)

(2)

(2)

(2)

< qx,+,, @I = Cd+,,d, ,...,d$l’,B, ,...,!$‘))*,02 = (q& ,9, ,...,$J~~,B, , = {yt, xi,+ i = l,..., k;t = 1, . . . , n}. The parameters of GTAR model to be esti-

. . . . /?f’)*, and 9 mated are @1,@2,+;, r, and d. Let Y; = (yn,+d, .%2+dr...,_h,+d)* and Y; = (y,~+,+d,...,yn._,+d)* be the respective observations generated by regime I and regime II in the order of occurrence; XT =

(x~tj_&x~~ldY...y~~$)Ty and xl

$2) =

=

(l,Yn,+d-I,

XT = (x~,:‘,+d9.. .y~~f!,+d)Tp where 4”’

Yn,+d-1,.

. .,Yn,+d-ply..

.,Xk,n+d)T

This is an arranged autoregression with the first s cases of in the first regime, and the remaining n - p - s cases in the second

, . . . , xk,n+d)*.

. . ..Yn.+d--pz

*

(~l,n,+d,~l,n~+d,..,Xl,n._p+d)

= (1,

regime. To implement the proposed analysis, it remains for us to derive the conditional posterior distribution of an unknown parameter given all the others. To this end, we choose the priors as follows. We take 81, 632 to be independent N(&, Vi’), and o:, U; independent IG(ri/2,rilZi/2), i = 1,2, where IG denotes the inverse gamma distribution and the hyper-parameters are assumed to be known. Moreover, we assume that r follows a uniform distribution on U(a, b) and d follows a discrete uniform distribution on the integers 1,2,. . . , do. Our interest lies in the marginal posterior distributions of 81, 02, o:, OS, r, and d. Using standard Bayesian techniques, we obtain the following conditional posterior distributions: 1. The conditional posterior distribution of @i is independent of @j for i # j.

where

with ki = (Xi*TXi*)-‘Xi*TYi*. 2. The conditional posterior distribution of c’ is independent of c$, for

Uiizi +

i.e.

niSf

0’

N x;,+~,,

i # j.

i = L2,

where n--P

n’ =

cqx,.., Qr),

t=1

with Pi =X,*@i.

n--P

n2 =

x r{x,.,>+ t=1

sf= n[‘(y,’

- fi)T(q* - pi)

C W.S. Chen I Statistics & Probability Letters 40 (1998) 15-22

18

3. The conditional posterior probability function of r is

p(rlOl, 02, 0-~,a 2, d, 9 ) c< exp

- Z i=1

(3)

2-~..2(Yi* - X i• O i ) T(Yi * - X i * O i )

.2 /0-1., 0"2"

t

Notice that nl and n2 are functions of r. 4. The conditional posterior probability function of d is a multinomial distribution with probability

p(dlOl,02,0"~,0"Z,r,~) = L(Ol,02,0"~,0"Z,r,d[N)

L(Ol,O2,0"~,0"Z,r,d[~),

(4)

--d=l

where d = 1,2 ..... do and

L(Ot,O2,0"~,0"2,r,d[~) = exp

~_2(y~ - - X i, Oi) T(Y/ • - X i• Oi)

--

n, 0"~2. /0"1

i=1

All conditional densities are identified with the exception of r. With regard to r, we will employ the Metropolis algorithm (Metropolis et al., 1953). Detailed discussions regarding Gibbs sampling using Metropolis algorithm can be found in Chen and Lee (1995). Note that the joint posterior of r and d, as given in Geweke and Terui (1993) for TAR model, need not be employed in our procedure. We utilized the Gibbs sampler which is especially useful in extracting marginal distributions from fully conditional distributions when the joint distribution is not easily obtained. Our approach avoids sophisticated analytical and numerical multiple integration.

3. Illustrative example To illustrate the usefulness of GTAR models, we consider two series of quarterly nominal gross national product (GNP) and money supply (M1) of Taiwan from the first quarter of 1965 to the first quarter of 1996, a total of 125 observations. Economists are not in agreement as to what specific items constitute the economy's money supply. Narrowly defined (and designated as M1), the money supply is composed of the following items of all commercial and cooperative banks: (1) coins, (2) paper money, and (3) demand deposits. Fig. 1(a) and (b) give time plots of GNP and Ml. Note that GNP is the current market value in dollars of all final goods and services produced in the economy in a given period. The data are obtained from the AREMOS data base. It is well known that money stock and current dollar measures of economic activity are positively correlated, and movements in monetary base (M1) precede nominal GNP. In our analysis, we focus on the growth rate, namely

Yt = Iog(GNPt ) - log(GNPt- l ) and

xt = log(Mr, t ) - log(Ml.t_ 1)Consequently, there are 124 data points each in the two series. The growth series of GNP and M1 are given in Fig. 2(a) and (b). We consider GTAR model given as follows:

Y'

=

{ (a]l)yt_l+fll')xt+al 1)

if

xt-a
¢~2)yt_,+fl~')xt+al 2)

if

Xt-d > r.

(5)

C.W.S. Chen I Statistics & Probability Letters 40 (1998) 15-22

19

xlO0 19200

'

'

'

.

.

.

.

.

.

.

14400

9600

4800

0 i

i

i

I

i

12

24

36

48

60

(a)

7~2

i

i

i

i

84

96

108

120

4I

I 96

I 108

t 120

time

xlO0 16800

12600

8400

4200

0

..........

I 12

" 24

3i6

I 48

~b)

I 60

| 72

8

time

Fig. 1. Part (a): Time plot of the GNP. Part (b): Time plot of the M1. Table 1 The parameter estimates of Example 1 Regime 1 2

fll (std)

t~l (std)

0.2 (std)

r (std)

d

0.5298 (0.114) 0.3942 (0.169)

-0.2424 (0.123) 0.3888 (0.193)

0.0018 (0.0004) 0.0020 (0.0006)

0.05490 (0.0199)

2

The hyper-parameters used are O0i = 0, V i = 0.05, /)i = 3, and 2i = #2/3 for i = 1 and 2, where #2 is the residual mean squared error of fitting an AR(1) regression time series model to the data. Moreover, we choose do = 4, a = gals, and b = ga85, where gak stands for the kth percentile of the data. The Gibbs sampler is run for 45 000 iterations. We record every 10th value in the sequence of the last 3000 in order to have more nearly independent contributions. The posterior mode of d and the posterior means together with standard deviations (std) of ill, q~l, a, and r are given in Table 1. Fig. 3 displays histogram and estimated marginal posterior density for fll and 4)1 in each regime. The results indicate that the growth rate of GNP is strongly related to growth rate of Ml and lagged growth rate of GNP. There may be question as to whether the sampler can get bogged down with moving between r

20

C W,S. Chen / Statistics & Probability Letters 40 (1998) 15-22 1.62

)

I

I

I

I

I

I

.94

.26

-0.42

-i.i0 I

I

I

I

I

I

I

I

I

I

12

24

36

48

60

72

84

96

108

120

I 72

I 84

I 96

108

I 120

(a)

time I

i. 62

.94

.26

-0.42

-I.i0

(b)

I 12

I 24

i 36

I 48

I 60 time

Fig. 2. Part (a): Time plot o f the growth series o f GNP. Part (b): Time plot o f the growth series o f MI.

and d. A 3D histogram of r and d is provided in Fig. 4. It is clear to see that all realizations of d have mass of 2, due to the fact that the conditional posterior probability function of d in (4) is a multinomial distribution and that the numerator in the computed posterior probability for a given value of d induces huge difference for each probability. Therefore, the probability of d=2 in (4) is approximately 1.0 and the probability of all other values of d is approximately 0.0. To check the adequacy of GTAR model, we consider the standardized series ali)/ai. The autocorrelation and partial autocorrelation functions of these standardized residuals are all within their (asymptotic) two standard error limits, suggesting that there is no significant residual serial correlation. The GTAR model in (5) fits the data very well. The most striking feature of the GTAR model in (5) is that it reveals through the sign of ~bl in Table 1 and Fig. 3 a negative growth in lagged growth rate of GNP when the growth rate of M1,t-2 is below 5.49% and a positive growth when the growth rate of Ml,t-2 is above 5.49%. The results also indicate that it is appropriate to make use of GTAR model rather than the usual linear time series model. The increase in the money supply tends to lower the interest rate and stimulate investment spending. Subsequently, the higher investment spending contributes to the total GNP. Moreover, the expansionary monetary policy may have significant effect on the output, but the effect is augmented when the growth rate of Ml,t-2 is above 5.49%. The quantity of money exerts a powerful influence on the price level when the growth rate of Ml,t-2

21

C.W.S. Chen / Statistics & Probability Letters 40 (1998) 15-22 . . . .

!

.

.

.

.

.

.

.

.

~

. . . .

I

. . . .

60

58

48

•

3e

3e

...................

20

48

2e

.............................................

10

10

: . . . . . . . . . . . . . :.

I

-1

I

i

-e.5

i

i

e

.~

~Bi___ . . . . . . . . . . . . . . . . . . . . . .I

,

i

e.5

.

.

.

.

I

1

.

.

.

.

I

~.5

4e

4e

38

30

211

211

18

18

-1

-8.5

8

e.6

1

1.5

-1

-0.5

e

8.5

1

1.5

-1

-e.S

e

e.5

1

1.5

Fig. 3. Histogram and estimated marginal posterior distribution for El and q~l in each regime by using the simulated realizations.

is above 5.49%. In the long run, money is the cause of inflation and inflation is a factor of economic growth. Hence, monetary policy is important for economic growth as well as economic fluctuation.

Acknowledgements The author thank the anonymous referee and my colleague, Professor Hsianghoo Steve Ching of Department of Public Finance, for constructive comments and valuable discussion. The author would like to thank research support from National Science Council of Taiwan (grant NSC85-2121-M-035-007).

22

C W.S. Chen / Statistics & Probability Letters 40 (1998) 15-22

68 58

28 4 le 9

2

3

0°1 0°2 r

Fig. 4. 3D histogram of r and d.

References Broemeling, L.D., Cook, P., 1992. Bayesian analysis of threshold autoregressions. Commun. Statist. Theory Methods 21, 2459-2482. Chen, C.W.S., Lee, J.C., 1995. Bayesian inference of threshold autoregressive models. J. Time Ser. Anal. 16, 483-492. Deutsch, M., Granger, C.W.J., Tefiisvirta, 1994. The combination of forecasts using changing weights. Intemat. J. Forecasting 10, 47-57. Gelfand, A.E., Hills, S.E., Racine-Poon, A., Smith, A.F.M., 1990. Illustration of Bayesian inference in normal data models using Gibbs sampling. J. Amer. Statist. Assoc. 85, 972-985. Gelfand, A.E., Smith, A.F.M., 1990. Sampling-hased approach to calculating marginal densities. J. Amer. Statist. Assoc. 85, 398--409. Geman, S., Geman, D., 1984. Stochastic relaxation, Gibbs distributions and the Bayesian restoration of images. IEEE Trans. Pattern Anal. Machine Intell. 6, 721-741. Geweke, J., Terui, N., 1993. Bayesian threshold autoregressive models for nonlinear time series. J. Time Ser. Anal. 14, 441-454. Goldfeld, S.M., Quandt, R.E., 1972. Nonlinear Methods is Econometrics. North-Holland, Amsterdam. Granger, C.W.J., Tef,isvirta, 1993. Modelling Dynamic Nonlinear Relationships. Oxford University Press, Oxford. McCulloch, R.E., Tsay, R.S., 1993. Bayesian analysis of threshold autoregressive processes with a random number regimes. Proc. 25th Symp. on the Interface, pp. 253-262. Metropolis, N., Rosenbluth, A.W., Rosenbluth, M.N., Teller, A.H., 1953. Equations of state calculations by fast computing machines. J. Chem. Phys. 21, 1087-1091. Quandt, R.E., 1958. The estimation of the parameters of a linear regression system obeying two separate regimes. J. Amer. Statist. Assoc. 53, 873-880. Tanner, M.A., Wong, W.H., 1987. The calculation of posterior distributions by data augmentation (with discussion). J. Amer. Statist. Assoc. 82, 528-550. Tong, H., 1978. In: Chen, C.H. (Ed.), On a Threshold Model in Pattern Recognition and Signal Processing, Sijhoff & Noordhoff, Amsterdam. Tong, H., 1983. Threshold Models in Non-linear Time Series Analysis, vol. 21. In: Krickegerg, K. (ed.), Lecture Notes in Statistics. Springer, New York. Tong, H., 1990. Nonlinear Time Series: A Dynamical System Approach. Oxford University Press, Oxford. Tong, H., Lim, K.S., 1980. Threshold autoregression, limit cycles and cyclical data (with discussion). J. Roy. Statist. Soc. Ser. B 42, 245-292. Tsay, R.S., 1989. Testing and modeling threshold autoregressive process. J. Amer. Statist. Assoc. 84, 231-240.