You can not select more than 25 topics Topics must start with a letter or number, can include dashes ('-') and can be up to 35 characters long.
ML-For-Beginners/translations/nl/7-TimeSeries/2-ARIMA
leestott 98900475bd
🌐 Update translations via Co-op Translator
3 weeks ago
..
solution 🌐 Update translations via Co-op Translator 3 weeks ago
working 🌐 Update translations via Co-op Translator 3 weeks ago
README.md 🌐 Update translations via Co-op Translator 3 weeks ago
assignment.md 🌐 Update translations via Co-op Translator 3 weeks ago

README.md

Tijdreeksvoorspelling met ARIMA

In de vorige les heb je wat geleerd over tijdreeksvoorspelling en een dataset geladen die de fluctuaties van het elektriciteitsverbruik over een bepaalde periode laat zien.

Introductie tot ARIMA

🎥 Klik op de afbeelding hierboven voor een video: Een korte introductie tot ARIMA-modellen. Het voorbeeld is gedaan in R, maar de concepten zijn universeel.

Quiz voorafgaand aan de les

Introductie

In deze les ontdek je een specifieke manier om modellen te bouwen met ARIMA: AutoRegressive Integrated Moving Average. ARIMA-modellen zijn bijzonder geschikt om data te modelleren die niet-stationair is.

Algemene concepten

Om met ARIMA te kunnen werken, zijn er enkele concepten die je moet begrijpen:

  • 🎓 Stationariteit. Vanuit een statistisch perspectief verwijst stationariteit naar data waarvan de verdeling niet verandert wanneer deze in de tijd wordt verschoven. Niet-stationaire data vertoont daarentegen fluctuaties door trends die moeten worden getransformeerd om te kunnen worden geanalyseerd. Seizoensgebondenheid kan bijvoorbeeld fluctuaties in data veroorzaken en kan worden geëlimineerd door een proces van 'seizoensverschillen'.

  • 🎓 Verschillen nemen. Het nemen van verschillen in data, wederom vanuit een statistisch perspectief, verwijst naar het proces van het transformeren van niet-stationaire data om deze stationair te maken door de niet-constante trend te verwijderen. "Het nemen van verschillen verwijdert de veranderingen in het niveau van een tijdreeks, elimineert trends en seizoensgebondenheid en stabiliseert daardoor het gemiddelde van de tijdreeks." Paper van Shixiong et al

ARIMA in de context van tijdreeksen

Laten we de onderdelen van ARIMA uitpakken om beter te begrijpen hoe het ons helpt tijdreeksen te modelleren en voorspellingen te maken.

  • AR - voor AutoRegressief. Autoregressieve modellen, zoals de naam al aangeeft, kijken 'terug' in de tijd om eerdere waarden in je data te analyseren en aannames over deze waarden te maken. Deze eerdere waarden worden 'lags' genoemd. Een voorbeeld zou data zijn die de maandelijkse verkoop van potloden toont. De verkoopcijfers van elke maand worden beschouwd als een 'evoluerende variabele' in de dataset. Dit model wordt gebouwd door de "evoluerende variabele van interesse te regresseren op zijn eigen vertraagde (d.w.z. eerdere) waarden." wikipedia

  • I - voor Geïntegreerd. In tegenstelling tot de vergelijkbare 'ARMA'-modellen verwijst de 'I' in ARIMA naar het geïntegreerde aspect. De data wordt 'geïntegreerd' wanneer stappen van verschillen nemen worden toegepast om niet-stationariteit te elimineren.

  • MA - voor Moving Average. Het moving-average-aspect van dit model verwijst naar de outputvariabele die wordt bepaald door de huidige en eerdere waarden van lags te observeren.

Kort gezegd: ARIMA wordt gebruikt om een model zo nauwkeurig mogelijk te laten aansluiten bij de speciale vorm van tijdreeksdata.

Oefening - bouw een ARIMA-model

Open de /working-map in deze les en vind het bestand notebook.ipynb.

  1. Voer de notebook uit om de statsmodels Python-bibliotheek te laden; je hebt deze nodig voor ARIMA-modellen.

  2. Laad de benodigde bibliotheken.

  3. Laad nu nog enkele andere bibliotheken die nuttig zijn voor het plotten van data:

    import os
    import warnings
    import matplotlib.pyplot as plt
    import numpy as np
    import pandas as pd
    import datetime as dt
    import math
    
    from pandas.plotting import autocorrelation_plot
    from statsmodels.tsa.statespace.sarimax import SARIMAX
    from sklearn.preprocessing import MinMaxScaler
    from common.utils import load_data, mape
    from IPython.display import Image
    
    %matplotlib inline
    pd.options.display.float_format = '{:,.2f}'.format
    np.set_printoptions(precision=2)
    warnings.filterwarnings("ignore") # specify to ignore warning messages
    
  4. Laad de data uit het bestand /data/energy.csv in een Pandas dataframe en bekijk het:

    energy = load_data('./data')[['load']]
    energy.head(10)
    
  5. Plot alle beschikbare energiedata van januari 2012 tot december 2014. Er zouden geen verrassingen moeten zijn, aangezien we deze data in de vorige les hebben gezien:

    energy.plot(y='load', subplots=True, figsize=(15, 8), fontsize=12)
    plt.xlabel('timestamp', fontsize=12)
    plt.ylabel('load', fontsize=12)
    plt.show()
    

    Nu gaan we een model bouwen!

Maak trainings- en testdatasets

Nu je data is geladen, kun je deze scheiden in trainings- en testsets. Je traint je model op de trainingsset. Zoals gebruikelijk evalueer je na het trainen de nauwkeurigheid van het model met behulp van de testset. Je moet ervoor zorgen dat de testset een latere periode in de tijd beslaat dan de trainingsset, zodat het model geen informatie uit toekomstige tijdsperioden verkrijgt.

  1. Wijs een periode van twee maanden toe van 1 september tot 31 oktober 2014 aan de trainingsset. De testset omvat de periode van twee maanden van 1 november tot 31 december 2014:

    train_start_dt = '2014-11-01 00:00:00'
    test_start_dt = '2014-12-30 00:00:00'
    

    Aangezien deze data het dagelijkse energieverbruik weerspiegelt, is er een sterk seizoensgebonden patroon, maar het verbruik is het meest vergelijkbaar met het verbruik in meer recente dagen.

  2. Visualiseer de verschillen:

    energy[(energy.index < test_start_dt) & (energy.index >= train_start_dt)][['load']].rename(columns={'load':'train'}) \
        .join(energy[test_start_dt:][['load']].rename(columns={'load':'test'}), how='outer') \
        .plot(y=['train', 'test'], figsize=(15, 8), fontsize=12)
    plt.xlabel('timestamp', fontsize=12)
    plt.ylabel('load', fontsize=12)
    plt.show()
    

    trainings- en testdata

    Daarom zou het gebruik van een relatief klein tijdsvenster voor het trainen van de data voldoende moeten zijn.

    Opmerking: Aangezien de functie die we gebruiken om het ARIMA-model te passen in-sample validatie gebruikt tijdens het passen, zullen we validatiedata weglaten.

Bereid de data voor op training

Nu moet je de data voorbereiden op training door filtering en schaling van je data uit te voeren. Filter je dataset om alleen de benodigde tijdsperioden en kolommen op te nemen, en schaal de data zodat deze wordt geprojecteerd in het interval 0,1.

  1. Filter de originele dataset om alleen de eerder genoemde tijdsperioden per set en alleen de benodigde kolom 'load' plus de datum op te nemen:

    train = energy.copy()[(energy.index >= train_start_dt) & (energy.index < test_start_dt)][['load']]
    test = energy.copy()[energy.index >= test_start_dt][['load']]
    
    print('Training data shape: ', train.shape)
    print('Test data shape: ', test.shape)
    

    Je kunt de vorm van de data zien:

    Training data shape:  (1416, 1)
    Test data shape:  (48, 1)
    
  2. Schaal de data naar het bereik (0, 1).

    scaler = MinMaxScaler()
    train['load'] = scaler.fit_transform(train)
    train.head(10)
    
  3. Visualiseer de originele versus geschaalde data:

    energy[(energy.index >= train_start_dt) & (energy.index < test_start_dt)][['load']].rename(columns={'load':'original load'}).plot.hist(bins=100, fontsize=12)
    train.rename(columns={'load':'scaled load'}).plot.hist(bins=100, fontsize=12)
    plt.show()
    

    origineel

    De originele data

    geschaald

    De geschaalde data

  4. Nu je de geschaalde data hebt gekalibreerd, kun je de testdata schalen:

    test['load'] = scaler.transform(test)
    test.head()
    

Implementeer ARIMA

Het is tijd om ARIMA te implementeren! Je gaat nu de statsmodels-bibliotheek gebruiken die je eerder hebt geïnstalleerd.

Nu moet je een aantal stappen volgen:

  1. Definieer het model door SARIMAX() aan te roepen en de modelparameters door te geven: p-, d- en q-parameters, en P-, D- en Q-parameters.
  2. Bereid het model voor op de trainingsdata door de functie fit() aan te roepen.
  3. Maak voorspellingen door de functie forecast() aan te roepen en het aantal stappen (de horizon) te specificeren om te voorspellen.

🎓 Waar zijn al deze parameters voor? In een ARIMA-model zijn er 3 parameters die worden gebruikt om de belangrijkste aspecten van een tijdreeks te modelleren: seizoensgebondenheid, trend en ruis. Deze parameters zijn:

p: de parameter die verband houdt met het autoregressieve aspect van het model, dat verleden waarden incorporeert. d: de parameter die verband houdt met het geïntegreerde deel van het model, dat de hoeveelheid verschillen nemen (🎓 herinner je verschillen nemen 👆?) beïnvloedt die op een tijdreeks wordt toegepast. q: de parameter die verband houdt met het moving-average deel van het model.

Opmerking: Als je data een seizoensgebonden aspect heeft - wat bij deze data het geval is - , gebruiken we een seizoensgebonden ARIMA-model (SARIMA). In dat geval moet je een andere set parameters gebruiken: P, D en Q, die dezelfde associaties beschrijven als p, d en q, maar betrekking hebben op de seizoensgebonden componenten van het model.

  1. Begin met het instellen van je gewenste horizonwaarde. Laten we 3 uur proberen:

    # Specify the number of steps to forecast ahead
    HORIZON = 3
    print('Forecasting horizon:', HORIZON, 'hours')
    

    Het selecteren van de beste waarden voor de parameters van een ARIMA-model kan uitdagend zijn, omdat het enigszins subjectief en tijdrovend is. Je kunt overwegen een auto_arima()-functie te gebruiken uit de pyramid-bibliotheek.

  2. Probeer voorlopig enkele handmatige selecties om een goed model te vinden.

    order = (4, 1, 0)
    seasonal_order = (1, 1, 0, 24)
    
    model = SARIMAX(endog=train, order=order, seasonal_order=seasonal_order)
    results = model.fit()
    
    print(results.summary())
    

    Een tabel met resultaten wordt afgedrukt.

Je hebt je eerste model gebouwd! Nu moeten we een manier vinden om het te evalueren.

Evalueer je model

Om je model te evalueren, kun je de zogenaamde walk forward-validatie uitvoeren. In de praktijk worden tijdreeksmodellen elke keer opnieuw getraind wanneer nieuwe data beschikbaar komt. Dit stelt het model in staat om de beste voorspelling te maken op elk tijdstip.

Begin aan het begin van de tijdreeks met deze techniek, train het model op de trainingsdataset. Maak vervolgens een voorspelling voor de volgende tijdstap. De voorspelling wordt geëvalueerd aan de hand van de bekende waarde. De trainingsset wordt vervolgens uitgebreid met de bekende waarde en het proces wordt herhaald.

Opmerking: Je moet het venster van de trainingsset vast houden voor efficiënter trainen, zodat elke keer dat je een nieuwe observatie aan de trainingsset toevoegt, je de observatie uit het begin van de set verwijdert.

Dit proces biedt een robuustere schatting van hoe het model in de praktijk zal presteren. Het komt echter met de rekenkundige kosten van het maken van zoveel modellen. Dit is acceptabel als de data klein is of als het model eenvoudig is, maar kan een probleem zijn op grotere schaal.

Walk-forward validatie is de gouden standaard voor tijdreeksmodel-evaluatie en wordt aanbevolen voor je eigen projecten.

  1. Maak eerst een testdatapunt voor elke HORIZON-stap.

    test_shifted = test.copy()
    
    for t in range(1, HORIZON+1):
        test_shifted['load+'+str(t)] = test_shifted['load'].shift(-t, freq='H')
    
    test_shifted = test_shifted.dropna(how='any')
    test_shifted.head(5)
    
    load load+1 load+2
    2014-12-30 00:00:00 0.33 0.29 0.27
    2014-12-30 01:00:00 0.29 0.27 0.27
    2014-12-30 02:00:00 0.27 0.27 0.30
    2014-12-30 03:00:00 0.27 0.30 0.41
    2014-12-30 04:00:00 0.30 0.41 0.57

    De data wordt horizontaal verschoven volgens het horizonpunt.

  2. Maak voorspellingen op je testdata met behulp van deze schuifvensterbenadering in een lus ter grootte van de testdatalengte:

    %%time
    training_window = 720 # dedicate 30 days (720 hours) for training
    
    train_ts = train['load']
    test_ts = test_shifted
    
    history = [x for x in train_ts]
    history = history[(-training_window):]
    
    predictions = list()
    
    order = (2, 1, 0)
    seasonal_order = (1, 1, 0, 24)
    
    for t in range(test_ts.shape[0]):
        model = SARIMAX(endog=history, order=order, seasonal_order=seasonal_order)
        model_fit = model.fit()
        yhat = model_fit.forecast(steps = HORIZON)
        predictions.append(yhat)
        obs = list(test_ts.iloc[t])
        # move the training window
        history.append(obs[0])
        history.pop(0)
        print(test_ts.index[t])
        print(t+1, ': predicted =', yhat, 'expected =', obs)
    

    Je kunt het trainen zien gebeuren:

    2014-12-30 00:00:00
    1 : predicted = [0.32 0.29 0.28] expected = [0.32945389435989236, 0.2900626678603402, 0.2739480752014323]
    
    2014-12-30 01:00:00
    2 : predicted = [0.3  0.29 0.3 ] expected = [0.2900626678603402, 0.2739480752014323, 0.26812891674127126]
    
    2014-12-30 02:00:00
    3 : predicted = [0.27 0.28 0.32] expected = [0.2739480752014323, 0.26812891674127126, 0.3025962399283795]
    
  3. Vergelijk de voorspellingen met de daadwerkelijke load:

    eval_df = pd.DataFrame(predictions, columns=['t+'+str(t) for t in range(1, HORIZON+1)])
    eval_df['timestamp'] = test.index[0:len(test.index)-HORIZON+1]
    eval_df = pd.melt(eval_df, id_vars='timestamp', value_name='prediction', var_name='h')
    eval_df['actual'] = np.array(np.transpose(test_ts)).ravel()
    eval_df[['prediction', 'actual']] = scaler.inverse_transform(eval_df[['prediction', 'actual']])
    eval_df.head()
    

    Output

    timestamp h prediction actual
    0 2014-12-30 00:00:00 t+1 3,008.74 3,023.00
    1 2014-12-30 01:00:00 t+1 2,955.53 2,935.00
    2 2014-12-30 02:00:00 t+1 2,900.17 2,899.00
    3 2014-12-30 03:00:00 t+1 2,917.69 2,886.00
    4 2014-12-30 04:00:00 t+1 2,946.99 2,963.00

    Observeer de voorspelling van de uurlijkse data, vergeleken met de daadwerkelijke load. Hoe nauwkeurig is dit?

Controleer de nauwkeurigheid van je model

Controleer de nauwkeurigheid van je model door de gemiddelde absolute procentuele fout (MAPE) over alle voorspellingen te testen.

🧮 Laat me de wiskunde zien

MAPE

MAPE wordt gebruikt om de nauwkeurigheid van voorspellingen weer te geven als een verhouding, gedefinieerd door de bovenstaande formule. Het verschil tussen de werkelijke waarde en de voorspelde waarde wordt gedeeld door de werkelijke waarde.

"De absolute waarde in deze berekening wordt opgeteld voor elk voorspeld tijdstip en gedeeld door het aantal passende punten n." wikipedia

  1. Druk de vergelijking uit in code:

    if(HORIZON > 1):
        eval_df['APE'] = (eval_df['prediction'] - eval_df['actual']).abs() / eval_df['actual']
        print(eval_df.groupby('h')['APE'].mean())
    
  2. Bereken de MAPE van één stap:

    print('One step forecast MAPE: ', (mape(eval_df[eval_df['h'] == 't+1']['prediction'], eval_df[eval_df['h'] == 't+1']['actual']))*100, '%')
    

    MAPE van één stap voorspelling: 0.5570581332313952 %

  3. Print de MAPE van de meerstapsvoorspelling:

    print('Multi-step forecast MAPE: ', mape(eval_df['prediction'], eval_df['actual'])*100, '%')
    
    Multi-step forecast MAPE:  1.1460048657704118 %
    

    Een mooi laag getal is het beste: bedenk dat een voorspelling met een MAPE van 10 een afwijking van 10% heeft.

  4. Maar zoals altijd is het makkelijker om dit soort nauwkeurigheidsmetingen visueel te zien, dus laten we het plotten:

     if(HORIZON == 1):
        ## Plotting single step forecast
        eval_df.plot(x='timestamp', y=['actual', 'prediction'], style=['r', 'b'], figsize=(15, 8))
    
    else:
        ## Plotting multi step forecast
        plot_df = eval_df[(eval_df.h=='t+1')][['timestamp', 'actual']]
        for t in range(1, HORIZON+1):
            plot_df['t+'+str(t)] = eval_df[(eval_df.h=='t+'+str(t))]['prediction'].values
    
        fig = plt.figure(figsize=(15, 8))
        ax = plt.plot(plot_df['timestamp'], plot_df['actual'], color='red', linewidth=4.0)
        ax = fig.add_subplot(111)
        for t in range(1, HORIZON+1):
            x = plot_df['timestamp'][(t-1):]
            y = plot_df['t+'+str(t)][0:len(x)]
            ax.plot(x, y, color='blue', linewidth=4*math.pow(.9,t), alpha=math.pow(0.8,t))
    
        ax.legend(loc='best')
    
    plt.xlabel('timestamp', fontsize=12)
    plt.ylabel('load', fontsize=12)
    plt.show()
    

    een tijdreeksmodel

🏆 Een zeer mooie grafiek, die een model met goede nauwkeurigheid laat zien. Goed gedaan!


🚀Uitdaging

Verdiep je in de manieren om de nauwkeurigheid van een tijdreeksmodel te testen. We behandelen MAPE in deze les, maar zijn er andere methoden die je kunt gebruiken? Onderzoek ze en annoteer ze. Een nuttig document is te vinden hier

Quiz na de les

Review & Zelfstudie

Deze les behandelt alleen de basis van tijdreeksvoorspelling met ARIMA. Neem de tijd om je kennis te verdiepen door deze repository en de verschillende modeltypes te verkennen om andere manieren te leren om tijdreeksmodellen te bouwen.

Opdracht

Een nieuw ARIMA-model


Disclaimer:
Dit document is vertaald met behulp van de AI-vertalingsservice Co-op Translator. Hoewel we streven naar nauwkeurigheid, dient u zich ervan bewust te zijn dat geautomatiseerde vertalingen fouten of onnauwkeurigheden kunnen bevatten. Het originele document in zijn oorspronkelijke taal moet worden beschouwd als de gezaghebbende bron. Voor cruciale informatie wordt professionele menselijke vertaling aanbevolen. Wij zijn niet aansprakelijk voor eventuele misverstanden of verkeerde interpretaties die voortvloeien uit het gebruik van deze vertaling.