I am sorry that I haven’t posted any new blog entries for the last couple months. I have been busy with my work and my side project which I do not want to discuss openly. Anyway, I would like to talk about Markov-switching regression today.

When we run time series regressions, it is usually the case where our estimated coefficients may not be stable and be subject to regime shifts. (This can be checked by CUSUM squared tests and other break point tests). To model this kind of regime shift , several types of econometric models were developed by Econometricians. One of them is the Markov-switching type regressions, which has gained recognition ever since the influential work by Hamilton (1989) . Nowadays, we have many statistical packages to run markov switching univariate or multivariate regressions. However, I never found anyone trying to post the manual code for the markov-switching regressions. So, let us write the code together! 🙂

The code I am going to present is the Quasi-likelihood Estimation method for a dynamic 2-state markov-switchhing AR(1) model:

Let be the observed variables. Assume our model has the following form:

where signals 2 different states. and are subject to regime shift depending on the state changes

First, Assuming asymptotic normality, and let , we can write the density of conditional on and as:

Second, we have:

Third, we have:

The fourth equation we have is:

where and are the transition probabilities from previous state to the next state. For further details of the Markov-switching model, please refer to most papers by Hamilton (1989,1994, 2005) etc.. Here, I only briefly outlined the key equations to set up the likelihood equation

Now if we are given initial values of , we iterate the 4 equations to obtain the sum of the likelihood. Then the quasi-likelihood equation is:

To maximize the equation above, we can use the optimx package in R or the nlm function in R

Here I post my R code for setting up the quasi-likelihood function for the example above:

likeli_fun <-function(param){ y <- data[,1] #Here, use your own dataset, dependent variable, put it in column 1 alpha<-param[1] #Intercept in regime 0 alpha1 <- param[2] #intercept in regime 1 var <-param[3]^2 # variance of the error p00 <- param[4] #The transition probability p11 <-param[5] beta <-param[6] #AR(1) slope in regime 0 beta1 <-param[7] #AR(1) slope in regime 1 #Set initial value for prob1 and prob0: A <- matrix(c(1-p00,-1+p11,1,-1+p00,1-p11,1),nrow=3,ncol=2) #This is the A matrix mentioned in Hamilton 1994 inv_A <-solve(t(A)%*%A) sol_A <- inv_A%*%t(A) #Why we do this? Please check out Hamilton (1994) prob0<-sol_A[1,3] #Recommended initial values by Hamilton (1994) prob1<-sol_A[2,3] #Start iterating to get the likelihood T <- length(y) likv <-0 j_iter <-2 while (j_iter <= N) { eps0 <- y[j_iter]-alpha-beta*y[j_iter-1] eps1 <- y[j_iter]-alpha1-beta1*y[j_iter-1] f0 <- (1/sqrt(2*pi*var))*exp(-0.5*eps0*eps0/var) #equation 1 for regime0 f1 <- (1/sqrt(2*pi*var))*exp(-0.5*eps1*eps1/var) #equation 1 for regime1 con_fz <- prob0*f0+prob1*f1 #Equation 2 fil_prob0 <- (prob0*f0)/con_fz #Equation 3 fil_prob1 <- (prob1*f1)/con_fz prob0 <- p00*fil_prob0+(1-p00)*fil_prob1 #Equation 4 prob1 <- (1-p11)*fil_prob0+p11*fil_prob1 lik <- log(con_fz) # Take the log likv <- likv-lik # add to the last likelihood j_iter <- j_iter+1 #Start another round of iteration } return(likv/T) }

Next step, all you need to do in R is to fire up nlm or optimx in R to maximize the likelihood to compute the estimates of the parameters.

In my next post, let us start using the packages in R and command in Eviews to do some interesting analysis using Markov-switching regressions.