Python: analyse des séries temporelles: création d'un modèle SARIMA

Construction du modèle SARIMA (1)

Examen du modèle ARIMA

Le modèle ARIMA est un modèle influencé par les valeurs précédentes Avec le modèle AR AR (p) qui est en corrélation avec les dernières valeurs p Un modèle MA MA (q) qui est affecté par les valeurs q précédentes dans un modèle affecté par l'erreur précédente a été synthétisé. ARMA (p, q) a été adapté à la série de différences avant le point d.

Et Quel est le modèle SARIMA? Le modèle ARIMA peut également être utilisé pour les données de séries chronologiques avec un cycle saisonnier C'est un modèle qui peut être étendu.

Le modèle SARIMA a des paramètres (p, d, q) en plus de Il a également les paramètres (sp, sd, sq, s).

Que sont sp, sd, sq?

sp :Autocorrélation saisonnière
sd :Dérivation de la saisonnalité
sq :Moyenne mobile saisonnière
C'est appelé.

Au fait, ARIMA(p,d,q)Du modèle

・ P:C'est ce qu'on appelle le degré d'autocorrélation et si le modèle est prédit en utilisant les p dernières valeurs.

· RÉ:C'est ce qu'on appelle l'induction, et le fait que la différence d'ordre d-ième était nécessaire pour stabiliser les données de la série chronologique.

・ Q:On l'appelle moyenne mobile, et le fait que le modèle soit affecté par les dernières valeurs q

Cela représentait.

Les significations de base de sp, sd et sq ne changent pas. pourtant Pour sp, sd, sq, les données actuelles sont influencées par les données historiques qui ont dépassé une ou plusieurs périodes saisonnières.

Par exemple, dans le cas de données avec fluctuations saisonnières d'un cycle de 12 mois, le paramètre de s représente le cycle. Cela devient s = 12.

Si sq = 1, données pour exactement 12 mois (il y a 1 cycle) Si sq = 2, cela signifie que le modèle est affecté par les données d'il y a 12 mois et les données d'il y a 24 mois.

Si c'est difficile à comprendre, remplacez simplement q par sq.

Déterminer les paramètres

Paramètres du modèle SARIMA, (p, d, q) (sp, sd, sq, s) pour Python Il n'y a pas de fonction qui en fait automatiquement la plus appropriée.

Par conséquent, la quantité d'informations standard(BIC dans ce cas(Norme de quantité d'information bayésienne) )Par
Vous devez écrire un programme pour savoir quelle valeur est la plus appropriée.

Je n'aborderai pas cette fois la norme de quantité d'informations Pour BIC, gardez à l'esprit que plus cette valeur est basse, plus la valeur du paramètre est appropriée.

#Données par ordre chronologique:DATA,Paramètres s(période):Entrez s pour afficher les meilleurs paramètres et leur BIC.

def selectparameter(DATA,s):
    p = d = q = range(0, 1)
    pdq = list(itertools.product(p, d, q))
    seasonal_pdq = [(x[0], x[1], x[2], 12) for x in list(itertools.product(p, d, q))]
    parameters = []
    BICs = np.array([])
    for param in pdq:
        for param_seasonal in seasonal_pdq:
            try:
                mod = sm.tsa.statespace.SARIMAX(DATA,
                                            order=param,
                seasonal_order=param_seasonal)
                results = mod.fit()
                parameters.append([param, param_seasonal, results.bic])
                BICs = np.append(y,results.bic)
            except:
                continue
    return print(parameters[np.argmin(BICs)])

Pour élaborer sur le contenu de la fonction, chaque paramètre est Pour 0,1 (cette fois, la limite supérieure du paramètre est fixée à 1 pour plus de simplicité) Le BIC du modèle SARIMA est calculé et le cas où le BIC est le plus petit est affiché.

Cependant, concernant les paramètres s, les données de séries chronologiques et Examinons-le en visualisant l'autocorrélation partielle décrite ci-dessous.

import warnings
import itertools
import pandas as pd
import numpy as np
import statsmodels.api as sm

#Chargement et organisation des données
sales_sparkling = pd.read_csv("./5060_tsa_data/monthly-australian-wine-sales-th-sparkling.csv")
index = pd.date_range("1980-01-31","1995-07-31",freq="M")
sales_sparkling.index=index
del sales_sparkling["Month"]

#Ce sont les données pour un an
sales_sparkling = sales_sparkling[:12]

#Définition des fonctions
def selectparameter(DATA,s):
    p = d = q = range(0, 2)
    pdq = list(itertools.product(p, d, q))
    seasonal_pdq = [(x[0], x[1], x[2], s) for x in list(itertools.product(p, d, q))]
    parameters = []
    BICs = np.array([])
    for param in pdq:
        for param_seasonal in seasonal_pdq:
            try:
                mod = sm.tsa.statespace.SARIMAX(DATA,
                                            order=param,
                                            seasonal_order=param_seasonal)
                results = mod.fit()
                parameters.append([param, param_seasonal, results.bic])
                BICs = np.append(BICs,results.bic)
            except:
                continue
    return parameters[np.argmin(BICs)]

#Mettez le cycle dans le deuxième argument
selectparameter(sales_sparkling,12)

Coefficient d'autocorrélation / coefficient d'autocorrélation partielle et sa visualisation

L'autocorrélation était la corrélation avec les données passées. Une autre façon d'analyser les données de séries chronologiques

La valeur de l'autocorrélation partielle est importante.

En parlant d'autocorrélation d'ordre k, il représente la corrélation entre yt et yt − k. L'autocorrélation partielle d'ordre k est les données entre yt et yt − k, c'est-à-dire La corrélation est obtenue en supprimant l'influence de yt − k + 1 de yt − 1.

Je vais vous expliquer en détail. Supposons que vous trouviez le coefficient d'autocorrélation pour une différence de 7 jours. Cependant, si les données d'un jour sont en corrélation avec les données de la veille, il y a 7 jours → il y a 6 jours → il y a 5 jours ... Il peut être corrélé à travers les données d'il y a 1 jour à aujourd'hui.

La corrélation obtenue en supprimant l'influence pendant cette période est appelée autocorrélation partielle.

Visualisons le coefficient d'autocorrélation et le coefficient d'autocorrélation partielle en utilisant Python. En général, si les données mensuelles sont saisonnières, le cycle sera de 12.

#Le graphique du coefficient d'autocorrélation est
sm.graphics.tsa.plot_acf(DATA)

#Le graphique de l'autocorrélation partielle est
sm.graphics.tsa.plot_pacf(DATA)

#Vous pouvez sortir avec.

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import statsmodels.api as sm
from pandas import datetime


#Chargement et organisation des données
sales_sparkling = pd.read_csv("./5060_tsa_data/monthly-australian-wine-sales-th-sparkling.csv")
index = pd.date_range("1980-01-31", "1995-07-31", freq="M")
sales_sparkling.index=index
del sales_sparkling["Month"]
#Visualisation de l'autocorrélation / autocorrélation partielle
fig=plt.figure(figsize=(12, 8))
#Produit un graphique du coefficient d'autocorrélation
ax1 = fig.add_subplot(211)
fig = sm.graphics.tsa.plot_acf(sales_sparkling, lags=30, ax=ax1)
#Produit un graphique du coefficient d'autocorrélation partielle
ax2 = fig.add_subplot(212)
fig = sm.graphics.tsa.plot_pacf(sales_sparkling, lags=30, ax=ax2)
plt.show()

Construction du modèle SARIMA (2)

Construire un modèle

Ici, je vais résumer la procédure d'analyse des données de séries chronologiques à l'aide du modèle SARIMA.

1. Lire les données

2. Organisez les données

3. Visualisation des données

4. Comprendre le cycle des données (déterminer les paramètres)

5. Détermination des paramètres autres que s

6. Construisez un modèle

7. Prédiction avec des données et leur visualisation

Cela se fait dans le flux.

Enfin, dans cette section, nous apprendrons de la construction du modèle à la prédiction. Une fois que vous connaissez l'optimum (p, d, q), (sp, sd, sq, s), il est temps de construire le modèle. Cette fois, nous utilisons des données sur 15 ans (180 mois).

#Pour construire un modèle
sm.tsa.statespace.SARIMAX(DATA,order=(p, d, q),seasonal_order=(sp, sd, sq, s)).fit()
#Est utilisé.

Construisons un modèle SARIMA de données sur les ventes de vins mousseux.

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import statsmodels.api as sm
from pandas import datetime


#Chargement et organisation des données
sales_sparkling = pd.read_csv("./5060_tsa_data/monthly-australian-wine-sales-th-sparkling.csv")
index = pd.date_range("1980-01-31", "1995-07-31", freq="M")
sales_sparkling.index=index
del sales_sparkling["Month"]

#Coupe du modèle
SARIMA_sparkring_sales = sm.tsa.statespace.SARIMAX(sales_sparkling,order=(0, 0, 0),seasonal_order=(0, 1, 1, 12)).fit()


#Sortie du BIC du modèle SARIMA construit
print(SARIMA_sparkring_sales.bic)

Exécution des prévisions et visualisation des données de prévision

Obtenir des données prédictives pour le modèle est facile

Nom du modèle.predict("Au début de la prédiction","À la fin de la prévision")Est utilisé.

Cependant, "au début de la prédiction" doit être l'heure dans les données de la série chronologique d'origine. Par exemple, les données sur les ventes de vins mousseux utilisées cette fois doivent être antérieures au 31/07/1995.

import pandas as pd
import matplotlib.pyplot as plt
import statsmodels.api as sm
from pandas import datetime


#Chargement et organisation des données
sales_sparkling = pd.read_csv("./5060_tsa_data/monthly-australian-wine-sales-th-sparkling.csv")
print(sales_sparkling.head())
print(sales_sparkling.tail())
index = pd.date_range("1980-01-31", "1995-07-31", freq="M")
sales_sparkling.index=index
del sales_sparkling["Month"]

#Coupe du modèle
SARIMA_sparkring_sales = sm.tsa.statespace.SARIMAX(sales_sparkling,order=(0, 0, 0),seasonal_order=(0, 1, 1, 12)).fit()

#Substituer les données de prédiction pour pred
pred = SARIMA_sparkring_sales.predict("1994-07-31","1997-12-31")

#Visualisation des données préd
plt.plot(pred)
plt.show()

Comparaison des données réelles et des données prévues

Ensuite, au lieu de produire uniquement des données de prévision Comparons les données de la série chronologique d'origine et les données de prévision en les générant en même temps.

import pandas as pd
import matplotlib.pyplot as plt
import statsmodels.api as sm
from pandas import datetime
import numpy as np


#Chargement et organisation des données
sales_sparkling = pd.read_csv("./5060_tsa_data/monthly-australian-wine-sales-th-sparkling.csv")
index = pd.date_range("1980-01-31", "1995-07-31", freq="M")
sales_sparkling.index=index
del sales_sparkling["Month"]

#Coupe du modèle
SARIMA_sparkring_sales = sm.tsa.statespace.SARIMAX(sales_sparkling,order=(0, 0, 0),seasonal_order=(0, 1, 1, 12)).fit()

#Substituer les données de prédiction pour pred
pred = SARIMA_sparkring_sales.predict("1994-7-31", "1997-12-31")

#Visualisation des données prédéfinies et des données de séries chronologiques originales
plt.plot(sales_sparkling)
plt.plot(pred, color="r")
plt.show()