1 Introduction

For this report, I will be analayzing Monthly Natural Gas Prices. Specifically I will be analyzing gas prices from January 2008 to August of 2020 To do this, I will build a time series with the data and then build and then compare various smoothing models.

2 Time Series Analysis

2.1 Training and Testing Data

To start, I will split the data into training and testing sets. I will withhold the last twelve observations for testing the various models that i will build.

Gas <- read.csv("MonthlyGasPrices.csv")
n.row = dim(Gas)[1]
data.Gas = Gas[(n.row-151):n.row, ]
test.Gas = data.Gas$Price[141:152]
train.Gas = data.Gas$Price[1:140]
Gas.ts=ts(Gas$Price[1:152], start=2008, frequency = 12)

2.2 Smoothing Models

Next, I will build the various smoothing models and then get their accuracy measures.

fit1 = ses(Gas.ts, h=12)
fit2 = holt(Gas.ts, initial="optimal", h=12)             ## optimal alpha and beta
fit3 = holt(Gas.ts,damped=TRUE, h=12 )                  ## additive damping
fit4 = holt(Gas.ts,exponential=TRUE, damped=TRUE, h =12) ## multiplicative damp
fit5 = hw(Gas.ts,h=12, seasonal="additive")              ## default h = 10
fit6 = hw(Gas.ts,h=12, seasonal="multiplicative")
fit7 = hw(Gas.ts,h=12, seasonal="additive",damped=TRUE)
fit8 = hw(Gas.ts,h=12, seasonal="multiplicative",damped=TRUE)

accuracy.table = round(rbind(accuracy(fit1), accuracy(fit2), accuracy(fit3), accuracy(fit4),
                             accuracy(fit5), accuracy(fit6), accuracy(fit7), accuracy(fit8)),4)
row.names(accuracy.table)=c("SES","Holt Linear","Holt Add. Damped", "Holt Exp. Damped",
                            "HW Add.","HW Exp.","HW Add. Damp", "HW Exp. Damp")
kable(accuracy.table, caption = "The accuracy measures of various exponential smoothing models 
      based on the training data")
The accuracy measures of various exponential smoothing models based on the training data
ME RMSE MAE MPE MAPE MASE ACF1
SES -0.0018 0.9308 0.5938 -1.1318 10.9038 0.3000 -0.0013
Holt Linear 0.0001 0.9343 0.6019 -1.0737 11.1148 0.3041 0.0048
Holt Add. Damped 0.0105 0.9293 0.5951 -0.5711 10.9418 0.3007 -0.0017
Holt Exp. Damped -0.0030 0.9340 0.6013 -1.1707 11.1011 0.3038 0.0041
HW Add. 0.0337 0.8992 0.5973 -0.0425 11.8137 0.3018 0.0026
HW Exp. 0.0584 0.8597 0.6054 0.6680 11.9053 0.3059 0.0604
HW Add. Damp 0.0126 0.8996 0.6069 -0.4287 12.2262 0.3067 0.0249
HW Exp. Damp 0.0006 0.8636 0.6060 -0.7841 11.9056 0.3062 0.0424

According to the table, the Holt Winters Exponential Damped model is the best so far.

Now to compare nonseasonal vs seasonal smoothing models.

par(mfrow=c(2,1), mar=c(3,4,3,1))
###### plot the original data
pred.id = 141:152
plot(1:140, train.Gas, lwd=2,type="o", ylab="Prices (Dollars)", xlab="", 
     xlim=c(1, 152), ylim=c(0, 20), cex=0.3,
     main="Non-seasonal Smoothing Models")
lines(pred.id, fit1$mean, col="red")
lines(pred.id, fit2$mean, col="blue")
lines(pred.id, fit3$mean, col="purple")
lines(pred.id, fit4$mean, col="navy")
##
points(pred.id, fit1$mean, pch=16, col="red", cex = 0.5)
points(pred.id, fit2$mean, pch=17, col="blue", cex = 0.5)
points(pred.id, fit3$mean, pch=19, col="purple", cex = 0.5)
points(pred.id, fit4$mean, pch=21, col="navy", cex = 0.5)
#points(fit0, col="black", pch=1)
legend("bottomright", lty=1, col=c("red","blue","purple", "navy"),pch=c(16,17,19,21),
   c("SES","Holt Linear","Holt Linear Damped", "Holt Multiplicative Damped"), 
   cex = 0.7, bty="n")
###########
plot(1:140, train.Gas, lwd=2,type="o", ylab="Beverage", xlab="", 
     xlim=c(1,200), ylim=c(0, 20), cex=0.3,
     main="Holt-Winterd Teend and Seasonal Smoothing Models")
lines(pred.id, fit5$mean, col="red")
lines(pred.id, fit6$mean, col="blue")
lines(pred.id, fit7$mean, col="purple")
lines(pred.id, fit8$mean, col="navy")
##
points(pred.id, fit5$mean, pch=16, col="red", cex = 0.5)
points(pred.id, fit6$mean, pch=17, col="blue", cex = 0.5)
points(pred.id, fit7$mean, pch=19, col="purple", cex = 0.5)
points(pred.id, fit8$mean, pch=21, col="navy", cex = 0.5)
###
legend("bottomright", lty=1, col=c("red","blue","purple", "navy"),pch=c(16,17,19,21),
   c("HW Additive","HW Multiplicative","HW Additive Damped", "HW Multiplicative Damped"), 
   cex = 0.7, bty="n")

It’s hard to tell from the graph which is the best. But we at least can see that Holt-Winters is doing much better than the non seasonal smoothing models.

acc.fun = function(test.data, mod.obj){
  PE=100*(test.data-mod.obj$mean)/mod.obj$mean
  MAPE = mean(abs(PE))
  ###
  E=test.data-mod.obj$mean
  MSE=mean(E^2)
  ###
  accuracy.metric=c(MSE=MSE, MAPE=MAPE)
  accuracy.metric
}

pred.accuracy = rbind(SES =acc.fun(test.data=test.Gas, mod.obj=fit1),
                      Holt.Add =acc.fun(test.data=test.Gas, mod.obj=fit2),
                      Holt.Add.Damp =acc.fun(test.data=test.Gas, mod.obj=fit3),
                      Holt.Exp =acc.fun(test.data=test.Gas, mod.obj=fit4),
                      HW.Add =acc.fun(test.data=test.Gas, mod.obj=fit5),
                      HW.Exp =acc.fun(test.data=test.Gas, mod.obj=fit6),
                      HW.Add.Damp =acc.fun(test.data=test.Gas, mod.obj=fit7),
                      HW.Exp.Damp =acc.fun(test.data=test.Gas, mod.obj=fit8))
kable(pred.accuracy, caption="The accuracy measures of various exponential smoothing models 
      based on the testing data")
The accuracy measures of various exponential smoothing models based on the testing data
MSE MAPE
SES 1.3032866 34.69953
Holt.Add 1.3848653 35.35283
Holt.Add.Damp 1.3018345 34.68724
Holt.Exp 1.3236138 34.87465
HW.Add 1.6139454 36.46042
HW.Exp 0.6821635 27.59715
HW.Add.Damp 1.8197499 37.95936
HW.Exp.Damp 1.0441864 32.01474

According to this table, the Holt Winters Additive Damped Model is the best. However, earlier the Holt winters exponential smoothing method worked best to minimize the various errors. These two methods actually did very similar the second time around.

2.3 Parameters

We can now calculate our alpha, beta, and gamma values. To do so, we must refit the model with the entire data set.

Gas <- read.csv("MonthlyGasPrices.csv")
n.row = dim(Gas)[1]
data.Gas = Gas[(n.row-151):n.row, ]
Gas.ts = ts(data.Gas$Price[1:152], start=2008, frequency = 12)
final.model = hw(Gas.ts, h=12, seasonal="additive") 
smoothing.parameter = final.model$model$par[1:3]
kable(smoothing.parameter, caption="Estimated values of the smoothing parameters in
      Holt-Winters linear trend with additive seasonality")
Estimated values of the smoothing parameters in Holt-Winters linear trend with additive seasonality
x
alpha 0.9998141
beta 0.0669256
gamma 0.0001844

3 Summary and Conclusion

This case study focused on Natural Gas Prices from January 2008 to August 2020. A time series was built, and various smoothing methods were used. Holt-Winters Additive and Exponential worked very similarly, so it was hard to tell which is better. Finally, our three parameter values were obtained, with alpha = 0.9998, beta = 0.0669, and gamma = 0.0002.