Regime Switching Regression Models

Readings and References:

• S. M. Goldfeld and R. E. Quandt, "Techniques for Estimating Switching Regressions," in Nonlinear Methods in Econometrics, North-Holland, 1972, Chapter 1, 3-35.
• J. D. Hamilton, "Estimation, Inference and Forecasting of Time Series Subject to Changes in Regime," in Handbook of Statistics, Vol. 11, ed. by G. S. Maddala, C. R. Rao, and H. D. Vinod, North-Holland, 1993, 231-260.
• R. Garcia and P. Perron, "An Analysis of Real Interest Rate under Regime Shift," Review of Economics and Statistics, 1996, 78, 115-125. (Paper)
• C.-J. Kim and C. R. Nelson, State-Space Models with Regime Switching, Classical and Gibbs-Sampling Approaches with Applications, The MIT Press, 1999, Chapter 4: Markov-Switching Models, 59-95.

Consider the M-state (regime or structure) linear regression equation:

Y_t = X_tb⁽ⁱ⁾ + e_t⁽ⁱ⁾
e_t⁽ⁱ⁾ ~ normal(0,s²⁽ⁱ⁾)

where [Y_t,X_t] is the observed sample data, t=1,2,...,N, and (i) is the state (regime) indicator, i=1,2,...,M. In particular, assuming normal distribution, the conditional probability density function of Y_t occured in the regime i is written as:

f(Y_t;b⁽ⁱ⁾,s²⁽ⁱ⁾) = 1/Ö(2ps²⁽ⁱ⁾) exp(-½(Y_t-X_tb⁽ⁱ⁾)²/s²⁽ⁱ⁾)

In general, the regression parameters could change several times (that is, switching backward and forward) within M possible regimes. Typically, M=2.

Deterministic Switching

Dummy Variable Approach

Dit =	1, if t in N(i)
	0 otherwise

Since D_it is observed, given the conditional probability density function f(Y_t;b⁽ⁱ⁾,s²⁽ⁱ⁾) of Y_t, the log-likelihood function for model estimation is:

ll(q) = å_t=1,2,...,Nln å_i=1,2,...,M D_itf(Y_t;b⁽ⁱ⁾,s²⁽ⁱ⁾)

where q = (b⁽ⁱ⁾,s²⁽ⁱ⁾, i=1,2,...,M)' is the parameter vector of the model.

Stochastic Switching

Dit =	1, if St = i (the regime is i at time t)
	0 otherwise

D_it is a random discrete variable taking value either 1 or 0, which is not directly observed. Let D_it^* be the latent variable for D_it such that D_it^* > 0 if D_it = 1, and D_it^* £ 0 otherwise (D_it = 0). We write:

D_it^* = a⁽ⁱ⁾ + u_t

Or, in a more general framework, D_it^* = Z_ta⁽ⁱ⁾ + u_t where Z_t is a set of exogenous or predetermined variables determining the regime structure.

When D_it^* > 0, u_t > - a⁽ⁱ⁾. Assuming u_t ~ normal(0,1), the probability of D_it = 1 is represented by the cumulative normal density evaluated from -¥ to a⁽ⁱ⁾. Then the probability that S_t=i is expressed as:

Pr(St=i;a(i))	= Pr(Dit=1) = Pr(Dit*>0) = Pr(ut>-a(i))
	= ò-¥ a(i) 1/Ö(2p) exp(-½ z2)dz

Of course, å_i=1,2,...,M Pr(S_t=i;a⁽ⁱ⁾) = 1. Then the conditional joint probability density function of Y_t and S_t=i is

f(Y_t,S_t=i;b⁽ⁱ⁾,s²⁽ⁱ⁾,a⁽ⁱ⁾)
= Pr(S_t=i;a⁽ⁱ⁾) f(Y_t|S_t=i;b⁽ⁱ⁾,s²⁽ⁱ⁾)
= ò_-¥ ^a(i)1/Ö(2p) exp(-½ z²)dz [1/Ö(2ps²⁽ⁱ⁾) exp(-½(Y_t-X_tb⁽ⁱ⁾)²/s²⁽ⁱ⁾)]

Let q⁽ⁱ⁾ = (b⁽ⁱ⁾,s²⁽ⁱ⁾,a⁽ⁱ⁾)' and q =(q⁽ⁱ⁾, i=1,2,...,M)'. Then, the unconditional probability density function of Y_t is

f(Y_t;q) = å_i=1,2,...,M f(Y_t,S_t=i;q⁽ⁱ⁾)

The resulting log-likelihood function is:

ll(q)	= åt=1,2,...,N ln f(Yt;q)
	= åt=1,2,...,N ln åi=1,2,...,M f(Yt,St=i;q(i))

As a by-product of the maximum likelihood estimation of the parameter vector q, the estimated conditional probability for the observation t occurred in regime i is computed as follows:

Pr(St=i|Yt;q)	= f(Yt,St=i;q(i)) / f(Yt;q)
	= f(Yt,St=i;q(i)) / åi=1,2,...,M f(Yt,St=i;q(i))

Mixture of Probability Distributions Revisited

f(Y_t|S_t=i;m⁽ⁱ⁾,s²⁽ⁱ⁾) = 1/Ö(2ps²⁽ⁱ⁾) exp(-½(Y_t-m⁽ⁱ⁾)²/s²⁽ⁱ⁾) with probability Pr(S_t=i), i=1,2.

Let Pr(S_t=1)=p, then Pr(S_t=2) = 1-p. Here p>0 must be estimated together with the parameters m⁽ⁱ⁾ and s²⁽ⁱ⁾ for i=1,2. The log-likelihood function is:

ll(p,m⁽¹⁾,m⁽²⁾,s²⁽¹⁾,s²⁽²⁾) = å_t=1,2,...,N ln [p f(Y_t|S_t=1;m⁽¹⁾,s²⁽¹⁾) + (1-p) f(Y_t|S_t=2;m⁽²⁾,s²⁽²⁾)]

Finally, the conditional probability of INCOME variable Y_t drawn from regime 1 is

p f(Y_t|S_t=1;m⁽¹⁾,s²⁽¹⁾) / [p f(Y_t|S_t=1;m⁽¹⁾,s²⁽¹⁾) + (1-p) f(Y_t|S_t=2;m⁽²⁾,s²⁽²⁾)]

Markov Switching

Pr(S_t=i|S_t-1=j) = p_ij > 0
å_i=1,2,...,M p_ij = 1

Then the joint probability of S_t=i and S_t-1=j is defined by

Pr(St=i,St-1=j)	= Pr(St=i|St-1=j)Pr(St-1=j) = pijPr(St-1=j)
and Pr(St=i)	= åj=1,2,...,M Pr(St=i,St-1=j)
	= åj=1,2,...,M pijPr(St-1=j)

For a two-state Markov switching process (M=2), we write:

p = p11 = Pr(St=1|St-1=1)	1-p = p21 = Pr(St=2|St-1=1)
q = p22 = Pr(St=2|St-1=2)	1-q = p12 = Pr(St=1|St-1=2)

Then,

Pr(S_t=1) = p Pr(S_t-1=1) + (1-q) Pr(S_t-1=2)
Pr(S_t=2) = (1-p) Pr(S_t-1=1) + q Pr(S_t-1=2)

The probabilities of a steady state (that is, S = S_t = S_t-1) are:

Pr(S=1) = (1-p)/(2-p-q)
Pr(S=2) = (1-q)/(2-p-q)

Let Pr(S_t=1) = Pr(S_t^*£0) and Pr(S_t=2) = Pr(S_t^*>0), where S_t^* is a latent variable for the discrete state indicator S_t defined by:

S_t^* = a₀ + a₁ (S_t-1-1) + u_t, u_t ~ normal(0,1)

Therefore,

p = Pr(St=1|St-1=1)	= Pr(St*£0|St-1=1) = Pr(ut£-a0)
	= ò-¥ -a01/Ö(2p) exp(-½ z2)dz
q = Pr(St=2|St-1=2)	= Pr(St*>0) = Pr(ut>-a0-a1)
	= ò-a0-a1 ¥1/Ö(2p) exp(-½ z2)dz
	= 1 - ò-¥ -a0-a11/Ö(2p) exp(-½ z2)dz

It clear that through probability transformation, 0 < p < 1 and 0 < q < 1. An alternative approach to ensure that 0 < p < 1 and 0 < q < 1 is to use logistic transformation as follows:

p = exp(p₀)/(1+exp(p₀))
q = exp(q₀)/(1+exp(q₀))

where p₀ and q₀ are unrestricted real numbers.

Dynamic Model

Let H_t-1 = {Y_t-1,Y_t-2,...} be the historical information up to time t-1. Suppose the conditional probability density function of Y_t (conditional to S_t=i, S_t-1=j and H_t-1, ignoring the parameter vector q for the moment) is f(Y_t|S_t=i,S_t-1=j,H_t-1). Then the conditional joint probability density function of Y_t, S_t, and S_t-1 (conditional to H_t-1) is:

f(Y_t,S_t=i,S_t-1=j|H_t-1) = Pr(S_t=i,S_t-1=j|H_t-1) f(Y_t|S_t=i,S_t-1=j,H_t-1)

and

f(Y_t|H_t-1;q) = å_i=1,2å_j=1,2 f(Y_t,S_t=i,S_t-1=j|H_t-1;q)

where q is the parameter vector. Finally, the conditional log-likelihood function is:

ll(q) = å_t=1,2,...,N ln f(Y_t|H_t-1;q)

Hamilton's Algorithm

(1)	Pr(St=i,St-1=j|Ht-1) = Pr(St=i|St-1=j) Pr(St-1=j|Ht-1)

Starting t=1, let Pr(S₀=1|H₀) = (1-p)/(2-p-q) and Pr(S₀=2|H₀) = (1-q)/(2-p-q). These are the steady-state or unconditional probability of S_t=1 and S_t=2, respectively.

From the computed Pr(S_t=i,S_t-1=j|H_t-1) and given f(Y_t|S_t=i,S_t-1=j,H_t-1), compute the joint pdfs f(Y_t,S_t=i,S_t-1=j|H_t-1) for all combination of i,j=1,2, and therefore the conditional likelihood at t is:

(2)	f(Yt|Ht-1)	= åi=1,2åj=1,2 f(Yt,St=i,St-1=j|Ht-1)
		= åi=1,2åj=1,2 Pr(St=i,St-1=j|Ht-1) f(Yt|St=i,St-1=j,Ht-1)

By definition, Pr(St=i,St-1=j|Ht)	= Pr(St=i,St-1=j|Yt,Ht-1)
	= f(Yt,St=i,St-1=j|Ht-1) / f(Yt|Ht-1)
	= f(Yt|St=i,St-1=j,Ht-1)Pr(St=i,St-1=j|Ht-1) / f(Yt|Ht-1)
and Pr(St=i|Ht)	= åj=1,2 Pr(St=i,St-1=j|Ht)

To iterate from t=1 to N, use Pr(S_t=i|H_t) to compute Pr(S_t+1=k,S_t=i|H_t) as (1) and the conditional likelihood f(Y_t+1|H_t) as (2). Finally, we sum all the log of component log-likelihood for maximization with respect to the parameter vector q.

The above algorithm may be generalized straightforwardly to consider M cases (M>2) and L lags (L>1), in which the probabilities and pdfs are evaluated over a large number of M^L+1 cases.

For parameter estimation, we apply standard numerical maximization of the log-likelihood function for the entire sample. The alternative method is the EM (Expectation-Maximization) Algorithm.

To make inference about the state variable S_t, we can use filtered probability: Pr(S_t=i|H_t), or smoothed probability: Pr(S_t=i|H_N). The later refers to the probability of S_t=i conditional on all the information in the sample. The smoothed probability is obtained by updating the filtered probability using all information.

Example: Hamilton's Model of Business Fluctuations

Y_t = 100*(ln(GDP_t)-ln(GDP_t-1))

Then 2-state ("boom" and "burst", or expansion and recession) growth in real GDP is modeled as an AR(4) process:

Y_t = m⁽ⁱ⁾ + e_t, i=1,2
e_t = r₁e_t-1 + r₂e_t-2 + r₃e_t-3 + r₄e_t-4 + u_t

where u_t ~ normal(0,s²). By assuming the first-order Markov switching process, the transition probabilities between two consecutive states S_t and S_t-1 are Pr(S_t=1|S_t-1=1)=p and Pr(S_t=2|S_t-1=2)=q. By definition, Pr(S_t=1|S_t-1=2)=1-q and Pr(S_t=2|S_t-1=1)=1-p.

The model can be represented equivalently as:

(Y_t-m⁽ⁱ⁾) = r₁(Y_t-1-m^(j)) + r₂(Y_t-2-m^(k)) + r₃(Y_t-3-m^(l)) + r₄(Y_t-4-m^(m)) + u_t

with i,j,k,l,m = 1,2. In total, there are 32 (=2⁴⁺¹) cases of probabilities and pdfs forming the likelihood function for maximization with respect to the parameter vector. In addition to the model estimation, both filtered and smoothed probabilities are computed for inference about the state variable S_t: Pr(S_t|H_t) and Pr(S_t|H_N) for each observation t. (See Program, require GPE2 application module MARKOV1.GPE [draft version]).

Example: Three-State Markov Switching Model of the Real Interest Rate

(Y_t-m⁽ⁱ⁾) = r₁(Y_t-1-m^(j)) + r₂(Y_t-2-m^(k)) + u_t
u_t ~ normal(0,s²⁽ⁱ⁾)

with i,j,k = 1,2,3, where Y_t is the ex-post real interest rate calculated by subtracting the CPI-based inflation rate from the nominal interest rate (three-month treasury bill rate). See Data.

Transition probability for the 3-state 1st-order Markov process is defined by p_ij = Pr(S_t=i|S_t-1=j) > 0, where å_i=1,2,3 p_ij = 1 for j=1,2,3.

In total, there are 27(=3²⁺¹) cases of probabilities and density functions forming the likelihood function for maximization with respect to 14 parameters in the model:

m⁽¹⁾, m⁽²⁾, m⁽³⁾
r₁, r₂
s²⁽¹⁾, s2⁽²⁾, s²⁽³⁾
p₁₁, p₁₂, p₁₃, p₂₁, p₂₂, p₂₃

Formulate and estimate the log-likelihood function for Garcia-Perron Model of real interest rate. (Program)