Day #325 - Time Series Forecasting Models

Models learn below features from input dataset and predict accordingly

• Trend
• Seasonality
• Smoothing Factor

Models

• Moving Average - Leverage recent weeks data
• Weighted Moving Average - Give varying weights for recent data
• Single Exponential - Smoothing factor added in the cost function
• Double Exponential - Smoothing & Trend
• Triple Exponential - Smoothing, Trend & Seasonality
• ARIMA - Auto Regression Integrated Moving Average
• Linear & Polynomial Models - Based on Features collected

	#Sample Data
	#Generate sample data of 100 records with mean = 850 and standard deviation = 900
	myvector <- rnorm(1000, 850, 900)
	myvector
	#Split it into timeseries weekly data for 104 weeks
	myts = ts(myvector, start=c(2014,1),end=c(2016,1),frequency=104)
	myts
	myvector <- rnorm(1000, 850, 900)
	myvector
	myts = ts(myvector, start=c(2014,1),end=c(2016,1),frequency=104)
	library(forecast)
	fit = HoltWinters(myts,beta = FALSE, gamma = FALSE)
	skirtsseriesforecasts2 <- forecast.HoltWinters(fit, h=19)
	plot.forecast(skirtsseriesforecasts2)
	skirtsseriesforecasts2
	hist(skirtsseriesforecasts2$residuals)
	fit$alpha
	fit$beta
	fit$gamma
	fit$seasonal
	#Data Plotting Analysis
	plot(myts)
	#Moving Average
	sm = ma(myts, order=4)#4 week average
	plot(sm)
	#Moving average does not reflect all the dips / heights. It averages out the pattern
	myvector <- rnorm(1000, 850, 900)
	myvector
	myts = ts(myvector, start=c(2014,1),end=c(2016,1),frequency=104)
	sm = ma(myts, order=4)#4 week average
	p <- predict(sm, 50, prediction.interval = TRUE)
	plot(fit, p, typ1 = 'l', xlab = 'Weeks', ylab = 'Replenishment Count', main='Replenishmely Weekly Analysis - MA')
	#Plotting with Single Exponential
	#It is common in simple exponential smoothing to use the first value in the time series as the initial value for the level.
	model = HoltWinters(myts,beta = FALSE, gamma = FALSE, l.start=23.56)
	plot(model)
	p <- predict(model, 50, prediction.interval = TRUE)
	plot(fit, p, typ1 = 'l', xlab = 'Weeks', ylab = 'Replenishment Count', main='Replenishmely Weekly Analysis')
	#Plotting with Double Exponential
	#Trend is considered in the predictions
	#This has higher accuracy than previous single exponential smoothing plot
	myvector <- rnorm(1000, 850, 900)
	myvector
	myts = ts(myvector, start=c(2014,1),end=c(2016,1),frequency=104)
	myts
	fit = HoltWinters(myts, gamma = FALSE)
	plot(fit)
	p <- predict(fit, 10, prediction.interval = TRUE)
	plot(fit, p, typ1 = 'l', xlab = 'Weeks', ylab = 'Replenishment Count', main='Replenishmely Weekly Analysis')
	p
	#Plotting with Triple Exponential
	fit = HoltWinters(myts, seasonal = "additive")
	myvector <- rnorm(1000, 850, 900)
	myvector
	myts = ts(myvector, start=c(2014,1),end=c(2016,1),frequency=104)
	fit = HoltWinters(myts, seasonal = "additive")
	p <- predict(fit, 10, prediction.interval = TRUE)
	plot(fit, p, typ1 = 'l', xlab = 'Weeks', ylab = 'Replenishment Count', main='Replenishmely Weekly Analysis')
	myvector <- rnorm(1000, 850, 900)
	myvector
	myts = ts(myvector, start=c(2014,1),end=c(2016,1),frequency=104)
	fit = HoltWinters(myts, seasonal = "multiplicative")
	p <- predict(fit, 10, prediction.interval = TRUE)
	plot(fit, p, typ1 = 'l', xlab = 'Weeks', ylab = 'Replenishment Count', main='Replenishmely Weekly Analysis')
	#Trend, Seasonality and smoothing factors are considered in Triple Exponential smoothing
	#Sample data predictions in red
	#ARIMA
	arimafit = arima(myts, order = c(0,0,0))
	summary(arimafit)
	plot(myts, xlim=c(2014,2016),lw=2,col="blue")
	lines(predict(arimafit,n.ahead=50)$pred,lw=2,col="red")
	#Auto-Regression
	arfit = ar(myts)
	summary(arfit)
	pred = predict(arfit,n.ahead=30)
	plot(myts,type="l",xlim=c(2014,2016),ylim=c(100,2200),xlab="weeks",ylab="forecast")
	lines(pred$pred,col="red")
	#http://www.statmethods.net/advstats/timeseries.html
	#numeric vector
	#rnorm(n, mean = , sd = )
	myvector <- rnorm(1000, 850, 900)
	myvector
	myts = ts(myvector, start=c(2009,1),end=c(2014,12),frequency=12)
	#The ts() function will convert a numeric vector into an R time series object.
	#The format is ts(vector, start=, end=, frequency=) where
	#start and end are the times of the first and last observation and
	#frequency is the number of observations per unit time (1=annual, 4=quartly, 12=monthly, etc.).
	myts
	myts2 = window(myts,start=c(2014,6), end=c(2014,12))
	#window(x, start = NULL, end = NULL,
	#frequency = NULL, deltat = NULL, extend = FALSE, ...)
	#start - the start time of the period of interest.
	#end - the end time of the period of interest.
	#window is a generic function which extracts the subset of the object x observed between the times start and end
	#plot(myts)
	acf(myts)
	pacf(myts)
	#seasonal decomposition
	fit = stl(myts,s.window = "period")
	#plot(fit)
	monthplot(myts)
	library(forecast)
	seasonplot(myts)
	#simple moving average
	#http://r-statistics.co/Time-Series-Forecasting-With-R.html
	sm = ma(myts, order=12)#12 month average
	lines(sm,col="red")
	#simple exponential mode
	#Single Exponential smoothing weights past observations with exponentially decreasing weights to forecast future values
	#Does not do well when there is a trend in the data
	fit = HoltWinters(myts,beta = FALSE, gamma = FALSE)
	plot(fit)
	#double exponential
	#Double exponential smoothing uses two constants and is better at handling trends
	#Smoothing is to introduce a term to take into account the possibility of a series exhibiting some form of trend
	#This slope component is itself updated via exponential smoothing
	fit = HoltWinters(myts, gamma = FALSE)
	plot(fit)
	#triple exponential
	#Triple exponential smoothing takes into account seasonal changes as well as trends
	#There are different types of seasonality: 'multiplicative' and 'additive' in nature
	fit = HoltWinters(myts)
	fit1 = HoltWinters(myts, seasonal = "multiplicative")
	fit2 = HoltWinters(myts, seasonal = "additive")
	plot(fit)
	plot(fit1)
	plot(fit2)
	library(forecast)
	forecast(fit,3)
	plot(forecast(fit,3))
	#http://stats.stackexchange.com/questions/144158/daily-time-series-analysis
	ets(myts)
	# ETS(Error,Trend,Seasonal)
	fit = tbats(myts)
	seasonal = !is.null(fit$seasonal)
	seasonal
	#weekly seasonality
	#http://stats.stackexchange.com/questions/144158/daily-time-series-analysis
	x.msts = msts(myts,seasonal.periods = c(7,365.25))
	model = tbats(x.msts)
	plot(forecast(model,h=100))
	#You should not use arima() or auto.arima(),
	#since these can only handle a single type of seasonality: either weekly or yearly

view raw TimeseriesExamples.R hosted with ❤ by GitHub

Happy Learning!!!