#Markmið Markmið okkar í þessu er að búa til spálíkan til þess að spá um sölu næstu mánaðar, við prufum nokkur mismunandi líkön og metum þau út frá ME MAD RMSE og MAPE.

#uppsetning gagna setjum upp pakka

library(zoo)
## 
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
## 
##     as.Date, as.Date.numeric
# Load the package
library(forecast)
## Registered S3 method overwritten by 'quantmod':
##   method            from
##   as.zoo.data.frame zoo

aðlögum svo gögn

Solugogn <- read.csv(file = '/Users/hilmir/Downloads/gögn_heimadæmi_i.csv', sep = ";")
names(Solugogn) <- c("Timabil", "Sala")

Solugogn$Timabil <- as.Date(Solugogn$Timabil, format="%Y-%m-%d")

Sala_timaseria <- ts(Solugogn$Sala, start=c(2014, 1), frequency=12)

#Moving average

Við byrjum á því að nota moving average aðferðina þar sem meðaltal gilda er tekið á ákvðenutímabili, í okkar tilfelli 12,6 og 3, og notað til þess að spá fyrir næsta gildi.

Það sem við gerum hér er að við búaum til for-lykkju sem rúllar í gegnum alla mánuði og spáir svo fyrir næsta gildi. Geymir það og heldur áfram fyrir alla mánuðina og reyknum svo skekkju, ME,MAD,RMSE,MAPE

##Hér er má sjá moving average fyrir 12 mánuði

spagildi <- rep(NA, length(Sala_timaseria))
raungildi <- rep(NA, length(Sala_timaseria))

Ma_manudir <- 12


for(i in (Ma_manudir + 1): length(Sala_timaseria)) {
  
  spagildi[i] <- mean(Sala_timaseria[(i-Ma_manudir):(i-1)])
  
  
  raungildi[i] <- Sala_timaseria[i]
}

skekkja <- spagildi - raungildi
me_ma12 <- mean(skekkja, na.rm = TRUE)
mad_ma12 <- mean(abs(skekkja), na.rm = TRUE)
rmse_ma12 <- sqrt(mean(skekkja^2, na.rm = TRUE))
mape_ma12 <- mean(abs(skekkja/raungildi), na.rm = TRUE) * 100


cat("ME:", me_ma12, "\nMAD:", mad_ma12, "\nRMSE:", rmse_ma12, "\nMAPE:", mape_ma12, "\n")
## ME: -3849365 
## MAD: 9289034 
## RMSE: 12571466 
## MAPE: 21.74893

ME og MAD er hátt en MAPE(Mean Absolute Percentage Error) er 21% sem er nokkuð gott meðað við önnur líkön hér

# Plot the actual vs forecasted values
plot(raungildi, type="l", col="blue", main="MA 12 spá vs Raungildi", xlab="Mánuður nr.", ylab="Value")
lines(spagildi, type="l", col="green")
legend("topright", legend=c("Raungidi", "Spá"), col=c("blue", "green"), lty=1)

Sjáum að spáinn fylgir vel breytingum yfir árið en er ekki gott í að mæla árstíðabundnar sveiflur eins og mátti gera ráð fyrir vegna þess að hún tekur heilsársgögn

##Hér er má sjá moving average fyrir 6 mánuði

spagildi <- rep(NA, length(Sala_timaseria))
raungildi <- rep(NA, length(Sala_timaseria))

Ma_manudir <- 6
for(i in (Ma_manudir + 1): length(Sala_timaseria)) {
  spagildi[i] <- mean(Sala_timaseria[(i-Ma_manudir):(i-1)])
  raungildi[i] <- Sala_timaseria[i]
}


skekkja <- spagildi - raungildi
me_ma6 <- mean(skekkja, na.rm = TRUE)
mad_ma6 <- mean(abs(skekkja), na.rm = TRUE)
rmse_ma6 <- sqrt(mean(skekkja^2, na.rm = TRUE))
mape_ma6 <- mean(abs(skekkja/raungildi), na.rm = TRUE) * 100


cat("ME:", me_ma6, "\nMAD:", mad_ma6, "\nRMSE:", rmse_ma6, "\nMAPE:", mape_ma6, "\n")
## ME: -1809023 
## MAD: 11694844 
## RMSE: 14524792 
## MAPE: 33.0032

Hér er ME er minna sem er áhugavert vegna þess að bæði MAD og MAPE er töluvert hærra en í MA 12, og úr því má draga að módelið er jafnviltlaust alltaf. Sem kannski passar vegna þess að summar trend í Júlí fylgja ennþá í desember.

plot(raungildi, type="l", col="blue", main="MA 6 spá vs Raungildi", xlab="Mánuður nr.", ylab="Value")
lines(spagildi, type="l", col="green")
legend("topright", legend=c("Raungidi", "Spá"), col=c("blue", "green"), lty=1)

Myndin staðfestir kenningu okkar að þegar spáin er háð þá er gildið hátt þegar salan er lág og gildið er lágt þegar salan er há

##Hér má sjá moving average fyrir 3 mánuði

spagildi <- rep(NA, length(Sala_timaseria))
raungildi <- rep(NA, length(Sala_timaseria))


Ma_manudir <- 3


for(i in (Ma_manudir + 1): length(Sala_timaseria)) {
  
  spagildi[i] <- mean(Sala_timaseria[(i-Ma_manudir):(i-1)])
  
  
  raungildi[i] <- Sala_timaseria[i]
}


skekkja <- spagildi - raungildi
me_ma3 <- mean(skekkja, na.rm = TRUE)
mad_ma3 <- mean(abs(skekkja), na.rm = TRUE)
rmse_ma3 <- sqrt(mean(skekkja^2, na.rm = TRUE))
mape_ma3 <- mean(abs(skekkja/raungildi), na.rm = TRUE) * 100


cat("ME:", me_ma3, "\nMAD:", mad_ma3, "\nRMSE:", rmse_ma3, "\nMAPE:", mape_ma3, "\n")
## ME: -820884.5 
## MAD: 10502794 
## RMSE: 12647819 
## MAPE: 29.17368

Hér er má sjá að gildi ME hefur lækkað enn meira en MAD er ennþá hærra heldur fyrir 6 mánuði og MAPE-ið er líka töluvert hærra

plot(raungildi, type="l", col="blue", main="MA 3 vs Raungildi", xlab="Mánuður nr.", ylab="Value")
lines(spagildi, type="l", col="green")
legend("topright", legend=c("Raungidi", "Spá"), col=c("blue", "green"), lty=1)

Hér virðist vera smá seinkunn í spám er það vegna þess að þrátt fyrir það tekur smá tíma fyrir stórar sveiflur að ná því að hafa áhrif.

MA 12 er ótrúlegt en satt besta módelið til þess að spá fyrir um sölunæsta mánaðar. En þó er MA3 módelið sem er best að spá fyrir heildarþörf en það er bara seinkunn í sveiflum.

#Exponenstial smoothing Exponenstial smoothing er aðferð sem setur meiri áherslu á nýjusutugildi en tekur þó önnur gildi inn og spáir þannig fyrir næsta mánuði.

Hér notum við sömu aðferð og í MA að við rúllum í gegnum módelið og spáum fyrir næsta mánuði hverju sinni með uppfærðum gögnum. Forritið finnur þau gildi sem henta best.

spagildi <- rep(NA, length(Sala_timaseria))
raungildi <- rep(NA, length(Sala_timaseria))


for(i in 2:(length(Sala_timaseria)-1)) {

  model <- ets(Sala_timaseria[1:i])

  next_forecast <- forecast(model, h=1)
  spagildi[i+1] <- next_forecast$mean

  raungildi[i+1] <- Sala_timaseria[i+1]
}


skekkja <- spagildi - raungildi
me_ets <- mean(skekkja, na.rm = TRUE)
mad_ets <- mean(abs(skekkja), na.rm = TRUE)
rmse_ets <- sqrt(mean(skekkja^2, na.rm = TRUE))
mape_ets <- mean(abs(skekkja/raungildi), na.rm = TRUE) * 100


cat("ME:", me_ets, "\nMAD:", mad_ets, "\nRMSE:", rmse_ets, "\nMAPE:", mape_ets, "\n")
## ME: 563100.7 
## MAD: 6775657 
## RMSE: 8712341 
## MAPE: 19.63658

Hér lýst okkur strax betur á niðurstöðurnar ME er minna en það hefur nokkru sinni verið en spáinn spáir aðeins of miklli sölu yfir allt tímabilið. MAD-ið og MAPE-ið er einnig mun betra en MAPE-ið um 19.6% og er það því betra en öll hin að finna spá næsta mánaðar.

plot(raungildi, type="l", col="blue", main="Exponential Smoothing spá vs Raungildi", xlab="Mánuður nr.", ylab="Value")
lines(spagildi, type="l", col="green")
legend("topright", legend=c("Raungidi", "Spá"), col=c("blue", "green"), lty=1)

Sjáum hér að sveiflurnar fylgja mjög með seinkunn um einn mánuð eins og mátti gera ráð fyrir enda er áhersla lögð mest á síðasta mánuð.

#ets2

Prófum aðra aðferð sem við notum summary skipiun sem gefur okkur þau gildi sem við þurfum ég veit þó ekki hvernig þessi gildi eru fundin.

model <- ets(Sala_timaseria)
summary(model)
## ETS(A,N,A) 
## 
## Call:
##  ets(y = Sala_timaseria) 
## 
##   Smoothing parameters:
##     alpha = 0.9999 
##     gamma = 1e-04 
## 
##   Initial states:
##     l = 25564904.1281 
##     s = -10962425 -8245448 -2293852 3482559 21417551 19749770
##            9864210 -5503431 -8605103 -3327045 -6400110 -9176676
## 
##   sigma:  2513849
## 
##      AIC     AICc      BIC 
## 2028.198 2039.107 2059.613 
## 
## Training set error measures:
##                    ME    RMSE     MAE       MPE     MAPE      MASE       ACF1
## Training set 328447.5 2201114 1654355 0.8751075 5.902972 0.2337661 0.05911586

Þegar ég set exponenstial smoothing upp í summary fallinu þá eru öll gildi miklubetri en áður MAPE komið niður í 6 og RMSE og ME bæði minni en nokkru sinni fyrr. En eins og ég er þó ekki alveg klár á því hvernig þetta finnur gildin og er þetta nokkurs konar “Black box” svo við skulum hafa það í huga.

#Arima Hér notum við sambærilega aðferð og við höfum verið að gera þar sem keyrum Arima fyrir hvern mánuð og spáum fyrir þeim næsta og líkt og í exponential smoothing finnur ARIMA bestu niðurstöður frá sögulegum gögnum. ##Arima1

spagildi <- rep(NA, length(Sala_timaseria))
raungildi <- rep(NA, length(Sala_timaseria))


for(i in 2:(length(Sala_timaseria)-1)) {

  model_arima <- auto.arima(Sala_timaseria[1:i])
  

  next_forecast <- forecast(model_arima, h=1)
  spagildi[i+1] <- next_forecast$mean
  

  raungildi[i+1] <- Sala_timaseria[i+1]
}

skekkja <- spagildi - raungildi
me_arima <- mean(skekkja, na.rm = TRUE)
mad_arima <- mean(abs(skekkja), na.rm = TRUE)
rmse_arima <- sqrt(mean(skekkja^2, na.rm = TRUE))
mape_arima <- mean(abs(skekkja/raungildi), na.rm = TRUE) * 100


cat("ME:", me_arima, "\nMAD:", mad_arima, "\nRMSE:", rmse_arima, "\nMAPE:", mape_arima, "\n")
## ME: -1187649 
## MAD: 6374186 
## RMSE: 8350935 
## MAPE: 18.93714

Sjáum hér að MAPE og MAD hefur aldrei verið betra(þ.e.a.s. með rúlla í gegnum mánuði aðferð ) eða 18.9% sem þýðir að þau eru 18.9% að eðaltali frá raungildi en samtals heildarskekkja er þó meiri en fyrir exponentsial smoothing.

plot(raungildi, type="l", col="blue", main="ARIMA vs Raungildi", xlab="Mánuður nr.", ylab="Value")
lines(spagildi, type="l", col="red")
legend("topright", legend=c("Raungidi", "Spá"), col=c("blue", "red"), lty=1)

Hér má sjá hvernig exponential smoothing fylgir mjög vel raungögnum og virðist verð betra því meira sem líður á notkun módelsins og er síðasta gildið nánast það sama.

##arima2

Notum aftur summary aðferðina.

model <- auto.arima(Sala_timaseria)
summary(model)
## Series: Sala_timaseria 
## ARIMA(0,1,1)(0,1,0)[12] 
## 
## Coefficients:
##           ma1
##       -0.2930
## s.e.   0.1521
## 
## sigma^2 = 6.382e+12:  log likelihood = -759.11
## AIC=1522.23   AICc=1522.5   BIC=1525.93
## 
## Training set error measures:
##                     ME    RMSE     MAE        MPE     MAPE      MASE       ACF1
## Training set -71104.46 2211958 1459975 -0.5563043 4.158545 0.2062994 0.02780417

Þegar við notum summary fallið sést aftur að gildin eru rosalega góð MAPE komið niður í 4.16% og allt virðist vera mjög gott. En hafa má aftur í huga að ég veit ekki almennilega hvernig þetta virkar eða hvort það gæti verið villa í mínum útreikningum.

En ARIMA virðist vera best sama hvort aðferðina sem notum. Því notum við hana til þess að spá fyrir um næstu 12 mánuði.

spa <- forecast(model, h=12)

# Plot the forecast stored in spa along with Sala_timaseria
plot(spa, xlab = "Ár", ylab = "Values", main = "Spá + sala")
lines(Sala_timaseria, col = "blue")  # Add the Sala_timaseria data as a blue line
legend("topleft", legend = c("Forecast", "Sala_timaseria"), col = c("black", "blue"), lty = 1, cex = 0.8)

Sjáum hér hvernig spáin fylgir gömlum sveiflum horfum og færum okkur yfir í excel

point_forecasts <- spa$mean
point_forecasts
##           Jan      Feb      Mar      Apr      May      Jun      Jul      Aug
## 2019 37930650 41057900 44303650 37926150 42348400 59506900 70691650 73874400
##           Sep      Oct      Nov      Dec
## 2019 58958650 50944900 38552900 35345900