Exercice 1

set.seed(123)
?arima.sim

Question (a)

Si \(d=0\) on a un ARMA. Si de plus \(d=0\), resp. \(q=0\), on trouve un MA, resp. AR.

Question (b)

\((1-B)^d X_t\) est la différentiation discrète d’ordre \(d\) de \(X_t\), notée \(\Delta^d X_t\) dans le cours. Ainsi, ARIMA(p,d,q) est un modèle ARMA(p,q) pour la série temporelle \(Y_t:=\Delta^d X_t\).

Question (c)

AR(1) de paramètre \(\phi_1=0.3\) et sa fonction d’autocorrélation :

x <- arima.sim(model=list(ar=0.3),n=100)
plot(x)

acf(x)

AR(1) de paramètre \(\phi_1=0.8\) et sa fonction d’autocorrélation :

x <- arima.sim(model=list(ar=0.8),n=100)
plot(x)

acf(x)

On observe que le chronogramme est oscille un peu moins que quand \(\phi=0.3\) : Il y a une plus forte corrélation (\(\sim 1\)) entre \(X_t\) et \(X_{t-1}\).

AR(1) de paramètre \(\phi_1=-0.8\) et sa fonction d’autocorrélation :

x <- arima.sim(model=list(ar=-0.8),n=100)
plot(x)

acf(x)

A l’inverse, on a ici une forte anticorrélation \((\sim -1)\) entre \(X_t\) et \(X_{t-1}\) : Si \(X_{t-1}\) est positive, il y a de fortes chances que \(X_{t}\) soit alors négative et réciproquement \(\Rightarrow\) Oscillations plus fréquentes.

Question (d)

Par définition, un processus \(X_t\) MA(q) satisfait : \[X_t=\sum_{k=0}^q\theta_k\;\varepsilon_{t-k}\,,\] où on a posé \(\theta_0:=1\) pour simplifier. En posant également \(\theta_k:=0\) pour \(k>q\), on calcule pour tout \(h\geq 0\) : \[ Cov(X_t,X_{t+h})=\sum_{k,\ell=0}^q\theta_k\theta_\ell \,Cov(\varepsilon_{t-k},\varepsilon_{t+h-\ell})=\sigma^2\sum_{k=0}^q\theta_k\theta_{k+h} \] et donc : \[ \rho(h)=\frac{\sum_{k=0}^q\theta_k\theta_{k+h} }{\sum_{k=0}^q\theta_k^2 }. \] En particulier, comme \(\theta_k=0\) dès que \(k>q\) par convention, on voit que \(\rho(h)=0\) dès que \(h>q\).

Question (e)

MA3 <- arima.sim(list(ma=c(0.9,0.6,0.9)),n=100)
plot(MA3)

acf(MA3)

Question (f)

C’est un modèle ARMA(1,2) avec un bruit blanc \(BB(0,1.5)\) gaussien. On voit dans l’aide de la commande \(\verb'arima.sim'\) que le bruit utilisé est par défaut gaussien, et on spécifie l’écart-type comme suit :

X <- arima.sim(n=200,list(ar=c(0.8),ma=c(-0.3,0.6)),sd=sqrt(1.5))
plot(X)

acf(X, lag.max=80)

ARIMA <- arima(X,order = c(1,0,2))
acf(ARIMA$residuals)

EXERCICE 2

Question (a)

Pour indexer la série temporelle par les trimestres en commençant par le deuxième trimestre de 1947, on fait comme ça :

gnpus_dat <- read.table("/Users/adrien/R/Data/gnpus.txt")
gnpus <- ts(gnpus_dat, frequency=4, start=c(1947,2))
plot(gnpus)

Question (b) : Test de blancheur

Box.test(gnpus)
## 
##  Box-Pierce test
## 
## data:  gnpus
## X-squared = 24.998, df = 1, p-value = 5.74e-07

On rejète donc l’hypothèse \(H_0\) d’un bruit blanc vu la faible valeur de la p-valeur.

Question (c) : Fonction d’autocorrélation empirique

acf(gnpus)

Vu l’ACF, on peut penser que cette série temporelle ressemble à un modèle AR, voir un MA(2).

Question (d) : Modélisation avec un AR

gnpusAR <- ar(gnpus)
gnpusAR
## 
## Call:
## ar(x = gnpus)
## 
## Coefficients:
##       1        2        3  
##  0.3463   0.1770  -0.1421  
## 
## Order selected 3  sigma^2 estimated as  9.676e-05

Cette fonction nous propose donc un AR(3), de modèle : \[ X_t=0.3463X_{t-1}+0.177X_{t-2}-0.1421 X_{t-3}+\epsilon_t \] où \(\epsilon_t\) est un bruit blanc \(BB(0,\sigma^2)\) de variance \(\sigma^2=9.676\times 10^{-05}\).

Question (e) : Etude des résidus

Residus <- gnpusAR$resid
plot(Residus)
library(caschrono) # ce package contient la fonction Box.test.2
## Loading required package: zoo
## 
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
## 
##     as.Date, as.Date.numeric

Box.test.2(Residus,1:40) 
##       Retard   p-value
##  [1,]      1 0.9115127
##  [2,]      2 0.9818001
##  [3,]      3 0.9415294
##  [4,]      4 0.8416601
##  [5,]      5 0.8444775
##  [6,]      6 0.9165832
##  [7,]      7 0.9214457
##  [8,]      8 0.8904874
##  [9,]      9 0.7792795
## [10,]     10 0.7562754
## [11,]     11 0.8220595
## [12,]     12 0.7840702
## [13,]     13 0.8329639
## [14,]     14 0.8165586
## [15,]     15 0.8236042
## [16,]     16 0.8480766
## [17,]     17 0.8795128
## [18,]     18 0.8862396
## [19,]     19 0.9173663
## [20,]     20 0.9402917
## [21,]     21 0.9525794
## [22,]     22 0.9637501
## [23,]     23 0.9624040
## [24,]     24 0.9659288
## [25,]     25 0.9745159
## [26,]     26 0.9811298
## [27,]     27 0.9846522
## [28,]     28 0.9832167
## [29,]     29 0.9867008
## [30,]     30 0.9557774
## [31,]     31 0.9660820
## [32,]     32 0.9743399
## [33,]     33 0.9757968
## [34,]     34 0.9389340
## [35,]     35 0.9531035
## [36,]     36 0.9620963
## [37,]     37 0.9621763
## [38,]     38 0.9484182
## [39,]     39 0.9601293
## [40,]     40 0.9684611
acf(Residus,na.action=na.pass)

shapiro.test(Residus)
## 
##  Shapiro-Wilk normality test
## 
## data:  Residus
## W = 0.98906, p-value = 0.2031

Au vu de l’acf et du test de Box-Pierce effectué pour un lag \(h\) variant de \(1\) à \(40\), qui ne rejète jamais \(H_0\), il n’est pas aberrant de supposer que les résidus forment un bruit blanc. De plus, le test de Shapiro-Wilk ne rejète pas non plus l’hypothèse de normalité des résidus. On peut donc conclure que la modélisation AR(3) ci-dessus est raisonnable, et que l’on peut même supposer le bruit gaussien. Cette dernière hypothèse nous permet d’obtenir un intervalle de confiance pour de la prédiction. En effet, ayant mesuré \(x_t,x_{t-1},x_{t-2}\), on peut faire la prédiction : \[ \tilde x_{t+1}=0.3463x_{t}+0.177x_{t-1}-0.1421 x_{t-2} \] avec un intervalle de confiance à \(95\%\) pour la vraie valeur \(x_t\) donné par : \[ [\tilde x_{t+1}-1.96\sigma,\tilde x_{t+1}+1.96\sigma]\qquad \mbox{où}\qquad \sigma=\sqrt{9.676\times 10^{-05}}\simeq 0.0098. \]