Tag Archives: intercept

La maudite constante des modèles ARIMA

Dans les modèles ARIMA, autant le dire tout de suite, les constantes c’est pénible ! Pour comprendre un peu mieux ce qui se passe, considérons ici un modèle ARMA avec constante. Ou pour commencer, juste un AR(1)

Supposons que la série  soit stationnaire, de moyenne . Alors en prenant l’espérance de part et d’autre dans l’expression précédente, , i.e. , ou . Plus généralement, avec un ARMA, , en prenant la encore l’espérance des deux cotés, on en déduit .

Simulons un processus AR(1),

Pour simuler des processus ARIMA, on pourrait utiliser la commande

> X=arima.sim(list(order=c(1,0,0),ar=1/3),n=1000)+2
> mean(X)
[1] 1.931767

mais comme le but est de comprendre ce qui se passe, autant faire les choses calmement,

> X=rep(NA,1010)
> X[1]=0
> for(t in 2:1010){X[t]=4/3+X[t-1]/3+rnorm(1)}
> X=X[-(1:10)]
> mean(X)
[1] 2.03397

I.e. ici le processus est de moyenne qui vaut 2.


Regardons maintenant ce que donnerait l’estimation du processus AR(1),

> arima(X, order = c(1, 0, 0))

Call:
arima(x = X, order = c(1, 0, 0))

Coefficients:
ar1  intercept
0.3738     2.0334
s.e.  0.0294     0.0487

sigma^2 estimated as 0.9318:  log likelihood = -1383.68

De manière un peu surprenante, le coefficient appelé intercept n’est pas la constante dans le modèle AR(1), mais la moyenne du processus. Autrement dit, R n’ajuste pas un processus

comme le laisserait penser l’intuition, mais un processus

Ces deux formes sont bien entendu équivalentes. Mais les coefficients estimés ne sont pas tout à fait ce que l’on attendait…
Plaçons nous maintenant dans le cas d’un processus non-stationnaire. L’extension naturelle serait de considérer un processus ARIMA(1,1,0), ou bien un processus tel que la différence soit un processus AR(1). Un processus ARIMA(1,1,0) avec une constante s’écrirait

en utilisant l’opérateur retard. Ceci fait penser à un modèle avec une tendance linéaire. Posons  (afin de se débarrasser de cette tendance). Alors

i.e.

ou encore

On note que  ou encore  suit alors un processus ARIMA(1,1,0) sans constante cette fois.
Afin de visualiser ce que donnerait l’inférence, simulons le processus suivant, qui est un processus AR(1) avec constante que l’on intègre:  avec

Peut-on retrouver les différents paramètres du processus avec R ?
Commençons (là encore) par simuler un tel processus,

> U=rep(NA,1010)
> U[1]=0
> for(t in 2:1010){U[t]=4/3+U[t-1]/3+rnorm(1)}
> U=U[-(1:10)]
> X=cumsum(U)

Sous R, on obtient l’estimation suivant si l’on tente de calibrer un modèle ARIMA(1,1,0)

> arima(X, order = c(1, 1, 0))

Call:
arima(x = X, order = c(1, 1, 0))

Coefficients:
ar1
0.8616
s.e.  0.0160

sigma^2 estimated as1.343:loglikelihood = -1565.63Ça ne convient pas du tout…. On peut tenter un processus AR(1) (avec constante) sur la série différenciée…

> arima(diff(X), order = c(1, 0, 0))

Call:
arima(x = diff(X), order = c(1, 0, 0))

Coefficients:
ar1  intercept
0.3564     2.0200
s.e.  0.0295     0.0486

sigma^2 estimated as 0.9782:  log likelihood = -1406.6

On progresse, sauf que comme auparavant, le terme qui est donne n’est pas la constante dans le modèle ARIMA, mais la moyenne du processus différencié. Mais cette fois, on a un interprétation, c’est que la constante est la pente de la tendance ! Si on estime la pente associée a , on récupère la même valeur:

> arima(X, order = c(1, 1, 0), xreg=1:length(X))

Call:
arima(x = X, order = c(1, 1, 0), xreg = 1:length(X))

Coefficients:
ar1  1:length(X)
0.3566       2.0519
s.e.  0.0296       0.0487

sigma^2 estimated as 0.9787:  log likelihood = -1406.82

Sur la figure ci-dessous, on retrouve le fait qu’en enlevant la tendance linéaire à  donneune série intégrée, sans constante (qui fait penser à une marche aléatoire).

Autrement dit

  • dans le cas d’une série stationnaire, la constante estimée n’est pas du tout la constante, mais la moyenne de la série temporelle
  • dans le cas d’une série non-stationnaire, la constante estimée dans la série différenciée a du sens, au sens ou il s’agit de la pente de la tendance (linéaire) du processus

Avant de conclure, une petite remarque. Quid de la prévision ? Si on commencer par une révision à l’aide du premier processus ARIMA, on obtient une prédiction pour avec un énorme intervalle de confiance (sans aucun bon sens)

> ARIMA1=arima(X, order = c(1, 1, 0))
> ARIMA2=arima(X, order = c(1, 1, 0), xreg=1:length(X))
> Xp1=predict(ARIMA1,20)
> Xp2=predict(ARIMA2,20,newxreg=
+ (length(X)+1):(length(X)+20))
> plot(960:1000,X[960:1000],xlim=c(960,1020),type="l")
> polygon(c(1001:1020,rev(1001:1020)),
+ c(Xp1$pred+2*Xp1$se,rev(Xp1$pred-2*Xp1$se)),
+ col=CL[3],border=NA)
> lines(1001:1020,Xp1$pred,col="red",lwd=2)

intervalle qui se visualise sur le graphique ci-dessous,

Si on regarde pour l’autre modèle,

> lines(1001:1020,Xp2$pred,col="blue",lwd=2)

Tous ceux qui auront reconnu ici des problèmes qui se posent dans la modélisation de la série http://freakonometrics.blog.free.fr/public/maths/viekt.png dans le modèle de Lee & Carter (1992) auront compris qu’il existe une vraie application a tout ce que je viens de raconter (je pense au particulier au commentaire poste par @ClaudeT en début de semaine) car la série des http://freakonometrics.blog.free.fr/public/maths/viekt.png ressemble très fortement à ce genre de série, non-stationnaire avec une tendance linéaire. Et qu’un modèle mal spécifié laisse à penser que l’incertitude est beaucoup beaucoup plus grande que ce qu’elle est vraiment (en plus de donner une prédiction surprenante à long terme !). Donc avant de faire tourner rapidement les codes, il vaut mieux regarder attentivement ce qu’ils font réellement…