Create a gist now

Instantly share code, notes, and snippets.

What would you like to do?
# Forecasting membership in the NRA for June 2013
# Trey Causey (tcausey@uw.edu)
# 26 February 2013
require(forecast)
# forecast package by Rob Hyndman
# http://robjhyndman.com/hyndsight/
data <- read.csv("nra_members.csv")
# Convert to time series
data_ts <- ts(data$nra_subscribers, frequency = 12, start = c(1981, 7), end = c(2012, 12))
# Automatically fit Error, Trend, Seasonality exponential smoothing model
fitted_ets <- ets(data_ts)
fcast_ets <- forecast(fitted_ets, h = 6) # Only forecast through June 2013
png(file = "ets_plot_full.png")
plot(fcast_ets) # Full time series
dev.off()
png(file = "ets_plot_2009.png")
plot(fcast_ets, 48) # Only 2009-2013
dev.off()
fcast_ets
# Point Forecast Lo 80 Hi 80 Lo 95 Hi 95
# Jan 2013 3098650 3049788 3147511 3023923 3173376
# Feb 2013 3096501 3010395 3182607 2964813 3228189
# Mar 2013 3094515 2976525 3212504 2914065 3274964
# Apr 2013 3092679 2944467 3240891 2866008 3319350
# May 2013 3090982 2913277 3268686 2819206 3362757
# Jun 2013 3089413 2882621 3296206 2773152 3405675
# Automatically fit ARIMA model
fitted_arima <- auto.arima(data_ts)
fcast_arima <- forecast(fitted_arima, h = 6) # Only forecast through June 2013
plot(fcast_arima) # Plot with uncertainty included
png(file = "arima_plot_full.png")
plot(fcast_arima) # Full time series
dev.off()
png(file = "arima_plot_2009.png")
plot(fcast_arima, 48) # Only 2009-2013
dev.off()
fcast_arima
# Point Forecast Lo 80 Hi 80 Lo 95 Hi 95
# Jan 2013 3090481 3046323 3134639 3022947 3158015
# Feb 2013 3082030 3008425 3155635 2969460 3194599
# Mar 2013 3078289 2979934 3176644 2927868 3228711
# Apr 2013 3066746 2944699 3188793 2880092 3253401
# May 2013 3073790 2928591 3218988 2851728 3295852
# Jun 2013 3076516 2908441 3244592 2819467 3333566
# ARIMA has lower forecast and narrower 95% CI
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment