require("fpp2")
## Loading required package: fpp2
## Warning: package 'fpp2' was built under R version 3.4.4
## Loading required package: ggplot2
## Loading required package: forecast
## Warning: package 'forecast' was built under R version 3.4.4
## Loading required package: fma
## Loading required package: expsmooth
library(fpp2)
setwd("/Users/brendanobrien/Downloads")
BOS.data=read.csv("/Users/brendanobrien/Downloads/monthly-boston-armed-robberies-j.csv", header = T)
BOS.data$Month = NULL
BOS.ts<- ts(data = BOS.data, start = c(1966,1), end = c(1975,10), frequency = 12, ts.eps = getOption("ts.eps"), class ="ts")
as.ts(BOS.ts)
## Jan Feb Mar Apr May Jun Jul Aug Sep Oct Nov Dec
## 1966 41 39 50 40 43 38 44 35 39 35 29 49
## 1967 50 59 63 32 39 47 53 60 57 52 70 90
## 1968 74 62 55 84 94 70 108 139 120 97 126 149
## 1969 158 124 140 109 114 77 120 133 110 92 97 78
## 1970 99 107 112 90 98 125 155 190 236 189 174 178
## 1971 136 161 171 149 184 155 276 224 213 279 268 287
## 1972 238 213 257 293 212 246 353 339 308 247 257 322
## 1973 298 273 312 249 286 279 309 401 309 328 353 354
## 1974 327 324 285 243 241 287 355 460 364 487 452 391
## 1975 500 451 375 372 302 316 398 394 431 431
str(BOS.ts)
## Time-Series [1:118] from 1966 to 1976: 41 39 50 40 43 38 44 35 39 35 ...
#variability increases over time
plot(BOS.ts)
seasonplot(BOS.ts, col = seq(1:109), main = "Season Plot", year.labels = TRUE)
BOS.ets <- ets(BOS.ts)
summary(BOS.ets)
## ETS(M,A,M)
##
## Call:
## ets(y = BOS.ts)
##
## Smoothing parameters:
## alpha = 0.5602
## beta = 1e-04
## gamma = 1e-04
##
## Initial states:
## l = 35.2503
## b = 2.6283
## s = 1.1161 1.0213 1.0058 1.0534 1.1793 1.0679
## 0.8226 0.8787 0.8578 0.988 0.9697 1.0394
##
## sigma: 0.178
##
## AIC AICc BIC
## 1358.263 1364.383 1405.365
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE
## Training set 0.8429871 32.76421 23.26007 -2.485931 13.58244 0.4403586
## ACF1
## Training set -0.1596109
plot(forecast((BOS.ets), h = 14), xlab = "Year", ylab = "A-R")
BOS.ets2 <- ets(BOS.ts, model = "MAN")
summary(BOS.ets2)
## ETS(M,A,N)
##
## Call:
## ets(y = BOS.ts, model = "MAN")
##
## Smoothing parameters:
## alpha = 0.6354
## beta = 1e-04
##
## Initial states:
## l = 43.7709
## b = 2.7402
##
## sigma: 0.199
##
## AIC AICc BIC
## 1374.755 1375.291 1388.608
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE
## Training set 0.7948842 38.82547 28.32731 -3.450859 16.47059 0.5362912
## ACF1
## Training set 0.03327989
plot(forecast((BOS.ets2), h = 14), xlab = "Year", ylab = "A-R")
BOS.ets3 <- ets(BOS.ts, model = "MMM")
summary(BOS.ets3)
## ETS(M,M,M)
##
## Call:
## ets(y = BOS.ts, model = "MMM")
##
## Smoothing parameters:
## alpha = 0.5646
## beta = 2e-04
## gamma = 1e-04
##
## Initial states:
## l = 38.4311
## b = 1.0229
## s = 1.1305 1.0111 0.9992 1.0582 1.1723 1.0953
## 0.8278 0.8736 0.8626 0.9785 0.9622 1.0287
##
## sigma: 0.18
##
## AIC AICc BIC
## 1360.918 1367.038 1408.020
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE
## Training set -2.352028 33.82083 23.67282 -2.476002 13.37398 0.4481726
## ACF1
## Training set -0.1453586
plot(forecast((BOS.ets3), h = 14), xlab = "Year", ylab = "A-R")
#according to AIC and RMSE metrics, the first ETS model using "MAM" was the best moedel
BOS.forecast <- forecast((BOS.ets), h = 14)
forecasted.BOS <- as.numeric(BOS.forecast$mean)
forecasted.BOS
## [1] 422.6540 464.8401 435.6426 408.9856 419.3271 366.3357 377.5530
## [8] 355.6355 464.5082 516.0244 463.7390 445.4265 455.0097 500.1998