Archive

Tag Archives: hidden markov model

The final task we have set ourselves, and which also happens to be a welcome initiation to the excellent ggplot2 package for me, is concerned with summarising the previously saved plots and tables into a single (and hopefully useful) dashboard. I have attempted to reduce clutter as much as possible by omitting axis ticks,values and legends in those cases where interpretation is forthcoming. After exhausting the template themes available from the ggthemes package, I have settled on the economist theme for all the plots in the dashboard,finding it to be the most pleasing on the eyes. The function that collects and arranges all the outputs from previous steps is the custom DashboardPlot() function which accepts as arguments the following :[1] List object of dashboard elements,[2]Fund name of choice ,[3] Regime name of choice.
l

#[7.3.3] Dashboard
dashboard.list <- list(FundRet.list=edhec.ret,
	                     Turbulence.plot=Turbulence.plot,
	                     Turbulence1.tbl=Turbulence1.tbl,
	                     Turbulence2.tbl=Turbulence2.tbl,
	                     plot.list=plot.list,
	                     CompareCumu=CompareCumu,
	                     FacetPlot=FacetPlot.list,
	                     FacetObj=Facet.obj,
	                     MarkovPlots=MarkovEstPlots)

#Dashboard for Global Macro fund in Inflation regime
# DashboardPlot(Inputs=dashboard.list,
# 	            Fund='Global.Macro',
# 	            Regime='Inflation')

#Dashboard for all regimes/fund combinations
for(i in fund.names){
	for(j in regimes.names){
		DashboardPlot(Inputs=dashboard.list,
			            Fund=i,
		              Regime=j
		)
	}
}
########################################################################################

l
This results in n.funds(13)*n.regimes(4)= 52 pdf files, each representing a unique combination of a particular fund and a particular regime.

I have uploaded the pdf files into wordpress and they should be available for viewing by clicking on the links in the following table :

l

Funds | Regimes Equity Currency Inflation Growth
Convertible Arbitrage C&EQ C&FX C&Inf C&G
CTA Global CTA&EQ CTA&FX CTA&Inf CTA.&G
Distressed Securities D&EQ D&FX D&Inf D&G
Emerging Markets Em&EQ Em&FX Em&Inf Em&G
Equity Market Neutral Eq&EQ Eq&FX Eq&Inf Eq&G
Event Driven Ev&EQ Ev&FX Ev&Inf Ev&G
Fixed Income Arbitrage Fix&EQ Fix&FX Fix&Inf Fix&G
Funds of Funds Fun&EQ Fun&FX Fun&Inf Fun&G
Global Macro Glo&EQ Glo&FX Glo&Inf Glo&G
Long Short Equity Lon&EQ Lon&FX Lon&Inf Lon&G
Merger Arbitrage Me&EQ Me&FX Me&Inf Me&G
Relative Value Rel&EQ Rel&FX Rel&Inf Rel&G
Short Selling Sh&EQ Sh&FX Sh&Inf Sh&G

l

This html table was generated using The Tables Generator with Compact mode ticked and Do not Generate CSS unticked.

To give an example of the information contained in a dashboard,consider the case of  Global Macro x Inflation and relevant points of interest below :

l

mkl

l

Points of interest :

[Barcharts & Table of estimates in top rows]

  1. The Inflation economic regime variable begins in the normal state,with a high probability of remaining in this state (Persistence : 90.1%) and a low probability of transitioning to the event state (Transition : 9.9%).The estimated mean and sigma for this normal state are 0.2559 and 0.275 respectively.
  2. The Event state for all regime variable is characterised by higher estimated sigmas and means than corresponding values for the normal state.
  3. The persistence of normal states across all economic regime variables is higher than that that of event states. The corollary is found in transition estimates where we are more likely to transition from the event state to the normal state than from the normal state to the event state for all regime variables.
  4. The dodged barchart in the final column seeks to answer the question : What is the average return for fund x when regime y is in state a or state b? For our example, the mean return to the Global Macro strategy when the inflation regime is in the normal state exceeds the corresponding value in the event state.

 

[Selection & Performance plots in bottom rows]

  1. The list of funds and regime variables are given as strips of text in two separate ggplots with a green rectangle beneath the chosen fund and a red rectangle beneath the chosen regime. They are here mainly for eye candy, filling in space and horsing around with the ggplot package!
  2. The monthly,cumulative returns and drawdowns are also drawn.It seems that that for the Global Macro x Inflation case,knowledge of Markov states does not translate to higher cumulative returns.The drawdowns however seem more favourable.
  3. We also have a time series plot of the chosen economic regime variable along with an overlay of event versus non-event states.The thin red bars signify moments in time where the regime variable is in the event state.The blue bars show when that variable is in its normal state.There are some strange discrepancies between the pdf file and this zoomed image above.The pdf file is a better option here even though the colour is strangely off depending on the zoom.Also there is some sort of white taper in the pdf version of this plot which is probably due to the jagged lines I mentioned in the previous post.
  4. The final barchart in the bottom right corner shows the scaled difference in means across states for the chosen regime variable. If I understood this correctly,this should answer the question : By how many standard deviations is a particular fund’s performance higher/lower during the event/normal regime? In our case,the Global Macro strategy performs worse in the event regime than the normal regime by -0.10 standard deviations (?) I suppose this corroborates the previous finding that the mean performance of the chosen fund,when seen in the context of the inflation regime variable,is superior in the normal state versus the event state.

Whether these results are consistent or my interpretations correct I do not know.Any mistakes are of course mine.Although the code for the function is a bit of a mess and for the most part just arranges the plots in the desired layout,I include it here for the sake of completeness :
l

#########################################################################################
# Dashboard
#########################################################################################

DashboardPlot <- function(Inputs,Fund,Regime){

#Set up data
	regime.map <- data.frame(regime=c('Equity','Currency','Inflation','Growth'),idx=c(1:4),stringsAsFactors=F)
	fund.names <- names(Inputs$FundRet.list)
  reg.names <- regime.map[,1]

#regime strip
  lbl.font <- c(rep('bold',4))
	reg.df<-data.frame(lab.x=rep(20,4),lab.y=seq(2,26,length=4),lab.txt=reg.names,lab.font=lbl.font,stringsAsFactors=F)
	empty.df <- data.frame(x=-2:26,y=-2:26)
	ymin <- filter(reg.df,lab.txt==Regime)$lab.y-2
	ymax<-ymin+diff(reg.df$lab.y)[1]-1.5
	rect.df <- data.frame(xmin=-4,xmax=28,ymin=ymin,ymax=ymax)

	gg.reg.list <-  ggplot(empty.df)+geom_blank()+theme_economist(base_size=5)+
                      theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(),
   	 									panel.border = element_blank(), panel.background = element_blank(), axis.title.x = element_blank(),
    									axis.title.y = element_blank(), axis.text.x = element_blank(), axis.text.y = element_blank(),
    									axis.ticks = element_blank(),legend.position='none')+
                   		geom_rect(data=rect.df,aes(xmin = xmin, xmax = xmax, ymin = ymin, ymax = ymax),fill='red',alpha=0.4)+
		                  geom_text(data=reg.df,aes(x=lab.x,y=(lab.y),label=lab.txt,fontface=lab.font,hjust=1,vjust=0),size=rel(3))

#fund strip
  lbl.font <- c(rep('bold',13))
	fund.df<-data.frame(lab.x=rep(20,13),lab.y=seq(0,26,length=13),lab.txt=fund.names,lbl.font=lbl.font,stringsAsFactors=F)
	empty.df <- data.frame(x=-2:26,y=-2:26)
	ymin <- filter(fund.df,lab.txt==Fund)$lab.y-0.5
	ymax<-ymin+diff(fund.df$lab.y)[1]-0.5
	rect.df <- data.frame(xmin=-4,xmax=28,ymin=ymin,ymax=ymax)
	gg.fund.list <-  ggplot(empty.df)+geom_blank()+theme_economist(base_size=5)+
                      theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(),
   	 									panel.border = element_blank(), panel.background = element_blank(), axis.title.x = element_blank(),
    									axis.title.y = element_blank(), axis.text.x = element_blank(), axis.text.y = element_blank(),
    									axis.ticks = element_blank(),legend.position='none')+
                   		geom_rect(data=rect.df,aes(xmin = xmin, xmax = xmax, ymin = ymin, ymax = ymax),fill='green',alpha=0.4)+
		                  geom_text(data=fund.df,aes(x=lab.x,y=(lab.y),label=lab.txt,fontface=lbl.font,hjust=1,vjust=0),size=rel(3))

#Data table
	reg.str <- paste(filter(regime.map,regime==Regime)[2],']]',sep='')
	tbl1 <- Inputs$Turbulence1.tbl
	tbl2 <- Inputs$Turbulence2.tbl
	initial.state <- tbl1[unlist(filter(regime.map,regime==Regime)[2]),1]

	mean.s1 <- round(tbl1[unlist(filter(regime.map,regime==Regime)[2]),4],3)
        mean.s2 <- round(tbl2[unlist(filter(regime.map,regime==Regime)[2]),4],3)
	mean.df <- data.frame(State=c('Normal','Event'),Mean=c(mean.s1,mean.s2),stringsAsFactors=F)

	sd.s1 <- round(tbl1[unlist(filter(regime.map,regime==Regime)[2]),5],3)
        sd.s2 <- round(tbl2[unlist(filter(regime.map,regime==Regime)[2]),5],3)
	sd.df <- data.frame(State=c('Normal','Event'),StDev=c(sd.s1,sd.s2),stringsAsFactors=F)

	#Check State and store fitted information
	if(initial.state=='State 1'){
		temp <-Inputs$Turbulence1.tbl
		state.type <-'Normal'
	}else{
		temp<-Inputs$Turbulence2.tbl
    state.type <- 'Event'
	}

	pers <- temp[unlist(filter(regime.map,regime==Regime)[2]),2]
	trans <- temp[unlist(filter(regime.map,regime==Regime)[2]),3]
	mean.f <- round(temp[unlist(filter(regime.map,regime==Regime)[2]),4],3)
	stddev <- round(temp[unlist(filter(regime.map,regime==Regime)[2]),5],3)

#Create Plots

	gg.ret.plot <- eval(parse(text=paste('Inputs$FundRet.list$',Fund,sep='')))
	gg.turb.plot <- eval(parse(text=paste('Inputs$Turbulence.plot[[',reg.str,sep='')))

	gg.turb.map1 <- eval(parse(text=paste('Inputs$Turbulence.map2[[',reg.str,'[[1]]',sep='')))
	gg.turb.map2 <- eval(parse(text=paste('Inputs$Turbulence.map2[[',reg.str,'[[2]]',sep='')))

	gg.avg.plot <- eval(parse(text=paste('Inputs$plot.list[[',reg.str,sep='')))
	gg.cumu.plot1 <- eval(parse(text=paste('Inputs$CompareCumu$',Fund,'[[',reg.str,'$plot[[1]]',sep='')))
	gg.cumu.plot2 <- eval(parse(text=paste('Inputs$CompareCumu$',Fund,'[[',reg.str,'$plot[[2]]',sep=''))) 

	gg.facet.plot <- eval(parse(text=paste('Inputs$FacetPlot$',Fund,sep='')))
	gg.mkv.mean <- Inputs$MarkovPlots[[1]]
	gg.mkv.sd <- Inputs$MarkovPlots[[3]]
	gg.mkv.p <- Inputs$MarkovPlots[[5]]
	gg.mkv.t <- Inputs$MarkovPlots[[6]]

#Dashboard table
	lab1 <- paste("~Fund   :  ",Fund,sep=' ')
	lab2 <- paste('~Regime   :  ',Regime,sep=' ')
        lab3 <- paste('       |__Initial State :  ',state.type,sep=' ')
        lab4 <- paste('       |__Persistence :  ',pers,sep=' ')
        lab5 <- paste('       |__Transition :  ',trans,sep=' ')
        lab6 <- paste('       |__Mean  :  ',mean.f,sep=' ')
        lab7 <- paste('       |__StDev  :  ',stddev,sep=' ')
        lbl.cols <- c('green','red','white','white','white','white','white')
        lbl.font <- c('bold','bold','bold','bold','bold','bold','bold')
        lbl.df <-data.frame(lab.x=rep(-2,7),lab.y=seq(-1,8,length=7),lab.txt=c(lab1,lab2,lab3,lab4,lab5,lab6,lab7),lbl.font=lbl.font,lbl.cols=lbl.cols,stringsAsFactors=F)

	empty.df <- data.frame(x=-2:12,y=-2:12)
	gg.empty.plot <-  ggplot(empty.df,aes(x=x,y=y))+geom_blank()+geom_hline(y=11,colour='white',size=1)+theme_tufte(base_size=5)+
                      theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(),
   	 									panel.border = element_blank(), panel.background = element_blank(), axis.title.x = element_blank(),
    									axis.title.y = element_blank(), axis.text.x = element_blank(), axis.text.y = element_blank(),
    									axis.ticks = element_blank(),plot.title = element_text(size = rel(4),colour='white'),plot.background = element_rect(colour = "skyblue4", fill = "skyblue4"))+
		                  labs(title='\n.:: Dashboard ::.')+
		                  geom_text(data=lbl.df,aes(x=lab.x,y=rev(lab.y),label=lab.txt,hjust=0,vjust=1,fontface=lbl.font),size=rel(3),colour=lbl.cols)

	gg.empty2.plot <-  ggplot(empty.df,aes(x=x,y=y))+geom_blank()+theme_tufte(base_size=5)+
                      theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(),
   	 									panel.border = element_blank(), panel.background = element_blank(), axis.title.x = element_blank(),
    									axis.title.y = element_blank(), axis.text.x = element_blank(), axis.text.y = element_blank(),
    									axis.ticks = element_blank(),plot.title = element_text(size = rel(4)),plot.background = element_rect(colour = "skyblue4", fill = "skyblue4"))

