Holt-Winters forecast using ggplot2

July 16, 2012
By

(This article was first published on FishyOperationsCategory Archives: r, and kindly contributed to R-bloggers)

R has great support for Holt-Winter filtering and forecasting. I sometimes use this functionality, HoltWinter & predict.HoltWinter, to forecast demand figures based on historical data. Using the HoltWinter functions in R is pretty straightforward.

Let's say our dataset looks as follows;

demand <- ts(BJsales, start = c(2000, 1), frequency = 12)
plot(demand)

plot of chunk unnamed-chunk-1

Now I pass the timeseries object to HoltWinter and plot the fitted data.

hw <- HoltWinters(demand)
plot(hw)

plot of chunk unnamed-chunk-2

Next, we calculate the forecast for 12 months with a confidence interval of .95 and plot the forecast together with the actual and fitted values.

forecast <- predict(hw, n.ahead = 12, prediction.interval = T, level = 0.95)
plot(hw, forecast)

plot of chunk unnamed-chunk-3

As you can see, this is pretty easy to accomplish. However, as I use ggplot2 to visualize a lot of my analyses, I would like to be able to do this in ggplot2 in order to maintain a certain uniformity in terms of visualization.

Therefore, I wrote a little function which extracts some data from the HoltWinter and predict.HoltWinter objects and feeds this to ggplot2;

#HWplot.R

library(ggplot2)
library(reshape)


HWplot<-function(ts_object,  n.ahead=4,  CI=.95,  error.ribbon='green', line.size=1){
	
	hw_object<-HoltWinters(ts_object)
	
	forecast<-predict(hw_object,  n.ahead=n.ahead,  prediction.interval=T,  level=CI)
	
	
	for_values<-data.frame(time=round(time(forecast),  3),  value_forecast=as.data.frame(forecast)$fit,  dev=as.data.frame(forecast)$upr-as.data.frame(forecast)$fit)
	
	fitted_values<-data.frame(time=round(time(hw_object$fitted),  3),  value_fitted=as.data.frame(hw_object$fitted)$xhat)
	
	actual_values<-data.frame(time=round(time(hw_object$x),  3),  Actual=c(hw_object$x))
	
	
	graphset<-merge(actual_values,  fitted_values,  by='time',  all=TRUE)
	graphset<-merge(graphset,  for_values,  all=TRUE,  by='time')
	graphset[is.na(graphset$dev),  ]$dev<-0
	
	graphset$Fitted<-c(rep(NA,  NROW(graphset)-(NROW(for_values) + NROW(fitted_values))),  fitted_values$value_fitted,  for_values$value_forecast)
	
	
	graphset.melt<-melt(graphset[, c('time', 'Actual', 'Fitted')], id='time')
	
	p<-ggplot(graphset.melt,  aes(x=time,  y=value)) + geom_ribbon(data=graphset, aes(x=time, y=Fitted, ymin=Fitted-dev,  ymax=Fitted + dev),  alpha=.2,  fill=error.ribbon) + geom_line(aes(colour=variable), size=line.size) + geom_vline(x=max(actual_values$time),  lty=2) + xlab('Time') + ylab('Value') + opts(legend.position='bottom') + scale_colour_hue('')
	return(p)

}

The above script is saved in a file called HWplot.R. If I load this file from R – via source() – I can directly call the function HWplot. The HWplot can be called as follows:

HWplot(ts_object, n.ahead=4, CI=.95, error.ribbon='green',line.size=1)

HWplot takes the following arguments;

  • ts_object: the timeseries data
  • n.ahead: number of periods to forecast
  • CI: confidence interval
  • error.ribbon: colour of the error ribbon
  • line.size: size of the lines
source("HWplot.R")
demand <- ts(BJsales, start = c(2000, 1), frequency = 12)
HWplot(demand, n.ahead = 12)

plot of chunk unnamed-chunk-4

It's also very easy to adjust the graph after it is returned by the function;

graph <- HWplot(demand, n.ahead = 12, error.ribbon = "red")

# add a title
graph <- graph + opts(title = "An example Holt-Winters (gg)plot")

# change the x scale a little
graph <- graph + scale_x_continuous(breaks = seq(1998, 2015))

# change the y-axis title
graph <- graph + ylab("Demand ($)")

# change the colour of the lines
graph <- graph + scale_colour_brewer("Legend", palette = "Set1")

# the result:
graph

plot of chunk unnamed-chunk-5

The HWplot R code: HWplot.R

The post Holt-Winters forecast using ggplot2 appeared first on FishyOperations.

To leave a comment for the author, please follow the link and comment on his blog: FishyOperationsCategory Archives: r.

R-bloggers.com offers daily e-mail updates about R news and tutorials on topics such as: visualization (ggplot2, Boxplots, maps, animation), programming (RStudio, Sweave, LaTeX, SQL, Eclipse, git, hadoop, Web Scraping) statistics (regression, PCA, time series, trading) and more...



If you got this far, why not subscribe for updates from the site? Choose your flavor: e-mail, twitter, RSS, or facebook...

Tags: , , , , ,

Comments are closed.