.png "DIKW Intelligence"){kind=logo}
Forecast in R
Hoewel er steeds meer informatie op internet te vinden is over Forecasten, is dit nog weinig te vinden in het Nederlands. Ik heb een Nederlandse Tutorial gemaakt zodat ik mijn kennis kan delen.
Voor een grote media uitgeverij heb ik in R een sales forecast / omzetvoorspelling gemaakt voor de periode 2016/2017. In deze (Nederlandse) Forecast tutorial leg ik stap-voor-stap forecasten in R uit met behulp van een voorbeeld dataset en reproduceerbare R codes. Ook wordt uitleg gegeven over het gebruik van (auto)ARIMA.
Forecast package
Forecasten is het voorspellen/ schatten van een toekomstige (financiële) trend of gebeurtenissen. In dit geval gaat het om het voorspellen van toekomstige omzet op basis van beschikbare omzetdata uit het verleden en heden.
In deze blog wordt gebruik gemaakt van de R package forecast. Rob J. Hyndman, maker van deze package, is Professor of Statistics en schrijver van het online textboek Forecasting: principles and practice. Ik heb veel geleerd over forecasting dankzij hem. Ook heb ik het genoegen gehad hem te ontmoeten tijdens een conference in Eindhoven. Meer informatie over professor Hyndman is te vinden op zijn website.
Data
Gebruikte dataset is beschikbaar gesteld op datamarket.com, door Rob J. Hyndman en is te downloaden via deze link. Ga naar ‘Export‘, download het csv(,) bestand en zet de dataset in de juiste working directory.
In deze blog voorspellen we de sales op basis van Monthly car sales in Quebec 1960-1968.
Aangeraden wordt om de dataset te downloaden en iedere stap te volgen, voor het beter begrijpen van de procedure.
Alleen de eerste zeven jaar (96 rijen) van de dataset worden gebruikt, het laatste jaar (rij 97 tot en met 108) wordt bewaard zodat de forecast resultaten vergeleken kunnen worden met de echte waarden.
1 2 3 4 5 6 | #Laad data in R data.orig = read.csv(file="monthly-car-sales-in-quebec-1960.csv") #Gebruik alleen eerste 7 jaar data = data.orig[1:96,] summary(data) |
1 2 3 4 5 6 7 8 9 | ## Month Monthly.car.sales.in.Quebec.1960.1968 ## 1960-01: 1 Min. : 5568 ## 1960-02: 1 1st Qu.:10934 ## 1960-03: 1 Median :13734 ## 1960-04: 1 Mean :14141 ## 1960-05: 1 3rd Qu.:17040 ## 1960-06: 1 Max. :24081 ## (Other):90 |
In het voorbeeld wordt gebruik gemaakt van een geschikte dataset waarvan de distributie al te zien was op de website. Met real-life data wordt aangeraden om eerst de data te controleren. Dat kan door het uitvoeren van exploratieve analyses en het plotten van de data zodat afwijkingen, missende waarden en vreemde patronen in de data ontdekt kunnen worden.
Timeseries
Voordat we de salesdata kunnen gebruiken maken we er eerst een timeseries object van. Controleer en order indien nodig datum in eigen dataset van hoog naar laag. Specificeer het begin en einde en geef een frequentie aan. In dit geval is er 12 keer per jaar een meting gedaan.
Let erop dat de timeserie hetzelfde begin- en eindpunt heeft als de dataframe.
1 2 3 4 | myts = ts(data$Monthly.car.sales.in.Quebec.1960.1968 , start=c(1960, 1), end=c(1967,12), frequency=12) str(myts) |
1 2 | ## Time-Series [1:96] from 1960 to 1968: 6550 8728 12026 14395 14587 13791 9498 8251 ... |
Door het aangeven van de frequentie weet R dat het gaat om een maandelijkse meting. Door het start,- en eindpunt te geven weet R om welke maanden/jaren het gaat.
Voorbeeld tijdreeks
1 2 | plot(myts) |
Onderliggende patronen
Met de functie decompose is het mogelijk om de timeserie op te breken in 3 verschillende timeseries, namelijk de trend , seizoen effecten en als laatste de overgebleven random effecten.
Additive timeserie
In dit geval is er sprake van een ‘additive’ decompositie, een timeserie waarbij seasonality gelijkmatig blijft naarmate de trend toeneemt. Indien de seasonality versterkt wordt bij een sterk toenemende trend, spreken we van een ‘multiplicative’ decomposition.
1 2 3 | myts.dcmp=decompose(myts,type="add") plot(myts.dcmp) |
Arima
Met de juiste achtergrond is het makkelijker om het model te perfectioneren en de uitkomsten beter te begrijpen. In het voorbeeld wordt gebruik gemaakt van de functie auto.arima.
auto.arima
De R functie auto.arima, ontwikkeld door Hyndman and Khandakar is het automatische algoritme dat verschillende aspecten van modelselectie uitvoert. Dit model kan als basismodel gebruikt worden, hoewel dit niet altijd het beste model is en de meest nauwkeurige voorspellingen zal geven. Het uiteindelijke forecastmodel wordt over het algemeen beter indien er meer aan ‘gesleuteld’ wordt. In deze blog leg ik het basismodel uit.
Waar staat ARIMA voor?
ARiMA staat voor AutoRegressive integrated Moving Average.
AR ‘AutoRegressive'(p) is het procesgedeelte waar met behulp van een gewogen som van waarden geschat wordt wat de toekomstige waarden zullen zijn. Met andere woorden, dit gedeelte van het model geeft aan in welke mate vorige waarden een effect hebben op de huidige waarden.
I ‘Integrated’ (d) Voordat er een arima model wordt gemaakt, is het belangrijk dat de data stationary is. Met stationary wordt bedoeld dat de eigenschappen van de data op een bepaald moment onafhankelijk zijn van de tijd. Met andere woorden: er bevindt zich geen trend en seasonality in de tijd.
Het ‘integrated’ gedeelte van een ARiMA model zorgt ervoor dat de data stationary wordt, door te differentiëren. Bij differentiëren neem je het verschil van achtereenvolgende waarden. Hierdoor worden trends en seasonality verwijderd uit data. Indien de data nog steeds niet stationair zijn, kan dit nog een keer worden gedaan.
MA ‘Moving Average'(q) . Een moving average, in het Nederlands voortschrijdend gemiddelde, heeft over het algemeen als belangrijkste eigenschap een onderliggende trend te herkennen. In dit geval wordt de moving average gebruikt voor het voorspellen van toekomstige waarden. Dit gedeelte neemt de errors/residuals uit het model mee om een forecast te kunnen maken voor de toekomst.
Het arima model zal er als volgt uit komen te zien ARIMA(p,d,q)
Het automatische algoritme auto.arima geeft een goede indicatie van het model dat het beste gekozen kan worden, waardoor we niet meer handmatig alle verschillende modelopties hoeven te proberen. Aangeraden wordt wel om ‘omliggende’ modellen te bekijken voor de beste resultaten. Uiteindelijk wijst een zo laag mogelijke AIC/ BIC op een beter model.
Hoe werkt het auto.arima model?
Indien de data non-stationary zijn, dient er ‘gedifferentieerd’ ‘ te worden totdat de data stationary zijn.
Verschillende testen zijn er aanwezig om te controleren of een tijdreeks stationary is. Voorbeelden zijn de Augmented Dickey–Fuller test (ADF) of de Kwiatkowski-Phillips-Schmidt-Shin (KPSS) test. Het automatische algoritme van een auto.arima model combineert de twee verschillende testen en differentieert indien nodig, zodat de gebruiker zelf niet van te voren de data hoeft te differentiëren.
Na het vaststellen van het aantal keer dat er gedifferentieerd zal worden (meestal niet meer dan 2 keer) worden de waarden voor AR(p) en MA(q) stapsgewijs bepaald door het algoritme. Het algoritme zoekt op deze manier het beste model.
Het is ook mogelijk om zelf een arima model te maken.
Hoe kan ik forecasten met auto.arima?
Installeer en / of laad de forecast package
1 2 3 | #install.packages("forecast") require(forecast) |
Maak een auto.arima model aan en bekijk het model:
1 2 3 | m.arima = auto.arima(myts) m.arima |
1 2 3 4 5 6 7 8 9 10 11 | ## Series: myts ## ARIMA(1,0,1)(0,1,1)[12] with drift ## ## Coefficients: ## ar1 ma1 sma1 drift ## 0.7623 -0.4935 -0.5239 82.3198 ## s.e. 0.1825 0.2398 0.1460 16.2801 ## ## sigma^2 estimated as 2326014: log likelihood=-734.85 ## AIC=1479.71 AICc=1480.48 BIC=1491.86 |
Seizoenseffecten
Zoals hierboven is te zien, bestaat er ook een seasonal arima model. In sommige gevallen is het nodig om ‘seasonal differencing’ toe te passen, door het verschil te nemen van de huidige waarde en die van bijvoorbeeld een jaar geleden. Een Arima seasonal model ziet er als volgt uit
ARIMA(p,d,q)(P,D,Q)m de M staat voor het aantal perioden in een seizoen, in dit geval 12.
Forecast
Bovenstaande model gebruiken we voor het forecasten.
Met h kan de forecastperiode aangegeven worden. In dit voorbeeld gebruiken we 12 maanden.
1 2 3 | m.fc = forecast(m.arima,h=12) plot(m.fc) |
Uitkomst Forecast
De donkerblauwe lijn laat de voorspelling voor de komende maanden zien, met daaromheen een 80- en 95% Confidence Interval (CI, betrouwbaarheidsinterval).
Over het algemeen worden de CI breder naarmate de tijdshorizon breder is. Echter stationary modellen (waar d=0) hebben een gelijkmatige CI.
1 2 3 | ## uitkomst forecast m.fc |
1 2 3 4 5 6 7 8 9 10 11 12 13 14 | ## Point Forecast Lo 80 Hi 80 Lo 95 Hi 95 ## Jan 1968 13093.36 11138.80 15047.92 10104.12 16082.60 ## Feb 1968 13253.30 11229.36 15277.24 10157.96 16348.65 ## Mar 1968 21413.70 19350.52 23476.88 18258.33 24569.06 ## Apr 1968 22284.46 20198.81 24370.11 19094.73 25474.18 ## May 1968 24326.90 22228.30 26425.49 21117.37 27536.42 ## Jun 1968 22553.02 20446.94 24659.10 19332.04 25773.99 ## Jul 1968 16358.20 14247.78 18468.62 13130.59 19585.81 ## Aug 1968 15339.44 13226.50 17452.38 12107.98 18570.90 ## Sep 1968 14451.79 12337.40 16566.19 11218.10 17685.49 ## Oct 1968 18338.88 16223.64 20454.13 15103.89 21573.87 ## Nov 1968 18332.40 16216.67 20448.14 15096.66 21568.15 ## Dec 1968 16010.55 13894.53 18126.57 12774.37 19246.73 |
Plot forecast met echte waarden
Het laatste jaar van de gebruikte data is apart gehouden zodat de voorspelling vergeleken kan worden met echte uitkomsten.
1 2 3 4 5 6 7 8 | # Maak een validatie-timeseries df.vld = data.orig[97:108,] ts.vld = ts(df.vld$Monthly.car.sales.in.Quebec.1960.1968, start=c(1968, 1), end=c(1968, 12), frequency=12) # Plot de echte data in forecastplot m.fc = forecast(m.arima,h=12) plot(m.fc) lines(ts.vld, col="red") |
Model verkennen
Vaak loont het om de omliggende modellen van het ,door auto.arima gegeven, model te verkennen.
Aangeraden wordt om de omliggende parameters p en q (AR en MA) te veranderen door + of – 1 . Hiervoor kan Arima gebruikt worden. Het aantal maal dat er gedifferentieerd wordt blijft hetzelfde. Daarnaast kan een drift toegevoegd worden. Een drift zorgt ervoor dat de veranderingen in de tijd van de voorspelling vergelijkbaar zijn met de veranderingen in het verleden (bijvoorbeeld een opwaartse trend).
1 2 3 | fit = Arima(myts, order=c(2,0,0), seasonal=c(1,1,1),include.drift = TRUE) fit |
1 2 3 4 5 6 7 8 9 10 11 | ## Series: myts ## ARIMA(2,0,0)(1,1,1)[12] with drift ## ## Coefficients: ## ar1 ar2 sar1 sma1 drift ## 0.2249 0.2117 0.3819 -0.9901 81.8131 ## s.e. 0.1072 0.1110 0.1497 0.7087 11.6620 ## ## sigma^2 estimated as 1994176: log likelihood=-734.04 ## AIC=1480.09 AICc=1481.18 BIC=1494.67 |
In sommige gevallen wordt er in auto.arima een model gemaakt met een lagere AIC waarden als ‘stepwise’ uit wordt gezet. Het nadeel hiervan is dat de berekening een stuk trager is.
1 2 | auto.arima(myts,stepwise = FALSE) |
1 2 3 4 5 6 7 8 9 10 11 | ## Series: myts ## ARIMA(2,0,0)(0,1,1)[12] with drift ## ## Coefficients: ## ar1 ar2 sma1 drift ## 0.2662 0.2011 -0.4911 82.5901 ## s.e. 0.1079 0.1116 0.1522 15.1559 ## ## sigma^2 estimated as 2335120: log likelihood=-734.77 ## AIC=1479.53 AICc=1480.3 BIC=1491.69 |
In dit geval loont het om stepwise uit te zetten.
1 2 3 | fit = Arima(myts, order=c(2,0,0), seasonal=c(0,1,1),include.drift = TRUE) fit |
1 2 3 4 5 6 7 8 9 10 11 | ## Series: myts ## ARIMA(2,0,0)(0,1,1)[12] with drift ## ## Coefficients: ## ar1 ar2 sma1 drift ## 0.2662 0.2011 -0.4911 82.5901 ## s.e. 0.1079 0.1116 0.1522 15.1559 ## ## sigma^2 estimated as 2335120: log likelihood=-734.77 ## AIC=1479.53 AICc=1480.3 BIC=1491.69 |
1 2 3 4 5 | # Plot de echte data in forecastplot m.fc = forecast(fit,h=12) plot(m.fc) lines(ts.vld, col="red") |
Conclusie
In deze blog heb ik de basis van arima uitgelegd met een voorbeelddataset. Het basis model is nog niet perfect maar komt al wel in de buurt van een nauwkeurige forecast. In de blogs die volgen ga ik verder in op het handmatig bouwen van een forecastmodel, het toevoegen van extra voorspellers en het effect van de gekozen tijdreeksfrequentie.
Voor vragen en opmerkingen ben ik bereikbaar via email. Ook voor extra hulp bij het forecasten kan informatie worden gevraagd.
Meer weten over wat data voor uw organisatie kan betekenen?
Neem contact op met Hugo Koopmans
Telefoon: +31643106780
E-mail: hugo.koopmans@dikw.com
Blogs
-
Datagedreven werken Deel 3 — door Nick van de Venn — last modified 03-10-2023
- Welke data heeft eigenlijk waarde voor uw organisatie?
-
Differential Privacy — door Nick van de Venn — last modified 18-09-2023
- Gevoelige gegevens verwerken zonder dat de gevoelige informatie kan uitlekken
-
Datagedreven werken Deel 2 — door Nick van de Venn — last modified 03-10-2023
- Welke data heeft eigenlijk waarde voor uw organisatie?
-
ChatGPT for Business Intelligence — door Nick van de Venn — last modified 18-09-2023
- Chatten met je datawarehouse, utopie of werkelijkheid?
-
Intelligence Factory — door Nick van de Venn — last modified 05-07-2023
- Agile design thinking met een ML-ops sausje
-
Data gedreven werken. Deel 1 — door Nick van de Venn — last modified 03-10-2023
- Bij datagedreven werken staan mens én data centraal
-
Hoe gebruik je ChatGPT om je data pipeline te bouwen? — door Nick van de Venn — last modified 11-05-2023
- Marcel-Jan zocht het uit voor zijn grote hobby: astronomie.
-
Lead consultant en manager Business Analytics Patrick Meulstee — door Nick van de Venn — last modified 06-04-2023
- Over remote werken op Bonaire
-
Hoogwaardige data opleidingen bij DIKW Academy — door nick van de venn — last modified 28-02-2023
- Waarom wij geloven in de praktische opzet van onze cursussen
-
Lead data engineer Remy Lamberty — door Nick van de Venn — last modified 27-03-2023
- Over remote werken vanuit India
-
Van pizzakoerier naar business intelligence consultant — door Nick van de Venn — last modified 15-12-2022
- Het verhaal hoe Brian een vaste waarde werd binnen DIKW
-
Bayesiaanse Statistiek — door Marc Jacobs — last modified 25-07-2022
- Wiskundig raamwerk voor ouderwets leren
-
R vs Python — door Nick van de Venn — last modified 05-07-2022
- Samenwerken is de sleutel
-
Je fietsroutes eenvoudig in kaart brengen... — door Nick van de Venn — last modified 05-05-2022
- Marcel-Jan doet het eenvoudig met behulp van Python.
-
Het DIKW model — door marco — last modified 24-02-2022
- In vier stappen waarde creëren met data
-
DIKW Academy: Waar theorie en praktijk samen komen — door marco — last modified 21-02-2022
- DIKW docenten delen hun expertise uit de praktijk
-
Laat data voor u renderen! — door marco — last modified 08-02-2022
- Bluemine: Analytics as a service
-
30 jaar intelligence: Nieuwe uitdagingen om met data waarde toe te voegen — door marco — last modified 08-02-2022
- Van oude computerterminal naar smartphone
-
De fascinerende wereld van testen — door marco — last modified 30-12-2021
- Verwacht het onverwachte
-
Verzekeraar creëert meerwaarde met slimme data hub — door marco — last modified 20-12-2021
- Klant maakt met betere voorspellingen met data
-
Op naar een mooi data gedreven 2022! — door marco — last modified 14-12-2021
- Data gedreven organisatie dient blijvend te worden gevoed
-
Machine Learning: De gereedschapskist van de data scientist — door marco — last modified 21-12-2021
- Machine Learning algoritmes zijn de gereedschappen voor een data scientist
-
Met data de wind in de zeilen — door marco — last modified 28-01-2022
- Met data management kiest u de juiste koers
-
Hoe ethisch is Facebook? — door marco — last modified 05-11-2021
- Is regulering en wetgeving voor AI nodig?
-
Boekbespreking: Data Teams van Jesse Anderson — door marco — last modified 02-11-2021
- Voor succesvolle big data projecten zijn drie teams nodig
-
Smells like AI — door marco — last modified 01-11-2021
- Artificial Intelligence creëert nieuwe muziek
-
De waarde van data voor het MKB — door marco — last modified 09-12-2021
- Bluemine ontzorgt MKB door data beheer
-
Data & AI: Kans of bedreiging? — door marco — last modified 21-12-2021
- Waarde creëren met data en AI zorgt voor nieuwe business mogelijkheden
-
Data gedrevenheid is proces van lange adem — door marco — last modified 02-11-2021
- Data is een ingrediënt dat zorgt voor meerwaarde op lange termijn
Data Science recente blogs
-
- Chatten met je datawarehouse, utopie of werkelijkheid?
-
- Agile design thinking met een ML-ops sausje
-
- Wiskundig raamwerk voor ouderwets leren
Data Science Nieuws & Evenementen
-
Data science opleidingen gaan weer van start! — door marco — last modified 08-02-2022
- Vanaf 21 september start ons succesnummer weer! Twaalf weken data science in R, we hebben er weer zin in
-
Aedes data science workshop 2 van 3 — door marco — last modified 07-02-2022
- Voor Aedes organiseert DIKW drie workshops data science
-
DIKW in top 50 beste data science bedrijven — door marco — last modified 22-10-2021
- DIKW is één van snelst groeiende bedrijven volgens MKB Data Science rapport