#Data table
  temptbl <- eval(parse(text=paste('Inputs$FacetObj$',Fund,'$df',sep='')))
	gg.tbl.plot <- QuickTbl(temptbl,title='Averaging Returns')

#Markov estimate tables
	gg.mkv.mean.tbl <- Inputs$MarkovPlots[[2]]
   	gg.mkv.sd.tbl <- Inputs$MarkovPlots[[4]]
	gg.mkv.p.tbl <- Inputs$MarkovPlots[[7]]
        gg.mkv.t.tbl <- Inputs$MarkovPlots[[8]]

#Arrange multiple plots and save as pdf
pdf(file = paste(Fund,'&',Regime,'.pdf',sep=''), width = 30, height = 17)
  layOut(
  	     list(gg.empty.plot,1:2,1),
  	     list(gg.empty2.plot,3:6,1),
  	     list(gg.reg.list,3,2),
  	     list(gg.fund.list,4:6,2),
   	  	 list(gg.mkv.mean,1,2),
  	     list(gg.mkv.sd,1,3),
  	     list(gg.mkv.p,1,4),
  	     list(gg.mkv.t,1,5),
  	  	 list(gg.facet.plot,1,6),
  	     list(gg.mkv.mean.tbl,2,2),
  	     list(gg.mkv.sd.tbl,2,3),
         list(gg.mkv.p.tbl,2,4),
  	     list(gg.mkv.t.tbl,2,5),
  	     list(gg.tbl.plot,2,6),
  	     list(gg.ret.plot,3:4,3:4),
  	     list(gg.cumu.plot1,5,3:4),
   	     list(gg.cumu.plot2,6,3:4),
  	     list(gg.avg.plot,5:6,5:6),
  	  	 list(gg.turb.plot,3:4,5:6)
   	)
dev.off()
}
#########################################################################################

l
Wish I could make this collapsable but wordpress seems to need a plugin for that.

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.

In this series of posts I intend to replicate some of the methodologies applied by Kritzman et al.(2012) in their paper  Regime Shifts: Implications for Dynamic Strategies – , wherein :

  • a two-state hidden Markov model is fitted to observed (or transformed) measures of financial turbulence,inflation and economic growth to partition history into meaningful regimes .
  • the in-sample performance of a variety of risk premia is probed by comparing scaled differences in respective state means across regime variables.
  • the out-of-sample performance of risk premia is examined by [1] tilting portfolios dynamically according to the likelihood of a particular regime and by [2] comparing the results of asset allocation that is conscious of regime switching versus the static alternative.

The excellent (but sadly idle) Quantivity blog has replicated some parts of the paper with much more rigour than I am capable of. This is perhaps as good a time as any to mention some differences between my post,the original paper and Quantivity’s post (QP).

One point of departure between my and QP’s post relates to the library/package used to estimate the hidden markov model. While QP used the RHmm package, I can only use the depmixS4 package (RHmm seems to be incompatible with newer versions of R…at least for me). As far as I can tell,the depmixS4 package does not support forecasting which means that I will probably not be able to take a gander at the out-of-sample issues alluded to above.

Another difference relates to the risk premia used in the original paper versus my post.

ko

Instead of using these data series (shamelessly snipped from the paper), I simply downloaded a revised edhec data set (which comes with the PerformanceAnalytics package) from the edhec-risk-homepage : http://www.edhec-risk.com/indexes/pure_style

These are simply returns on various hedge fund strategies as opposed to deliberately constructed risk premia that are expected to respond to different states across economic/financial regime variables. A snapshot of the strategies follows, along with their most recent summary results :

nmop

 


 

[Data]

In terms of the data used in this project, it is useful to make a distinction at the outset between economic regime variables which are used to partition history into different states in a meaningful way) and risk premia (simple returns in my case) which are used to assess performance across states for each economic regime variable.

 

–Economic Regime Variables–

[1] Financial Market Turbulence

The paper defines this variable as a condition in which asset prices behave in an unusual fashion given historical patterns of behaviour and proposes to measure it using the Mahalanobis distance :

cxFinancial market turbulence comes in two flavours : Equity and Foreign exchange.

The Equity Turbulence Index (ETI) was constructed by [1] downloading sp500 daily returns,[2] estimating mean and covariances using equally weighted historical returns for the past 10 years and [3] averaging the daily turbulence scores within each month. The Currency Turbulence Index (CTI) was constructed by [1] downloading daily currency returns (vs.USD) for 8 of the G-10 countries, [2] estimating mean and covariances using equally weighted historical returns for the past 3 years and [3] averaging the daily turbulence scores within each month.

The code for downloading relevant time series follows :

######################[a] Economic Regime Variables######################################

#[a.1] Equity Turbulence
symbol <- c('^GSPC')
beg <- as.Date('1970/01/01')
end <- as.Date('2014/01/01')

getSymbols(symbol,from=beg,to=end)
SPX.ret <- ROC(eval(parse(text=paste(substr(symbol,start=which(strsplit(symbol,split='^')[[1]]=='^')+1,stop=nchar(symbol)),'$GSPC.Close',sep=''))),n=1)
attr(SPX.ret,'type')='Equity'

#[a.2] Currency Turbulence
G9curr <- c('DEXUSAL','DEXUSUK','DEXCAUS','DEXNOUS','DEXJPUS','DEXUSNZ','DEXSDUS','DEXSZUS')
col.xts <- c('AUS','GBP','CAD','NOR','JAP','NZD','SWD','CHF')
inverted.idx <- grep('DEXUS',G9curr,value=F)
normal.idx <- index(G9curr)[-inverted.idx]

getSymbols(G9curr,src='FRED')
fx.xts <- do.call(merge,lapply(G9curr,function(x) eval(parse(text=x))))
names(fx.xts) <- col.xts

for(i in inverted.idx) fx.xts[,i] <- 1/eval(parse(text=G9curr[i]))
fx.xts <- ROC(fx.xts)
attr(fx.xts,'type')='Currency'

lI

[2] Inflation

Inflation was measured using the Consumer Price Index For All Urban Consumers :All Items, available from the FRED website and directly accessible via r :

#[a.3] Inflation Turbulence
getSymbols('CPIAUCNS',src='FRED')
inf.xts <- ROC(CPIAUCNS['1946-01-01::',1],n=1)
attr(inf.xts,'type')='Inflation'

l
[3] Economic Growth

Economic growth was measured using quarterly percentage growth in Real Gross Domestic Product :

#[a.4] Economic growth Turbulence
getSymbols('GDPC1',src='FRED')
ecg.xts <- ROC(GDPC1['1947-01-01::',1],n=1)
attr(ecg.xts,'type')='Growth'

ll
Now that we have downloaded the necessary time series for our regime variables,let’s actually generate the Mahalanobis distances for the ETI and CTI. As far as I can tell, Inflation and economic growth variables are already in the desired form (simple returns), so there is no need for further transformations. After storing the regime variables in a list object, I apply a custom function to each element of that list :

#[1] Store return data in list object

	regime.variables<-list(Equity=SPX.ret,Currency=fx.xts,Inflation=inf.xts,Growth=ecg.xts)

#[2] Generate Turbulence indices where required (we focus on the final complete turbulence index where available)

	Turbulence.list <- lapply(regime.variables,function(x) GenMahalanobis(x)$complete.turbulence)

l
The GenMahalanobis function detects which regime variable is being passed to it,generates daily turbulence indices where appropriate and averages the result within each month :

#########################################################################################
# Rolling Mahalanobis function
#########################################################################################

GenMahalanobis <- function(inp.ret){
#Set up variables	

	#inp.ret <- inf.xts
	ret.xts <- na.omit(inp.ret)
	n.var <- dim(ret.xts)[2]
	map <- data.frame(Type=c('Equity','Currency'),Window=c(2500,750))
	estwind <- map[map$Type==attributes(ret.xts)$type,2]
	mahalanobis.list <- NULL

#Check which type of variable we have
	cond.eval <- eval(parse(text='attributes(ret.xts)$type'))
	if(cond.eval=='Equity'|| cond.eval=='Currency'){

		#Looping across variables and windows & storing results in list elements
			for(i in 1:n.var){
				window.beg <- 1
				window.idx <- 1
				window.end <- estwind

				mahalanobis.list[[i]] <- list()
				names(mahalanobis.list)[i] <- colnames(ret.xts)[i]

				temp <- NULL
				mahalanobis.complete <- NULL

				while(window.end<=dim(ret.xts)[1]){
		  		mahalanobis.list[[i]][[window.idx]] <- list()
					names(mahalanobis.list[[i]])[[window.idx]] <- sprintf('Window%i',window.idx)

					window.data <- eval(parse(text='window.beg:window.end'))
					temp <- mahalanobis(x=ret.xts[window.data,i],center=mean(ret.xts[window.data,i]),cov=cov(ret.xts[window.data,i]))

					mahalanobis.list[[i]][[window.idx]] <- xts(temp,order.by=as.Date(names(temp)))

    			window.beg = window.beg+estwind
					window.end = window.end+estwind
    			window.idx = window.idx+1
				}
			mahalanobis.list[[i]]$xts <- do.call(rbind,mahalanobis.list[[i]])
		}

	#Collect in single xts structure
		mahalanobis.list$distances <- NULL
		for(i in 1:n.var){
			mahalanobis.list$distances <- cbind(mahalanobis.list$distances,mahalanobis.list[[i]]$xts)
			colnames(mahalanobis.list$distances)[i] <-  colnames(ret.xts)[i]
		}

	#Average across variables & calculate mean within each month
		mahalanobis.list$daily.turbulence <- NULL
		mahalanobis.list$monthly.turbulence <- NULL
		mahalanobis.list$daily.turbulence <- xts(rowMeans(mahalanobis.list$distances),order.by=index(mahalanobis.list$distances))
		mahalanobis.list$complete.turbulence <- apply.monthly(x=mahalanobis.list$daily.turbulence,FUN=mean)
  	names(mahalanobis.list$daily.turbulence) <- names (mahalanobis.list$complete.turbulence) <- attributes(ret.xts)$type

}else{

	mahalanobis.list$complete.turbulence <- xts(mahalanobis(x=ret.xts,center=mean(ret.xts),cov=cov(ret.xts)),order.by=index(ret.xts))
  names(mahalanobis.list$complete.turbulence) <- attributes(ret.xts)$type

}

#return list object
return(mahalanobis.list)
}

#########################################################################################

l

— Risk Premia —

For this category of variables,I simply downloaded the updated csv file for the edhec data set,imported it into r and stored a ggplot of each strategy’s return as an element of a list object :

#[7.1] Load and date the updated edhec dataset & store fund returns in list

orig.csv <- read.csv('edhecfinal.csv',header=T,stringsAsFactors=F)
edhec.xts <- xts(order.by=seq(as.Date('1997-01-01'),as.Date('2014-07-01'),by='1 month'),orig.csv[,-1])

n.funds <- ncol(edhec.xts)
fund.names <- names(edhec.xts)
temp.df <- data.frame(TimeIdx=index(edhec.xts),coredata(edhec.xts))

edhec.ret <- list()
	for(i in 1:n.funds){
		edhec.ret[[i]] <- list()
		names(edhec.ret)[i] <- fund.names[i]
		edhec.ret[[i]] <-EdhecRet(edhec.xts[,i])
	}

l
The EdhecRet function simply returns a ggplot of time series returns for the chosen strategy. The reason for storing these plots as elements of a list object is to simplify access to fund/regime results in the final dashboard of results. For example, to access the monthly return plot for Short.Selling strategy,edhec.ret$Short.Selling yields :

kk

 

The code for that function is simply :

#########################################################################################
# Simple Fund Returns
#########################################################################################

EdhecRet <- function(xts.ret){
	loc <- environment()
	temp.df <- NULL
	temp.df <- data.frame(TimeIdx=as.Date(index(xts.ret)),coredata(xts.ret),stringsAsFactors=F)
	pp <- ggplot(data=temp.df,aes(x=TimeIdx,y=eval(parse(text=colnames(temp.df)))),environment=loc)+
						      theme_tufte(base_size=11)+
			            geom_line(size=0.5,color='skyblue4') +
			            labs(y='Values',title='Monthly Returns',x='')+
			        		theme(legend.position = "none")
	return(pp)
}
#########################################################################################

l
Now that we have downloaded or generated the necessary data series, the next posts shall deal with [1] fitting a 2 state hidden markov model to the economic regime variables and [2] assessing in-sample performance of different strategies across the hidden states for each regime variable. Ultimately I wish to present a dashboard of results for a chosen fund/regime combination.