Asset Pricing [9b : Regime Switching]

Continuing from the previous post where we downloaded,transformed (where necessary) and stored the data series for each of the two categories of variables (economic regime variables & risk premia),the remainder of this series of posts shall be concerned with the following tasks :

  1. Fit a two state hidden markov model to economic regime variables (Equity turbulence,Currency Turbulence,Inflation,Economic Growth).
  2. For each regime variable,overlay the resulting probability state map on a time series plot of the regime variable and save final graph as an element of a list.
  3. For each regime variable,extract from the fitted object a table of estimation results and store the table plot as an element of a list.
  4. For each hedge fund strategy,ensure time series are aligned across economic regime variables.
  5. Explore In-Sample Performance across funds :
    • Event-state Mean vs Normal-state Mean.
    • Scaled difference in state means.
    • Cumulative returns & Drawdowns when states are known vs unknown
  6. Create Dashboard of results for chosen Fund/Regime combination.

Although the code was not initially created with a dashboard in mind,it is probably for the best to just show the dashboard of plots and tables at the end as opposed to displaying the outputs of intermediate steps. 


 

[1] Fitting the hidden markov model

The Turbulence.list object holds our regime variables as elements of a list. To fit a two-state hidden markov model to each of the variables,I used the lapply function to pass each element of that list to a custom HMfit() function.
g

#[3] Fit the turbulence indices where available
	Turbulence.fit <- lapply(Turbulence.list,function(x) HMfit(x))

k
The HMfit() function calls the fitting functions in the depmixS4 package and returns a list object as given below :
j

#########################################################################################
# Hidden Markov Fit
#########################################################################################
HMfit <- function(inp.ret){

	ret.xts <- inp.ret
	states <- 2

	model.spec <- depmix(eval(parse(text=names(ret.xts)))~1,nstates=states,data=ret.xts)
	model.fit <- fit(model.spec)
	model.posterior <- posterior(model.fit)
        sum.frame <- summary(model.fit)
	normal.idx <- which(sum.frame[,1]>0 & sum.frame[,1]==min(sum.frame[,1]))
	event.idx <- index(sum.frame)[-normal.idx]

	HM.fit <- list(
		normal.state = normal.idx,
		event.state = event.idx,
		turb =ret.xts,
		spec =model.spec,
                fit = model.fit,
                summary = summary(model.fit),
		normal.prob = model.posterior[,normal.idx+1],
                normal.class = model.posterior[,'state']==normal.idx,
		event.prob = model.posterior[,event.idx+1],
		event.class = model.posterior[,'state']==event.idx
	)

	HM.fit$df = data.frame(Time.Idx=index(HM.fit$turb),
			       Turbulence=HM.fit$turb,
                               Normal.Prob=HM.fit$normal.prob,
			       Event.Prob=HM.fit$event.prob,
			       Normal.Class=ifelse(HM.fit$normal.class==T,1,0),
			       Event.Class=ifelse(HM.fit$event.class==T,1,0),
			       row.names=NULL)
  return(HM.fit)
}

l
The most useful structure this function returns is a data frame containing :a time stamp,the turbulence index,state probabilities and state maps.

 

[2] Plotting economic regime variables and overlay with state maps

For a reason I cannot understand, the state map which shows when a particular economic regime variable is in the event (or normal) state and which therefore should be a line plot of either 0 or 1, becomes jagged when it should be straight. To overlay the state map on top of a timeseries plot of an economic regime variable,I simply multiplied the event state map ( a time series of 1’s and 0s corresponding to when the regime variable is in an event or non-event state ) by the highest value of the regime variable and used the resulting index as an input to the geom_ribbon layer. As before I use the lapply function to pass the previously fitted objects to a custom function TurbPlot().

l

#[4] Plot the Turbulence indices and overlay state maps
	Turbulence.plot <- lapply(Turbulence.fit,function(x) TurbPlot(x))

j
The TurbPlot() function creates one ‘overlay’ data frame for each state (event,normal) from the fitted object and then passes these to the ggplot function ( and most importantly the geom_ribbon layer ):
l

#########################################################################################
# Turbulence and regime Plotting
#########################################################################################
TurbPlot <- function(fitted.obj){

  plot.df <- fitted.obj$df
  ymax <- max(plot.df[,names(fitted.obj$turb)])+1

  overlay.df <- mutate(select(plot.df,Time.Idx,Event.Class),ID=Event.Class*ymax)
    overlay.df$yvals <- fitted.obj$turb
    overlay2.df <- mutate(select(plot.df,Time.Idx,Normal.Class),ID=Normal.Class*ymax)
    overlay2.df$yvals <- fitted.obj$turb

  gg.s1 <- ggplot()+
            theme_tufte(base_size=11)+
            geom_ribbon(data=overlay.df,aes(x=Time.Idx,ymin=0,ymax=ID),alpha=0.6,fill='red')+
            geom_ribbon(data=overlay2.df,aes(x=Time.Idx,ymin=0,ymax=ID),alpha=0.6,fill='lightblue')+
            geom_area(data=plot.df,aes_string(x='Time.Idx',y=names(fitted.obj$turb)),fill='skyblue4',colour='darkblue',size=0.5)+
            labs(x='Time',y='Values',title=paste('Turbulence Index For',names(fitted.obj$turb),'Regime',sep=' '))

return(gg.s1)
}
#########################################################################################

l
The function returns a ggplot to be used in the final dashboard.

 

[3] Construct table of estimation results

Continuing with our list of fitted objects,let’s extract pertinent information from it across the event and normal states and make a simple text plot of mean and sigma estimates across all four regime variables.
;l

#[6] Table of Estimation results
  Turbulence1.tbl <- Extract.est(Turbulence.fit,state='Normal')
  Turbulence2.tbl <- Extract.est(Turbulence.fit,state='Event')
  MarkovEstPlots <- MarkovEstPlot(Turbulence1.tbl,Turbulence2.tbl)

ll
Since we want a table that collects information across regimes,we do not need lapply to act on every element of the list and instead pass the entire list to a custom Extract.est() function.

j

#########################################################################################
#Extracting Estimates from list object
#########################################################################################
Extract.est <- function(fitted.list,state){
 list.store <- list()
 idx.str<-paste('ifelse(state==','"Normal"',',','x$normal.state',',','x$event.state)',sep='')
 list.store <- lapply(fitted.list,function(x) data.frame(
                      Initial_State = paste('State',which(x$fit@init==1)),
  		      Persistence   = paste(round(x$fit@transition[[eval(parse(text=idx.str))]]@parameters$coefficients[eval(parse(text=idx.str))],3)*100,'%',sep=''),
  		      Transition    = paste(round(x$fit@transition[[eval(parse(text=idx.str))]]@parameters$coefficients[-eval(parse(text=idx.str))],3)*100,'%',sep=''),
  		      Mean          = x$summary[eval(parse(text=idx.str)),1],
                      Std_Dev	    = x$summary[eval(parse(text=idx.str)),2],
                      stringsAsFactors=F
                    ))
  est.tbl <- do.call(rbind,list.store)
  return(est.tbl)
}
#########################################################################################

f
This function returns a data frame that contains information on (1)The initial state of the regime variable,(2)The probability of staying in one state,(3) The probability of transitioning to another state,(3) The mean estimate and (4) The sigma estimate across all regimes.

Once we have extracted the desired data from our fitted objects,I use the MarkovEstPlot() function to make relevant stacked bar charts and associated value tables to be used in the final dashboard.

l

#########################################################################################
# Markov Estimates Plotter
#########################################################################################

MarkovEstPlot <- function (tbl.norm,tbl.event){

	tbl.norm <- coredata(tbl.norm)
	tbl.event <- coredata(tbl.event)
	
	mean.df <- data.frame(c(rownames(tbl.norm)),select(tbl.event,Mean),select(tbl.norm,Mean),stringsAsFactors=F)
	colnames(mean.df) <- c('Regime','Event','Normal')
	mean.tbl <- mean.df
	
	sd.df <-data.frame(c(rownames(tbl.event)),select(tbl.event,Std_Dev),select(tbl.norm,Std_Dev),stringsAsFactors=F)
	colnames(sd.df) <- c('Regime','Event','Normal')
	sd.tbl <- sd.df
	
	pers.df <- data.frame(c(rownames(tbl.norm)),select(tbl.event,Persistence),select(tbl.norm,Persistence),stringsAsFactors=F)
  colnames(pers.df) <- colnames(mean.df)
	pers.df <-data.frame(Regime=c(rownames(tbl.norm)),apply(select(pers.df,Event,Normal),2,function(x) round(as.numeric(sub('%','',x)),4)),stringsAsFactors=F)
	pers.tbl <- pers.df
	
	trans.df <- data.frame(c(rownames(tbl.norm)),select(tbl.event,Transition),select(tbl.norm,Transition),stringsAsFactors=F)
  colnames(trans.df) <- colnames(mean.df)
	trans.df <- data.frame(Regime=c(rownames(tbl.norm)),apply(select(trans.df,Event,Normal),2,function(x) round(as.numeric(sub('%','',x)),4)),stringsAsFactors=F)
	trans.tbl <- trans.df
	
	mean.df <- gather(mean.df,variable,value,-Regime)
	sd.df<-gather(sd.df,variable,value,-Regime)
	pers.df<-gather(pers.df,variable,value,-Regime)
	trans.df<-gather(trans.df,variable,value,-Regime)
	colnames(mean.df)[2] <- colnames(sd.df)[2] <-colnames(pers.df)[2] <- colnames(trans.df)[2] <-'State' 
	
	
	gg.m <- ggplot()+
		        theme_economist(base_size=5)+
		  			geom_bar(data = mean.df,aes(x=Regime,y=value,fill=State),stat = "identity",position='dodge')+
		  			scale_fill_brewer(palette='Set1')+
						labs(title='Mean Estimates Across\nRegimes & States')+
		  			coord_flip()+
      			theme(axis.ticks=element_blank(),legend.position='none',axis.title=element_blank(),axis.text.x=element_blank(),panel.grid.major = element_blank(), panel.grid.minor = element_blank())

	
	gg.sd <- ggplot()+
				    theme_economist(base_size=5)+
		  			geom_bar(data = sd.df,aes(x=Regime,y=value,fill=State),stat = "identity",position='dodge')+
		  			scale_fill_brewer(palette='Set1')+
						labs(title='Sigma Estimates Across\nRegimes & States')+
      			coord_flip()+
      			theme(axis.ticks=element_blank(),legend.position='none',axis.title=element_blank(),axis.text.x=element_blank(),panel.grid.major = element_blank(), panel.grid.minor = element_blank())

	  gg.m.tbl <- QuickTbl(mean.tbl,title='Estimates of Mean')
	  gg.sd.tbl <-QuickTbl(sd.tbl,title='Estimates of Sigma')
	
		gg.p <- ggplot()+
				    theme_economist(base_size=5)+
		  			geom_bar(data = pers.df,aes(x=Regime,y=value,fill=State),stat = "identity",position='dodge')+
		  			scale_fill_brewer(palette='Set1')+
					  labs(title='Persistence Estimates Across\nRegimes & States')+
      			coord_flip()+
      		  theme(axis.ticks=element_blank(),legend.position='none',axis.title=element_blank(),axis.text.x=element_blank(),panel.grid.major = element_blank(), panel.grid.minor = element_blank())


		gg.t <- ggplot()+
				    theme_economist(base_size=5)+
		  			geom_bar(data = trans.df,aes(x=Regime,y=value,fill=State),stat = "identity",position='dodge')+
		  			scale_fill_brewer(palette='Set1')+
						labs(title='Transition Estimates Across\nRegimes & States')+
      		  coord_flip()+
      			theme(axis.ticks=element_blank(),legend.position='none',axis.title=element_blank(),axis.text.x=element_blank(),panel.grid.major = element_blank(), panel.grid.minor = element_blank())

		gg.p.tbl <- QuickTbl(pers.tbl,title='Estimates of Persistence (%)')
	  gg.t.tbl <-QuickTbl(trans.tbl,title='Estimates of Transition (%)')
	
						
	return(list(gg.m,gg.m.tbl,gg.sd,gg.sd.tbl,gg.p,gg.t,gg.p.tbl,gg.t.tbl))
}
##########################################################################################

hk
This function returns a list object with 4 stacked bar chart plots and 4 text plots corresponding to mean,sigma,persistence and transition estimates across regimes.

 

[4] Aligning risk premia and regime variable time series

Since returns on various hedge fund strategy returns may have different time indices than those of our regime variables,it is first necessary to ensure that we are looking at data that is correctly aligned in time. This is where mistakes may easily creep up (at least for me) and if there is anything wrong in subsequent outputs,the culprit is likely to lurk in this step. The Idea is to align returns across 13 funds with economic regime variables and store the output in a list object named Aligned.

l

#[7.2] Create list object to hold Turbulence aligned edhec return data
num.turbs <- length(Turbulence.fit)
names.turbs <- names(Turbulence.fit)
wind.beg <- min(index(edhec.xts))
wind.end <- max(index(edhec.xts))

Aligned <- list()
  for(i in 1:num.turbs){
    	Aligned[[i]] <- list()
        	names(Aligned)[i] <- names.turbs[i]
		sub.turb <- filter(Turbulence.fit[[i]]$df,Time.Idx>=wind.beg & Time.Idx<=wind.end)
		min.t <- paste(format(min(sub.turb$Time.Idx),'%Y-%m'),'-01',sep='')
		max.t <- paste(format(max(sub.turb$Time.Idx),'%Y-%m'),'-01',sep='')
		if(names.turbs[i]=='Growth'){
         		edhec.sub <- edhec.xts[sub.turb$Time.Idx,]
		}else{
			edhec.sub <- window(edhec.xts,start=min.t,end=max.t)
		}
	Aligned[[i]] <- cbind(sub.turb[,-1],as.data.frame(edhec.sub))
	}

l

Each element of the Aligned list is named after an economic regime variable and contains  a data frame of aligned fund returns. To get an idea as to the contents of the Aligned list above, let us look at  head(Aligned$Equity)

 sdsd

 

The columns continue for the remaining strategies.Data frames for other aligned time series are accessed by Aligned$Currency,Aligned$Inflation,Aligned$Growth. Now that we have aligned fund returns,regime indices and state probabilities,we are in the position to explore the in-sample performance of each fund.

Advertisements
3 comments
  1. vonjd said:

    This is an excellent piece of work!

    You might not know that but because I was so impressed by Kritzman’s paper I brought it to Quantivity’s attention and he was kind enough to replicate parts of it. Now you follow up on that, which is just great.

    Yes, the RHmm has been deprecated altogether now. I tried to reach the authors but to no avail. I also posted a question here: http://stats.stackexchange.com/questions/116575/why-was-the-popular-rhmm-package-removed-from-cran

    We are also doing research along similar lines, perhaps it would be worthwhile to liaise.

    I am very much looking forward to your next posts! Keep up the good work!

    • Thank you very much for stopping by and also for the kind comment! I also wish that quantivity continues his great blog but for some reason it seems to have been idle for some time. I merely try to replicate things here and there…probably incorrectly/inefficiently which is why my posts are not really meant to be replicated :-)

Leave a Reply

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out / Change )

Twitter picture

You are commenting using your Twitter account. Log Out / Change )

Facebook photo

You are commenting using your Facebook account. Log Out / Change )

Google+ photo

You are commenting using your Google+ account. Log Out / Change )

Connecting to %s

%d bloggers like this: