In [1]:
import numpy as np
import numpy.random as npr
import matplotlib.pyplot as plt
import pandas as pd
from pandas import read_csv
from statsmodels.tsa.arima_process import ArmaProcess
from statsmodels.graphics.tsaplots import plot_acf
from statsmodels.tsa.arima_model import ARMA
from statsmodels.stats.diagnostic import acorr_ljungbox
import statsmodels.api as sm
# from statsmodels.tsa.seasonal import seasonal_decompose as decompose
# from statsmodels.tsa.stattools import acf
In [15]:
phi_1 = 0.5
ar_parms = [1,-phi_1] # Attention, les paramètres AR sont [1,-phi_1,...,-phi_p] et [1,theta_1,...,theta_q] pour MA
ma_parms = [1]
AR_t = ArmaProcess(ar=ar_parms,ma=ma_parms).generate_sample(300)
plt.plot(AR_t)
plt.show()
In [16]:
# ACF theorique
plot_acf(ArmaProcess(ar=ar_parms,ma=ma_parms).acf())
plt.show()
In [17]:
# ACF empirique
plot_acf(AR_t,lags=20)
plt.show()
Question (d)¶
In [33]:
ar_parms = [1] # C'est un MA(3)
ma_parms = [1,0.9,0.6,0.9]
MA_t = ArmaProcess(ar=ar_parms,ma=ma_parms).generate_sample(100)
plt.plot(MA_t)
plt.show()
In [34]:
# ACF theorique
plot_acf(ArmaProcess(ar=ar_parms,ma=ma_parms).acf())
plt.show()
In [35]:
# ACF empirique
plot_acf(MA_t,lags=20)
plt.show()
Question (e)¶
In [24]:
ar_parms = [1,-0.8] # C'est un ARMA(1,2)
ma_parms = [1,-0.3,0.6]
X_t = ArmaProcess(ar=ar_parms,ma=ma_parms).generate_sample(200)
plt.plot(X_t)
plt.show()
In [72]:
model = ARMA(X_t,order=(1,2)).fit(trend='nc')
model.summary()
Out[72]:
In [75]:
# Les coefficients estimés ne sont pas loin des valeurs théoriques : AR = 0.8, MA = -0.3, 0.6
model.params
Out[75]:
In [73]:
residus = model.resid
plt.plot(residus)
plot_acf(residus)
plt.show()
In [44]:
acorr_ljungbox(residus,lags=range(60),boxpierce=True)[3]
# p-valeurs des 60 premiers lags pour le test de Box-Pierce
Out[44]:
In [45]:
acorr_ljungbox(residus,lags=range(60),boxpierce=True)[1]
# p-valeurs des 60 premiers lags pour Ljung-Box test
Out[45]:
In [47]:
acorr_ljungbox(residus,lags=range(60),boxpierce=True)[3] < 0.01
Out[47]:
L'étude des résidus laisse penser que la modélisation est bonne : d'après le test de Box-Pierce (et Ljung-Box) l'hypothèse $H_0$ que les résidus sont gaussien n'est pas rejetée.
Exercice 2¶
In [76]:
GNP_US = read_csv('gnpus.csv')
In [77]:
type(GNP_US)
Out[77]:
In [78]:
plt.plot(GNP_US)
plt.show()
In [79]:
# Plot de l'ACF (avec statsmodels.graphics.tsaplots.plot_acf)
plot_acf(GNP_US,lags=12)
plt.show()
In [80]:
#Après tatonnements, on choisit p=4 pour que l'hypothèse {résidus=BB} ne soit pas rejetée
AR_GNPUS = ARMA(X_t,order=(4,0)).fit(trend='nc')
AR_GNPUS.summary()
Out[80]:
In [81]:
plot_acf(AR_GNPUS.resid)
plt.plot()
Out[81]:
In [84]:
acorr_ljungbox(AR_GNPUS.resid,lags=range(60),boxpierce=True)[3] < 0.01
Out[84]: