Séries chronologiques

October 1, 2007
Chapter 1

Introduction

nombre d’+observation+s N = ms = 40
nombre d’+observations N = ms = 40

On ne peut pas évider les séries temporelles. Exemples :

· Nombre de chômeurs en France

· Température dans le coeur d’un réacteur nucléaire.

· Précipitation quotidienne

· Fuel consommation trimestrielle

1.1 Décomposition d’une série temporelle

Les séries temporelle sont souvent composées de 3 parties

1. Tendance

2. Cycle saisonnier

3. Perturbation aléatoire

X(t) = T (t) + S(t) + R(t)

Identification d’un système : Identifier les composants T (t) et S(t). On examine dans la suite un exemple concret
Fuel consommation exemple La table 1.1 représente la quantité (en m3 ) de fuel consommée dans une ville
pendant 5 ans. La figure 1.2 trace le trajectoire (observation), on remarque facilement que il présente une période et

année 1 2 3 4
1965 874 679 616 816
1966 866 700 603 814
1967 843 719 594 819
1968 906 703 634 844
1969 952 745 635 871
moyenne 888.2 709.2 616.4 835.8

Table 1.1: Fuel consommation

une tendance et il est perturbée par une variable aléatoire.

1.2 Identification du système, première méthode

On cherche à estimer les 3 composants T (t), S(t) et R(t).
 1000

950

900

850

800

750

700

650

600

550
0 5 10 15 20

Figure 1.1: Fuel consommation

1.2.1 Identification du composant saisonnier

On désigne

nombre de cycle m = 5, période s = 4,

nombre d’observations N = ms = 40
N
1 X
la moyenne totale M = X(t) = 761.65
N t=1

Première méthode Une méthode très simple pour trouver S(1), S(2), S(3) et S(4) est de calculer la différence entre
moyenne totale et la moyenne partielle

S(1) = 882.2 − 761.65 = 126.55, S(2) = 709.2 − 761.65 = −52.4

S(3) = 616.4 − 761.65 = −145.25, S(4) = 832.8 − 761.65 = 71.15

Soit t = k + sj

X(k + sj) = T (k + sj) + S(k + sj) + R(k + sj) j = 0, 1, 2, 3, 4 k = 1, 2, 3, 4

En soustrayant le composant S(k + sj) on obtient la figure 1.3. Cette méthode marche bien en générale quand le
composant tendanciel est une fonction linéaire. Soit T (t) = a + bt (b: gradient pente, ici très faible, a, b coefficients
inconnus). On examine maintenant le défaut de cette méthode, en ignorant le composant aléatoire, la moyenne totale
est (Attention : en effet, la moyenne sur R(t) est proche de 0)

N N
1 X 1 N (N + 1) 1 X
M̄ = (T (t) + S(t)) = (N a + b )+ S(t)
N j=1 N 2 N j=1

N +1
= a+ b
2
 1000
900
"stat" using 1:2
"stat" using 1:3
"stat" using 1:4

800
700
600
500
400
300
200
100
0
0.5

0.45

0.4

0.35

0.3

0.25

0.2

0.15

0.1

0.05

0
Figure 1.2: Fuel consommation

PN
j=1 S(j) = 0. La moyenne partielle correspond à k est

m−1
1 X
M̄k = (T (k + sj) + S(k + sj))
m j=0

m−1
1 X
= (a + (k + sj)b + S(k + sj))
m j=0

m−1
m−1 1 X
= a + bk + bs + S(k + sj)
2 m j=0

La méthode précédente consiste donc de soustraire

m−1
m−1 1 X ms + 1
S̄(k) = M̄ (k) − M̄ = a + bk + bs + S(k + sj) − (a + b)
2 m j=0 2

Pm−1
= bk − b s+1
2 +
1
m j=0 S(k + sj)

Donc si la pente b est petite, alors la méthode marche bien.

Deuxième méthode On effectue des calculs successifs pour obtenir la table 1.2. Grâce à la dernière colonne (moyenne
mobile centrée) on peut estimer le composant saisonnier par

1
S(3) = ((616 − 745.25) + (603 − 741.75) + (594 − 750.) + (634 − 777.7)) = −142.00
4

et S(1) = 128.594 S(2) = −46.487, P S(4) = 65.00. On remarque que S(1) + S(2) + S(3) + S(4) = 5.107 (non zéro).
ceci peut être réajusté pour que S(i) = 0

5.107 5.107
S(1) = 128.594 − 4 = 127.74 S(2) = −46.487 − 4 = −47.763

5.107 5.107
S(3) = −142.00 − 4 = −143, 254 S(4) = 65.00 − 4 = 63.746

Le résultat est tracé dans la figure 1.4. Cette méthode est facile à implementer. A noter que dans la 3ème colonne il
manque des données CMA au début et à la fin de la série. Cette méthode est assez efficace dans le cas ou b est petit.
 1000

950

900

850

800

750

700

650

600

550
0 5 10 15 20

Figure 1.3: composant saisonnier éliminé

1.2.2 Identification du composant aléatoire

Nous avons dans cet exemple estimé le composant saisonnier (S(1), S(2), S(3), S(4)) et la tendance (sans calculer les
coefficients a et b). On peut alors, grâce à l’observation, calculer le composant aléatoire. Le résultat obtenu est
représenté par la figure 1.5. Dans les applications, il n’y a très peu d’intérêt à trouver le composant aléatoire, car
les perturbations aléatoires changent constamment, seuls T (t) et S(t) sont importants. On trace maintenant ces trois
composants sur la même figure 1.6.

1.2.3 Identification du composant tendanciel

On suppose que la tendance est linéaire
T (t) = a + bt

Par une méthode d’approximation 1.7 classique on obtient l’estimation

T̂ = 713.38 + 3.647t

1.3 Identification du système, deuxième méthode

On présente ici une méthode alternative, mais sans entrer en détails. Soit

Xt = β0 + β1 t + β2 d1 + β3 d2 + β4 d3 + β5 d4 + ε(t)

où β0 + β1 t représente le composant tendanciel, β2 d1 + β3 d2 + β4 d3 + β5 d4 le composant saisonnier avec


 1 ième saison
di =

0 sinon

En appliquant la méthode de régression, on obtient la solution suivante

Xt = 3.38t + 858.83d1 + 675.45d2 + 579.28d3 + 792.30d4 + ε(t)

Dans la figure 1.8 on compare cette méthode avec l’observation.

 saisons X(t) MA CMA
1 874

2 679
746.250
3 616 745.250
744.250
4 816 746.875
749.500
5 866 747.875
746.250
6 700 746.000
745.750
7 603 741.750
737.750
8 814 740.125
742.500
9 834 741.375
740.250
10 719 740.875
741.500
11 594 750.500
759.500
12 819 757.500
755.500
13 906 760.500
765.500
14 703 768.625
771.750
15 634 777.500
783.250
16 844 788.500
793.750
17 952 793.875
794.000
18 745 797.375
800.750
19 635

20 871

Table 1.2: calcul de moyenne mobile

1.4 Prédiction
Identification du système, pourquoi ? Une application concrète est la prédiction. En fait, d’après les calculs précédents
on a :
Xt = T (t) + S(t) + R(t)
et
T̂ (t) = 713.38 + 3.647t

Ŝ(1) = 127.34, Ŝ(2) = −46.763, Ŝ(3) = −143.254, Ŝ(4) = 63.746

Question : quelle est la quantité consommée pour l’année 1970 ? Pour répondre à cette question, on peut naturellement
utiliser la prédiction suivante, le résultat est dans table 1.3
X̂t = T̂ (t) + Ŝ(t), t = 21, 22, 23, 24
La qualité de cette prédiction est tracée dans la figure 1.9.

période tendance saisonnier prédiction observation

t T S X̂ X
21 789.963 127.340 917.303 981
22 793.610 -47.763 745.847 759
23 797.257 -143.254 654.003 674
24 800.904 63.746 864.650 900

Table 1.3: Prédiction

1.5 Études préliminaire d’une série

1.5.1 Les données

Xt , t = 1, 2, · · · , n n : taille échantillon
 1000

950

900

850

800

750

700

650

600

550
0 5 10 15 20

Figure 1.4: Fuel consommation et Moyenne mobile centrée

série observée. Ces observations sont faites à intervalle de temps régulière. A noter que cet échantillon n’est pas issu
d’un n-uplet de variables aléatoires indépendantes. Voici quelques exemples des séries temporelles

· Relevés mensuels d’un compte bancaire

· longueurs journalière de la tige d’une plante

· Évaluation d’heur en heur du stock pour chaque article d’un hypermarché.

1.5.2 Observation de graphique

On n’hésitera pas, avant et pendant l’analyse, à user de graphique. Par simple observation, on pourra s’apercevoir
par exemple

– Existence de valeurs manquantes

– Valeurs accidentelles
– Mouvement tendanciel
– mouvement périodique
– ruptures de comportement
– analogie de comportement de plusieurs séries représentées simultanément.

On peut aussi observer les données transformées.

Xtλ +a
– b où a ∈ R, b ∈ R, λ ∈ R − {0} (Box-Cox)
– log Xt choisit pour atténuer une croissance exponentielle.
Xt+1 −Xt
– Xt taux de croissance

1.5.3 Corrections des données

Modifier les données brutes avant de se lancer dans une étude. Quelques taches possibles.

– Évaluation préliminaire de données manquantes, remplacement de données accidentelles.

– Découpage en sous séries : rupture de comportement
– Standardisation(harmonisation), calcul de moyenne pour se ramener à un intervalle de temps d’amplitude
fixe entre les relevés successifs. Par exemple le mot mois) doit être normalisé selon le contexte.
 50

40

30

20

10

-10

-20

-30
0 5 10 15 20

Figure 1.5: Composant irrégulier

∗ en consommation 365/12 jours

∗ en production (365 − jours fériés)/12
– Transformations des données pour des raison divers
∗ Cartésienne vers polaire
∗ log.
∗ etc.
 800

600

400

200

-200
0 5 10 15 20

Figure 1.6: Décomposition

Figure 1.7: Estimation de la tendance

1000

950

900

850

800

750

700

650

600

550
0 5 10 15 20

Figure 1.8: Modèle de régression linéaire et observation réelle

1000

950

900

850

800

750

700

650

600

550
0 5 10 15 20 25

Figure 1.9: Prédiction et observation réelle

Chapter 2

Modèles des séries temporelles

La résolution des divers problèmes s’appuie sur des modèles décrivant la façon dont la dérive évolue. Il est commode
de distinguer trois types de modèles: les modèles d’ajustement, les modèles auto-projectifs et les modèles explicatifs.
Nous allons présenter brièvement ces trois modèles.

2.1 Modèles d’ajustement

D’une façon générale, les modèles d’ajustement sont sous la forme

Xt = f (t, εt )

où f est une fonction indexée par un nombre fini de paramètres inconnus et où εt est une v. a centrée sur la quelle
diverses hypothèses additionnelles peuvent être faites. Il existe en pratique 3 modèles généraux dans cette catégorie
de série temporelle, articulant les composants de façons différentes. Notations:

Xt : Chronique donnée observation

Tt : tendance à l’instantt
Vt : chronique de variation saisonnière
ε(t) : chronique de variation accidentelle(R(t))
St : cycle saisonnière.

2.1.1 Modèle Additif

Le modèle additif consiste à poser
X t = T t + St + ε t
Le modèle le plus simple suppose que le saisonnier de période p est constant

St+p = St , ∀t

Il est alors déterminé par p coefficients c1 , c2 , · · · , cp . On impose la contrainte suivante :

p
X
ci = 0
i=1

Example 2.1.1 Soit

Tt = 0.05t + 100

S1 = −1, S2 = 2, S3 = −2, S4 = 1

εt ∼ N (0, 0.52)

Voici la simulation de cet exemple présentée dans la figure 2.1.

De toutes façon, il y a une infinité de solutions au problème de décomposition, correspondant d’abord à l’infinité
de choix possible de tendance

Xt − T t = est la série détendancialisée, elle est égale à saisonnier + résidu

X t − St = est la série désaisonnalisée, elle est égale à tendance + résidu.

on la nomme aussi série corrigée des variations saisonnières
 107

106

105

104

103

102

101

100

99

98

97
0 10 20 30 40 50 60 70 80 90 100

Figure 2.1: Modèle additif

2.1.2 Premier modèle multiplicatif

le deuxième modèle est multiplicatif : il suppose que les variations saisonnières et le facteur accidentel ont un effet
multiplicatif sur la taille des observations. On le note
Xt = Tt (1 + St )(1 + Rt )

= Tt + Tt St + (Tt + Tt St )Rt = Tt + Vt + εt
Dans la mesure où les termes sont positifs, on peut calculer le logarithmes
ln Xt = ln Tt + ln(1 + St ) + ln(1 + Rt )
Nous prenons comme exemple le modèle suivant
Xt = ea+bt (1 + St )(1 + Rt )
avec
a = ln 100, b = 0.05

S1 = 0.06, S2 = −0.06, S3 = −0.04, S4 = 0.04

Rt ∼ N (0, 202 )
Le résultat de simulation est présenté dans la figure 2.2. On distingue encore nettement la périodicité des variations
saisonnières, mais la différence Xt − Tt s’accroı̂t avec le temps. On trace maintenant cet exemple dans la figure 2.3

2.1.3 Deuxième modèle multiplicatif

Le troisième modèle est
Xt = Tt (1 + St ) + εt
Par rapport au modèle précédent, seule la variation saisonnière Tt St est proportionnelle à la tendance, on peut alors
écrire
Xt = Tt (1 + St ) + εt
Voici un exemple dont le résultat simulé est présenté dans la figure 2.4
Tt = 0.05t + log(100)

S1 = 0.06, S2 = −0.06, S3 = −0.04, S4 = 0.04

εt ∼ N (0, 0.12)
 18000

16000

14000

12000

10000

8000

6000

4000

2000

0
0 10 20 30 40 50 60 70 80 90 100

Figure 2.2: Premier modèle multiplicatif

Est-ce que la variation accidentelle est proportionnelle ou non à la tendance ? Il est difficile de donner une réponse à
cette question uniquement à partir d’une simple représentation graphique. En effet on distingue entre les deux modèles
multiplicatifs suivant les besoins du modèle et de la technique de description employée. Comment alors définir un
modèle ?
1. Graphique pour déterminer les facteurs explicatifs de son évolution.
2. Méthode numérique (chapitres suivants)
Ces deux méthodes sont deux approches complémentaires.

2.2 Modèles auto-projectifs

Dans un modèle auto-projectif, on suppose que Xt est une fonction des ses valeurs passées et d’une perturbation
aléatoire εt
Xt = f (Xt−1 , Xt−2 , · · · , εt )
Une classe de tels modèles, particulièrement utiles pour la prévision, sont les modèles ARIMA. Comme on verra, cette
classe permet d’atteindre à l’aide d’un nombre de paramètres relativement limité une gamme de modèles très variés.
En outre, il est possible de proposer des méthodes dites méthodes d’identification, permettant de choisir dans cet
ensemble de modèles celui qui semble le mieux adapté aux données dont on dispose. Une fois ce modèle choisi, on
peut estimer les paramètres et déterminer les prévisions optimales.
Parmi les méthodes de prévision ainsi mises en évidence on trouve, en particulier, les méthodes de lissage exponen-
tiel, l’avantage essentiel de la méthode (identification, estimation, prévision) fondée sur les modèles ARIMA, appelée
aussi méthode de BOX et Jenkins, est qu’elle permet de sélectionner la méthode de prévision optimale dans un vaste
ensemble possibilités. tandis que les méthodes classiques (ajustement et lissage) comportent une part plus grande
d’arbitraire dans le choix de la fonction utilisée pour l’ajustement.

2.3 Modèles explicatifs

Dans cette dernière catégorie de modèles la variable Xt est exprimée en fonction d’un vecteur de variables observable,
dite exogènes, Yt et d’une perturbation aléatoire εt :

Xy = f (Yt , εt )

Où Yt est soit déterministe, soit aléatoire; dans ce dernier cas, les processus {Yt } et {ut } ont certaines propriétés
d’indépendance ou de non corrélation.
 10

4
0 10 20 30 40 50 60 70 80 90 100

Figure 2.3: Log du premier modèle multiplicatif

Ces modèles sont les méthodes de base de l’économétrie et nous les considérons essentiellement pour faire le lien
entre eux et les modèles auto-projectifs.

2.3.1 Modèle explicatif statique

Comme son nom indique, dans un modèle explicatif statique les variables Yt ne contiennent pas de valeurs passées de
Xt et les εt sont indépendants entre eux. Par exemple

Xt = a + bYt + εt , t = 1, · · · , T

où les Yt sont indépendantes de tous les εt et où les εt sont centrés et indépendants entre eux.

2.3.2 Modèle explicatif dynamique

Un modèle explicatif peut être dynamique soit parce que les εt sont auto-corrélées, soit parce que Yt contient des
variables passées de Xt , c’est à dire des variables dite endogènes retardées. On distingue donc deux classes différentes
· Perturbation auto-corrélées
· Variables endogènes retardées

2.3.3 TP: simulation des chroniques

Dans cette section, on simule quelques séries chronologiques en utilisant le langage C++.

TP: simulation
Les séries temporelles peuvent être simulées par ordinateur.
Exercice 2.3.1 Soient X ∼ U [0, 1] et Y ∼ U [0, 1] deux variables aléatoires indépendantes. Calculer la loi de la
variable aléatoire Z = X − Y . Montrer qu’il s’agit d’une v. a. centrée.

Exercice 2.3.2 Variable aléatoire N (0, 1)

Soit X1 , X2 , · · · , Xn · · · des variables aléatoires i.i.d. et Xk ∼ U (0, 1) (loi uniforme sur (0, 1)). Définissons
1
Y2n−1 = (−2 ln X2n−1 ) 2 cos 2πX2n ,
1
Y2n = (−2 ln X2n−1 ) 2 sin 2πX2n .

Montrer que Y1 , Y2 , · · · , Yn , · · · sont i.i.d. et Yk ∼ N (0, 1).

 16000

14000

12000

10000

8000

6000

4000

2000

0
0 10 20 30 40 50 60 70 80 90 100

Figure 2.4: Deuxième modèle multiplicatif

Solution : Le résultat de cette exercice nous donne une méthode pour simuler les suites des variables aléatoires
indépendantes gaussiennes N (0, 1). Il est donc très utile pour la simulation numérique. On cherche donc à montrer :
Soit X et Y indépendantes de loi U [0, 1], Alors le couple (U, V ) fabriqué par
1
U = (−2 ln X) 2 cos 2πY,
1
V = (−2 ln X) 2 sin 2πY.

est indépendantes, et suivent la même loi N (0, 1) Proof. On calcule d’abord la loi de (U, V ). Soit h(u, v) une fonction
mesurable. On introduit un changement de variable
  2 
  u + v2
1
 u = (−2 ln x) 2 cos 2πy,

 x = exp −
 2
⇐⇒
 1   
v = (−2 ln x) 2 sin 2πy.  y = 1 arctan v


2π u
On calcule Jacobien
 2   2 
u + v2 u + v2
∂x1 ∂x1 −u exp − −v exp −
∂u ∂v 2 2 1 u2 + v 2
|J| = = =− exp{− }
∂x2 ∂x2 1 v 1 u 2π 2
∂u ∂v −
2π u2 + v 2 2π u2 + v 2
Donc  
1 1
E h(U, V ) = E h (−2 ln X) 2 cos 2πY, (−2 ln X) 2 sin 2πY

Z 1 Z 1  
1 1
= h (−2 ln x) 2 cos 2πy, (−2 ln x) 2 sin 2πy dxdy
0 0
Z Z
1 u2 + v 2
= h(u, v) exp{− }dudv
R R 2π 2
Donc
1 u2 + v 2
p(u, v) = exp{− }
2π 2
   
1 u2 1 v2
= √ exp{− } × √ exp{− }
2π 2 2π 2
    
1 u2 1 v2
On vérifie facilement que √ exp{− } et √ exp{− } sont respectivement la densité de U et V . En
2π 2 2π 2
prenant dans la démonstration une fonction h(u, v) qui dépend uniquement de u. En effet
 1

E h(U ) = E h (−2 ln X) 2 cos 2πY

Z 1 Z 1  
1
= h (−2 ln x) 2 cos 2πy dxdy
0 0
Z Z
1 u2 + v 2
= h(u) exp{− }dudv
R R 2π 2
Z Z
1 v2 1 u2
= √ exp{− }dv h(u) √ exp{− }du
R 2π 2 R 2π 2
Z
1 u2
= h(u) √ exp{− }du
R 2π 2

Exercice 2.3.3 TP : simulation des chroniques. On peut simuler à l’aide de ordinateur les série suivantes


 −1

2
Xt = 0.05t + 100 + St + εt où St = εt ∼ N (0, 0.52 ) (2.1)

 −2

1


 0.3

a+bt −0.3
Xt = e (1 + St )(1 + Rt ) où St = Rt ∼ N (0, 202) (2.2)

 0.2

−0.2
avec a = ln 100 et b = 0.05
Tracer ln Xt . 

 −1

2
Xt = (0.05t + 100)(1 + St ) + εt où St = εt ∼ N (0, 1) (2.3)

 −2

1
2.3.4 Implémentation
On construit d’abord une classe BBG:
/* ------ BBG.h ---------------*/
#ifndef BBG_H_IS_DEFINE
#define BBG_H_IS_DEFINE
class BBG {
double W1, W2;
int seesaw;
public:
double gaussien();
BBG(int);
};
#endif
/* ------ fin BBG.h -----------*/
Voici les membres fonctions pour BBG.h
/* ------ BBG.C ---------------*/
#include <stdlib.h> /* drand48(), srand48(), sin(), cos() */
#include <math.h> /* M_PI, log() */
#include <iostream.h> /* cout */
#include "BBG.h"

BBG::BBG(int seed) {
seesaw = 0;
srand48(seed);
}

double BBG::gaussien()
{
double c1,c2;

if (seesaw==0)
{
c1=2.0*M_PI*drand48();
c2=-2.0*log(drand48());
c2=sqrt(c2);
W1 = c2*cos(c1);
W2 = c2*sin(c1);
seesaw=1;
// cout <<"seesaw = "<<seesaw<<" "<<W1<<"\n ";
return(W1);
}
else
{
seesaw=0;
// cout <<"seesaw = "<<seesaw<<" "<<W2<<" \n ";
return(W2);
}
}
/* ---- fin BBG.C -----------*/
Pour tester le générateur de nombre v.a. N (0, 1), on utilise le subroutine testBBG.C
/* ---- testBBG.C -------------- */
#include<iostream.h>
#include"BBG.h"

void testBBG(int seed, int n) {

BBG a(seed);
double sum =0.0;
double *X = new double[n];
for(int i = 0; i<n; ++i) {
 X[i] = a.gaussien();
sum += X[i];
}
sum = sum/n;
double var = 0.0;
for(int i = 0; i<n; ++i)
var += (X[i] - sum)*(X[i] - sum);
var = var/(n-1);
cout << " sum = " << sum/n << "\n";
cout << " var = " << var << "\n";
}
/* ----- Fin testBBG.C -------------- */

En C ou C++, tous les applications contiennent un seul main.C

/* ----- main.C -------------- */

#include<math.h> /* M_PI */
#include<stdlib.h>
/* #include<unistd.h> */ /* getpid */
#include<iostream.h>
#include"BBG.h"

extern void testBBG(int,int);

main()
{
int seed=getpid();
// int seed=100;
/* srand48(getpid()); */
int n = 100;
testBBG(seed,n);
}

Remarques

1. Dans les fichiers *.h, on peut uniquement mettre

Type de finition struct BBG {int seed; };
Template déclarations template<class T> class Z;
Fonction déclarations extern int strlen(const char*);
Inline fonction définition inline double dmax(double x, double y) {...};
Data déclaration extern int a;
constant définition const double pi = 3.14159;
Il faut jamais met-
Enumeration enum coulour {red, yellow, green};
Name déclaration class Matrix;
Include #include <std.h>
Macros définition #define VERSION 12
conditional compilation directives #ifdef _cpluscplus
comments /* tres important */
tre
Ordinary function definition double dmax (doube x, double y)...
Data definition int a;
2. Il faut connaı̂tre les options standartes de g++, par exemple -o -g -Wall -C -O -E lire attentivement le
manuelle man g++

3. Un bon Makefile est très important aussi, voir le fichier joint.

4. Comprendre les options -I -L

Voici le fichier Makefile

# ----- ZHANG Huilong Makefile for C++

#
# ------------ C++ compiler [g++] options on Unix (Sparc Linux HP)
#
# -g : produces additional symbol table information for debugging
# -c : suppresses loading and produces .o files
# -C : checks that subscripts are within declared array bounds
# -o : names the final output file
# -O : optimizes
# -S : stops after the stage of compilation proper
# -p : prepares object files for profiling (prof)
# -pg : prepares object files for profiling (gprof)
# -Wall : issues warnings for conditions

# ------------ directories

# ------------ libraries

LIBMATH = -lm

# ------------ include
INC = -I
# ------------ macros definition
SRC =
#
# ------------ object files
OBJS = main.o\
BBG.o \
testBBG.o \
# ------------ compilateur
COMPG++ = g++
COMPGCC = gcc
# ------------ options
RULE = -o -g $@ $?
CLEAN = @ /bin/rm $?
MESS = @ echo "compilation ... $@ done
# OPT = -g
OPT = -g -Wall
# OPT = -O
# ---------------
# OPT = -g compiler avec cette option vous permet de debuger
# a l’aide de xxgdb ou gdb
# OPT = -O programme compile avec cette option sera optimise, et
# tourne plus vite, mais on ne pourra pas debuger.
# ---------------

essai : $(OBJS)
$(COMPG++) $(OPT) $(OBJS) -o essai -lm
$(MESS)
clean:
@rm -f $(OBJS)
# ----- fin Makefile for C++
2.4 Recherche et estimation de la tendance
Le mot lissage : exprime le mouvement tendanciel d’une succession de points (t, Zt ) situés au milieu d’un nuage de
points (t, Xt ).
Si on désire exprimer la tendance avec une courbe très lisse, une méthode consiste à proposer une expression
analytique. on calcule souvent cette expression à partir de données lissées jugée convenable (non à partir des données
initiales). Alors c’est quoi une bonne série lissée ?
Il est difficile de définir. les motivations sont souvent très variées. Il ne s’agit pas le plus souvent d’atteindre le
cas limite : la droite, ni de s’ajuster au point de passer chaque point observé. Lisser le moins fort possible tout en
atteignant un objectif précis

· Se débarrasser des résidus.

· Se débarrasser des résidus et du saisonnier

de telle sorte que ces résidus (ou résidus + saisonnière) fluctuent autour de la courbe lisse.

2.4.1 Moyennes par intervalles

Pour simplifier la présentation, on considère le problème sans variations saisonnières

Xt = T t + εt

Le but consiste à mettre en évidence la tendance Tt . La méthodes par intervalle est utilisée pour résumer les séries
chronologiques comportant un grand nombre d’observations. La démarche consiste à approcher la chronique par une
fonction en escalier ou par une fonction affine par morceaux. Plus précisement, soient I1 , I2 , · · · , Ij , · · · intervalles de
temps (les semaines, les mois, les années, etc) on calcule les moyennes des observations àppartenant à ces intervalles
pour obtenir la série mj , où j est une indice centrée, pas toujours une indice coincidant avec t dans Xt . Cette méthode
très simple présente certains avantages :

· Mettre en évidence la tendance de la série

· La variance de la variation accidentelle est divisée par le nombre de terms considérés

· Les valeurs successives sont indépendantes si les εt le sont.

l’inconvenient de cette méthode est la diminuition importante du nombre de valeurs, elle est appliquée uniquement
dans des séries présentant un nombre important d’observations.

exemple 2.4.1 On considère la série suivante. La figure 2.4.1 compare sa simulation et la moyenne par intervalle,
la longueur de la moyenne est de 4.

Tt = 0.05t + 100

S1 = −1, S2 = 2, S3 = −2, S4 = 1

εt ∼ N (0, 0.52)

2.4.2 Moyennes mobiles

En pratique, le nombre d’observations est insuffisant pour que l’on puisse se satisfaire de moyennes par intervalle.
Deux cas se présentent

· le nombre q d’observations considérées est impaire

i=+k
1 X
mt = Xt+i , k vérifie q = 2k + 1
q
i=−k

Par exemple
1
mt = (Xt−2 + Xt−1 + Xt + Xt+1 + Xt+2 )
5
 107

106

105

104

103

102

101

100

99

98

97
0 10 20 30 40 50 60 70 80 90 100

Figure 2.5: Moyenne par intervalle

· Le nombre q d’observations considérées est paire, la formule est

" i=k−1 #
1 X (Xt−k + Xt+k )
mt = Xt+i + k vérifie q = 2k
q 2
i=−k+1

Par exemple, soit q = 4  

1 1
mt = Xt−1 + Xt + Xt+1 + (Xt−2 + Xt+2 )
4 2
On remarque facilement que dans l’exemple de Fuel consommation, chapitre 1, table 1.2, on a utilisé cette
moyenne mobile avec q = 4.

Remarque 2.4.2 On note que

1. Il s’agit d’une moyenne pondérée.
2. Pour une suite d’observations (t = 1, · · · , T ), les instants où on peut calculer mt sont (k + 1, k + 2, · · · , T − k),
soit seulement t − 2k instants.
3. Non independance entre les mt . Par exemple mt et mt+1 ne sont pas indépendants, car

Xt−k , · · · Xt+k =⇒ mt

Xt−k+1 ··· Xt+k+1 =⇒ mt+1

Xt−k+1 , · · · , Xt+k sont prises en compte à la fois dans le calcul de mt et dans le calcul de mt+1 .

Nous étudions dans la suite de cette section l’effet du filtre de moyenne mobile sur la tendance d’une chronique.

Chronique de modèle Xt = a + εt
Il s’agit d’un modèle d’ajustement simplifié : constant + perturbation aléatoire. On suppose E εt = 0, Var εt = σ 2 et
les {εt } sont indépendantes identiques distribuées (Bruit Blanc). Soit q = 2k + 1 on a alors ∀h ≥ 1
i=k
1 X
mt = εt+i + a
q
i=−k

i=k
1 X
mt+h = εt+h+i + a
q
i=−k
 σ2
Proposition 2.4.3 La Série mt ne possède pas de tendance et elle est de variance q

σ2
Preuve: E mt = a, Var mt = q

Remarque 2.4.4 On remarque que

1. On dit que la serie ne possède pas de composante tendancielle (où que la tendance est nulle) dans la mésure où
sa moyenne est constante.

2. La variance est divisée par q par rapport à la série initiale, les mt sont soumises à une variation accidentelle
d’amplitude beaucoup plus petite (lisse).

Chronique de modèle Xt = at + b + εt
Ce modèle fait partie aussi des chroniques à tendance polynomiale. On suppose que la longueur de moyenne mobile q
est impaire
i=k i=k
1 X 1 X
mt = [a(t + i) + b + εt+i ] = [at + ai + b + εt+i ]
q q
i=−k i=−k

"
i=k 0 i=k
#
1 X X X
= at + b + a i+a i+ εt+i
q i=0 i=−k i=−k

i=k
1 X
= at + b + εt+i
q
i=−k

Proposition 2.4.5 Si Xt est une chronique dont la composante tendantielle est linéaire, alors la série mt des
moyennes mobiles concervent la même tendance et a pour composante accidentelle la moyenne mobile des composantes
accidentelle de la série initiale.

On admet (en exercice) le même résultat pour une moyenne mobile de longueur q = paire.

Chronique de modèle Xt = at2 + bt + c + εt

Nous pouvons nous limiter à l’étude du terme de seconde dégré. (Pourquoi ? A réflechir...). Soit∗ q = 2k + 1, impaire
i=k i=k
1 X   1 X  2 
mt = a(t + i)2 = a(t + 2it + i2 )
q q
i=−k i=−k

" k k
# k
2 a X X a X 2
= at + 2t i+ i = at2 + 2
2
i
q q i=1
i=−k i=−k

k(k + 1)
= at2 + a
3
Proposition 2.4.6 Soit Xt une chronique dont la composante tendancielle est un polynôme de 2ième degré at2 +bt+c,
alors la composante tendancielle de la série des moyennes mobiles est égale à

ak(k + 1)
Tt = at2 + bt + c +
3
Remarque 2.4.7 Il est évident que

1. l’estimateur de la tendance par le filtre des moyennes mobiles est biaisé, le biais est constant.

2. On montre que dans le cas d’une tendance du 3ième degré, c’est un biais linéaire en t et dans le cas d’une
tendance du 4ieme degré, c’est un biais en t2 , etc.

Il faut donc se montrer prudent dans l’utilisation des moyennes mobiles et ramener le problème autant que possible à
une fonction linéaire en t. La figure 2.6 illustre un exemple d’une série à tendance polynomiale de 2ième degré et ces
moyennes mobiles.
k
X k(k + 1)(2k + 1)
∗ i2 =
1
6
 500

450

400

350

300

250

200

150

100

50

0
0 5 10 15 20 25 30

Figure 2.6: Moyenne mobile de longueur 12 pour une chronique de degré 2

Présence de variations saisonnières

Soit
X t = T t + St + ε t
avec St+p = St ∀t où p est la période de St . Nous nous limitons à l’étude de la série St , car nous savons que la
moyenne mobile d’une somme est la somme des moyennes mobiles. On considère une moyenne mobile de longueur
q = 2k + 1.
k
1 X
mt = St+i
q
i=−k

Remarque 2.4.8 On remarque que

1. Il est donc évident que si la longueur de la moyenne mobile q et la période de variation saisonnière sont égales
(p = q 2k + 1 = p) alors la moyenne mobile calculée sur St est constante.
2. Il est évident également que la moyenne mobile est constante si sa longueur est un multiple de la période q = rp

Proposition 2.4.9 Dans le cadre du modèle additif, le filtre de la moyenne mobile de longueur q élimine les variations
saisonnière dont la période p divise q.

En figure 2.7, nous donnons la représentation simultanée d’une série et sa moyenne mobile de longueur 4.

Quelques remarques supplémentaires sur les moyennes mobiles

Il existe d’autres moyennes mobiles plus génerales, non symétriques, elle sont utilisées pour la prévision à courte terme
pour conserver des tendances ou pour éliminer des variations saisonnières particulières. Par exemple, les moyennes de
Spencer
· Moyennes de Spencer de longueur 15
1
(−3, −6, −5, 3, 21, 46, 67, 74, 67, 46, 21, 3, −5, −6, −3)
320
Réduction de variance résiduelle est de 0.19.
· Moyennes de Spencer de longueur 21
1
(−1, −3, −5, −5, −2, 6, 18, 33, 47, 57, 60, 57, · · ·, −1)
350
Réduction de variance résiduelle est de 0.15.
 15

14

13

12

11

10

8
0 5 10 15 20 25 30 35 40

Figure 2.7: Moyenne mobile de longueur 4 pour une chronique de période 4

Les moyennes Spencer annulent les saisonnalités de période 4 en concervant les tendances polynomiales de degré ≤ 3
(sans biais).
Exercice 2.4.10 Soit la moyenne
1 1
mt = (Xt−1 + Xt + Xt+1 + (Xt−2 + Xt+2 ))
4 2
Montrer que cette moyenne mobile élimine les variations saisonnières de période 4.
Solution q = 4 paire. On sait que St−2 = S1 , St−1 = S2 , St = S3 , St+1 = S4 , St+2 = S5 = S1 .
 
1 1
St−1 + St + St+1 + (St−2 + St+2 )
4 2
 
1 1 1
= S2 + S3 + S4 + (S1 + S1 ) = (S1 + S2 + S3 + S4 ) = 0
4 2 4

Exercice 2.4.11 Soit St est évolutive linéaire de période 4. Montrer que la moyenne mobile pondérée
1
mt = (Xt−3 + 2Xt−2 + 3Xt−1 + 4Xt + 3Xt+1 + 2Xt+2 + Xt+3 )
16
élimine la variation saisonnière
Solution Soit S1 , S2 , S3 , S4 coefficients saisonniers. On pose

St−3 = S1 , St−2 = S2 , St−1 = S3 , St = S4

St+1 = S1 , St+2 = S2 , St+3 = S3

On a alors
1
(St−3 + 2St−2 + 3St−1 + 4St + 3St+1 + 2St+2 + St+3 )
16
1 1
= (S1 + 2S2 + 3S3 + 4S4 + 3S1 + 2S2 + S3 ) = (4S1 + 4S2 + 4S3 + 4S4 ) = 0
16 16

2.4.3 Différences successives

Les différence successives sont utilisées pour étudier la tendance et pour obtenir une chronique stationnaire (la notion
stationnaire sera étudiée dans les chapitres suivants).
Définition
Les différences successives d’une chronique sont obtenues par transformation linéaire des observations initiales.

DXt = Xt+1 − Xt , D2 Xt = DXt+1 − DXt = D(DXt )

d’une façon générale

Dk Xt = Dk−1 Xt+1 − Dk−1 Xt = D(Dk−1 Xt )
Donc la série Dk Xt est définie de t = 1 à t = T − k. On note que le filtre des différences successives conserve la
périodicité des variations saisonnière éventuelles, car St+1 − St est périodique et de période p. On remarque également.

DXt = Xt+1 − Xt = Tt+1 − Tt + St+1 − St + εt+1 − εt

C’est une série de modèle additif. On étudie dans la suite l’effet des différences successives sur les modèles de chronique

Effet sur bruit blanc

Soit Xt = εt alors DXt = εt+1 − εt , donc le filtre linéaires des différences premières n’est pas un lissage puisque
la variance de DXt est multipliée par 2. On obtient une suite de v.a. non indépendantes (DXt et DXt+1 sont
dépendantes). En figure 2.8 illustre les différences successives d’un bruit blanc. On remarque l’amplitude des variations
augmente avec l’ordre des différences successives.

6
"BB"
"DX"
"D2X"
4

-2

-4

-6

-8
0 5 10 15 20 25 30 35 40 45 50

Figure 2.8: Représentation simultanée d’un BBG, DXt et D2 Xt

Proposition 2.4.12 La série des différences premières d’une chronique de modèle Xt = a + εt où εt est un bruit
blanç est une suite de v. a. de moyenne nulle, de variance 2σ 2 .

Effet sur le modèle polynomial

On considère une chronique de modèle Xt = P (t) + εt où P (t) est un polynôme en t de degré n. εt un bruit blanc de
variance σ 2 .

· Soit Xt = at + b + εt alors DXt = Xt+1 − Xt = a + εt+1 − εt On a donc

E DXt = a, Var DXt = 2σ 2

c’est à dire DXt est une suite de v. a. de moyennes a et de variance 2σ 2 .

· Soit Xt = at2 + bt + c + εt alors

DXt = 2at + a + b + εt+1 − εt
 chronique composée d’une tendance linéaire et d’une variation accidentelle.
D 2 Xt = 2a + (εt+2 − εt+1 ) − (εt+1 − εt )
= 2a + εt+2 − 2εt+1 + εt
chronique composée d’une tendance constante et d’une variation accidentelle. On calcule facilement que
E D2 Xt = 2a Var D2 Xt = 6σ 2

La figure 2.9 donne la représentation d’une chronique simulée dont la tendance est un polynôme de degré 2, le DXt et
D2 Xt . Les différences premières presentent une tendance linéaire et les différences secondes une tendance constante,
on constate également que l’amplititude des variations augmente avec l’ordre des différences.

0 10 20 30 40 50 60 70 80 90 100

Figure 2.9: Représentation simultanée d’une chronique (échelles différentes), d’une chronique de terme Xt = at2 +bt+c+εt ,
DXt et D2 Xt

Proposition 2.4.13 Soit Xt une chronique dont le terme général de la tendance est défini par un polynôme de degré
n. Alors la série des différences Dn Xt est de la forme a + εt , où εt est une suite de v. a. de moyenne nulle et de
variance constante par rapport à t.

2.4.4 Autocorrélation et correlogramme

Nous avons vu précedemment que les moyennes mobiles d’une chronique sont corrélées.

Définition et propriété
Définition 2.4.14 On appelle coefficient d’autocorrélation de rang k de la chronique Xt , t = 1, · · · T le coefficient † :
T
X −k
(Xt − X̄)(Xt+k − X̄) T
t=1 1X
rk = X̄ = Xt
T
X T t=1
(Xt − X̄)2
t=1
† rappel: Soient ξ, η deux v. a. alors le coefficient de corrélation est défini par
E (ξ − E ξ)(η − E η)
r= √ √
Var ξ Var η
ce coefficient peut être estimé par Pn
i=1 (ξi− ξ̄)(ηi − η̄)
r̄ = qP
n
pPn
2 2
i=1 (ξi − ξ̄) i=1 (ηi − η̄)
où ξi et ηi sont deux échantillons de taille n. Il y a donc une légère différence avec la définition 2.4.14, ce choix est dicté par le fait que εt
est de moyenne µ et de variance σ2 , indépendantes de t.
Définition 2.4.15 On appelle corrélogramme la représentation graphique de l’application: k −→ rk .
Le coefficient rk est défini pour mesurer la dépendance entre Xt et Xt+k .
Théorème 2.4.16 On considère un bruit blanc εt . les coefficient d’autocorrélation rk de rang supérieur ou égale à 1
suivent approximativement une loi normale de moyenne − T1 et de variance T1 .
ce théorème permet de juger de la taille des coefficients d’autocorrélation: en négligeant la moyenne − T1 , on peut dire
que la probablilité qu’un coefficient soit supérieur en valeur absolut à 1.96
√
T
est égale à 0.05.
1.96
) = 0.05
P (|rk | ≥
T
Pour illustrer le théorème, nous avont simulé un bruit blanc de 500 termes dont les coefficients d’autocorrélation sont
présentés dans la table 2.1. (Voir le corrélogramme Fig. 2.10.)

r1 = −0.026720 r2 = 0.052195 r3 = −0.069710 r4 = −0.097432

r5 = 0.042616 r6 = −0.021875 r7 = 0.044295 r8 = −0.023381
r9 = 0.043714 r10 = −0.020600 r11 = −0.018972 r12 = −0.000013
r13 = −0.006135 r14 = 0.022235 r15 = 0.005640 r16 = −0.047127
r17 = 0.013550 r18 = −0.053561 r19 = 0.005232 r20 = 0.034764

Table 2.1: coefficient d’autocorrélation k = 1, · · · , 20 sur une série BBG

0.8

0.6

0.4

0.2

-0.2
0 5 10 15 20

Figure 2.10: corrélogramme d’un bruit blanc simulé

Analyse des variation accidentelle obtenues par l’ajustement d’un modèle à des données observées est fondamentale
pusiqu’elle permet de vérifier l’adéquation du modèle choisi. Cette vérification est faite en deux étapes:
1. Calcul des résidus, définis par la différence entre les observations et les valeurs estimées, forment un bruit blanc
(si le modèle est correct), à partir duquel on peut plus tirer d’information.
2. Calcul du corrélogramme, et vérifier l’hypothèse selons la quelle il s’agit bien un bruit blane avec le théoreme
2.4.16

Corrélogramme de chronique de tendance non nulle

Le corrélogramme de chronique de tendance non nulle est en général caractéristique : les coefficients d’autoccorrélation
restent relativement élevés jusqu’à un rang important si la série est monotone. Lorque la tendance change de sens au
cours de la période étudiée, certains coefficients d’autocorrélation prennent des valeurs négatives, ceci correspond aux
changement de régime. (ruptur de comportement). On peut simuler un exemple en prenant deux tendances différentes
(monotone croissant suivi par monotone décroissant figure 2.11) en mileux de la série et tracer le corrélogramme associé,
figure 2.12.
 30

25

20

15

10

-5

-10

-15
0 10 20 30 40 50 60 70 80 90 100

Figure 2.11: changement de régime : trajectoir

Autocorrélation et filtre des moyennes mobiles

Le filtre des moyennes mobiles transforme la chronique en une serie fortement corrélée. Considérons les moyenne
mobiles de longueur impaire (q = 2k + 1) mt et mt+h , appliquées sur un bruit blanc ,
k k
1 X 1 X
mt = εt+i mt+h = εt+h+j
q q
i=−k j=−k

On étudie maintenant les coefficients de corrélation théoriques entre mt et mt+h .

k k
1 X X
Cov (mt , mt+h ) = ( )2 E (εt+i , εt+h+j )
q
i=−k j=−k

Ces espérance sont tous nulles sauf si les indices sont les même
t+k
1 X 2k − h + 1 2
Cov (mt , mt+h ) = ( )2 σ2 = σ
q (2k + 1)2
i=t+h−k

où q = 2k + 1 et h ≤ 2k. Donc le coefficient théorique est

Cov (mt , mt+h ) 2k − h + 1 h
rh = √ p = =1−
Var mt Var mt+h 2k + 1 2k + 1
Ce coefficient est d’autant plus grand que la longueur de la moyenne mobile est grande, et d’autant plus petit que le
rang h est grand. pour h > 2k, il est nul. La figure 2.13 est une simulation. On notera que le coefficient de corrélation
est positive ou nul, cela signifie que les séries mt et mt+h varient dans le même sens.

Autocorrélation et filtre des différences successives

Nous avons vu que la présence d’une tendance croissante ou décroissante dans la série étudiée se traduit par de fortes
valeurs des coefficients d’autocorrélation est en lui même intéressant car on cherche souvent à éliminer la tendance de
la chronique pour n’en conserver que la part aléatoire. On sait que les différences successives d’ordre n transforment
les chroniques à tendance polynomiale de degré n en chronique à tendance constante. Au niveau des coefficients
d’autocorrélation, l’effet n’est pas négligeable. On étudie dans la suite

Autocorrélation des différences premières d’un modèle linéaire

Soit Xt = at + b + εt . On calcule facilement
DXt = a + εt+1 − εt
Cov (DXt+1 , DXt ) = E (εt+2 − εt+1 )(εt+1 − εt ) = −E(εt+1 )2
 1

0.8

0.6

0.4

0.2

-0.2

-0.4

-0.6
0 20 40 60 80 100

Figure 2.12: changement de régime : corrélogramme

Car εt est un bruit blanc. On ontient ainsi le coefficient autocorrélation d’ordre 1:

δ1 = √ Cov (DXt+1
√
,DXt )
Var DXt+1 Var DXt
σ2 1
= 2σ2 = − 2

Proposition 2.4.17 le coefficient d’autocorrélation de rang 1 des différences premières d’une chronique à tendance
linéaire est égal à −0.5 et les coefficients d’autocorrélation de rang supérieur sont nuls.

Autocorrélation des différences successives d’une chronique à tendance polynomiale

Etudions le cas où Xt = at2 + bt + c + εt . On sait que les différences secondes sont

D2 Xt = 2a + εt+2 − 2εt+1 + εt = 2a + Et

Les Et+h et Et sont indépendantes pour h ≥ 3, pour d’autre h on a

E (D2 Xt+1 − 2a)(D2 Xt − 2a) = E (Et+1 Et )

= E (εt+3 − 2εt+2 + εt+1 )(εt+2 − 2εt+1 + εt )
= −2E ε2t+2 − 2E ε2t+1 = −4σ 2

On en déduit que le coefficient de corrélation σ1 = − 32 . calculons maintenant σ2

E (D2 Xt+2 − 2a)(D2 Xt − 2a) = E (Et+2 Et )

= E (εt+4 − 2εt+3 + εt+2 )(εt+2 − 2εt+1 + εt )
= E ε2t+2 = σ 2
1
donc δ2 = 6

Proposition 2.4.18 Le coefficent autocorrélation de rang 1 des différences secondes d’une chronique à tendance
polynômiale de degré 2 égal à − 23 et le coefficent d’autocorrélation de rang 2 à 61 . Les coefficents d’autocorrélation de
rang supérieur sont nuls.
Voir en figure 2.14 On admet plus géneéralement
Proposition 2.4.19 Les coefficients d’autocorrélation de rang supérieur ou égal à n + 1 des différences successives
d’ordre n sont nuls pour une chronique de tendance polynômiale de degré n, sans saisonnalité.
Une façon d’éliminer la tendance d’une chronique est donc de calculer les différences successives jusquà l’obtention
de coefficients d’autocorrélation dimimuant rapidement. On demande comme exercice la simulation et le calcul d’un
exemple de Xt = at2 + bt + c + εt afin de vérifier cette proposition (figure 2.15). Les coefficients d’autocorrélation
 1

0.8

0.6

0.4

0.2

-0.2
0 5 10 15 20

Figure 2.13: corrélogramme des moyennes mobiles de longueur 7 d’un bruit blanc simulé 600 termes

calculés sur la série initiale sont tous presque égaux à 1, ce qui dénote une tendance monotone. Les coefficients
d’autocorrélation reste élevés sur la série des différences premières (compris entre 0.5 et 1), ce qui montre l’existence
d’une tendance non nulle. Par contre, les coeffcient d’autocorrélation sur les différences seconde mettent en évidence
l’absence de tendance. Il a donc fallu aller jusqu’aux différence du seconde ordre pour annuler la tendance existante
dans la série initiale.

2.4.5 Méthode des moindres carrés

Une autre méthode pour calculer la tendance est la méthode des moindres carrés, en régressant les observations par
le temps.

Tendance polynômiale
Nous avons (dans la section 2.4.2) étudié le cas de tendance polynômiale à l’aide de moyennes mobiles, il paraı̂t
normal, lorsque la forme de la tendance est fixée, d’utiliser le critère des moindres carrés pour estimer les coefficients
du polynôme supposé. Nous nous limitons ici au cas où la série ne subit pas de variations saisonnières.
Le critère des moindres carrés appliqué à un polynôme de degré 2 est donné par la formule suivante
T
X
S= {Xt − (at2 + bt + c)}2
t=1

La méthode consiste à déterminer les coefficients a, b et c tels que la quantité S soit minimale. Cette somme représente
la somme des carrés des écarts entre les observations Xt et les valeurs du polynôme. On sait que le minimum est
atteint pour les valeurs des coefficients qui annulent les dérivées partielles premières:
∂S X
= − 2{Xt − (at2 + bt + c)}t2 = 0
∂a
∂S X
= − 2{Xt − (at2 + bt + c)}t = 0
∂b
∂S X
= − 2{Xt − (at2 + bt + c)} = 0
∂c
On obtient donc le sytème de trois équations à trois inconnues suivant:
X X X X
a t4 + b t3 + c t2 = X t t2
X X X X
a t3 + b t2 + c t = Xt t
X X X
2
a t +b t + cT = Xt
 1

0.8

0.6

0.4

0.2

0.0

−0.2

−0.4

−0.6

4 8 12

Figure 2.14: corrélogramme des différences premières d’une chronique à tendance linéaire.

Ces trois équations, appelées équations normales, se limitent lorsque le polynôme est de degré 1 (tendance linéaire),
aux deux équations suivantes: P 2 P P
P b t +P c t = Xt t
b t + cT = Xt

Donc on en déduit facilement les formules suivantes‡

Cov (t, Xt )
b= c = E Xt − bE t
Var t
Rappelons, pour faciliter es calculs:
PT T (T +1)
t =
PTt=1 2 2
T (T +1)(2T +1
t =
PTt=1 3 6
t = ( T (T2+1) )2
PTt=1 4 T (T +1)(2T +1)(3T 2 +3T −1)
t=1 t = 30

Toutes les tendances ne sont évidemment pas linéaires ou polynomiales, la présence d’une limite vers laquelle la série
converge rend impossibble l’ajustement de la tendance par un polynôme. Nous présentons dans la suite que des
modèles les plus courants.
‡ Ici, on prolonge la définition d’espérance (variance) mathématique ainsi que son estimateur empirique à des fonctions déterministes:
T
1 X
E f (t) ∼ f (t)
T t=1
T T
1 X 1 X
Var f (t) ∼ (f (t) − f (t))2 = E f (t)2 − (E f (t))2
T t=1 T t=1
 1

0.8

0.6

0.4

0.2

0.0

−0.2

−0.4

−0.6

4 8 12

Figure 2.15: corrélogramme des différences successives d’une chronique à tendance polynômiale de degré 2

Tendance exponentielle

La tendance exponentielle, souvent constatée dans les série d’observation économiques, s’exprime de la façcon suivante
:
Tt = cebt ou log(Tt ) = bt + log(c)

On suppose c > 0, l’ajustement est réalisé en considérant le logarithme des observations et de la tendance : on cherche
les valeurs b et c telles que la somme S soit minimale :

T
X
S= {(log(Xt ) − (bt + log(c)))}2
t=1

On se ramène donc au cas l’ajustement linéaire, bien que le critère des moindres carrés, utilisé sur la série des
logarithmes, ne soit pas équivalent au critère des moindres carrés appliqué aux observations.

Modèle exponentiel à saturation

La tendance est sous la forme

Tt = a − bekt

où a est un paramètre connu (limite de la tendance), b est positif et k négatif. La tendance est croissante vers a pour
les grandes valeurs de t et l’accroissement marginal Tt est décroissant. L’ajustement est réalisé en se ramenant au cas
précedent par le changement de variable suivant:

Tt = Tt − Tt−1 = b(1 − e−k )ekt

Modèle logistique
Le modèle logistique est défini de la façon suivante :
1
Tt =
a + bekt
où b > a > 0 et k < 0. Cette fonction est encore croissante mais l’accroissement marginal n’est pas monotone :
il a tendance à augmenter jusqu’à la valeur t = log(a/b)
k , qui est mathématiquement l’abscisse du point d’inflexions,
puis à diminuer. Ce modèle est particulièrement utilisé pour ajuster les taux d’équipement des ménages en biens de
consommation durable. Pour effectuer l’ajustement, on se ramène au modèle précédent en posant :
1
Tt = = a + bekt
Tt

Voir la figure 2.16§

80

70

60

50

40

30

20

10

0
10 20 30 40 50 60 70 80 90 100

Figure 2.16: Tendance exponentielle, exponentielle à saturation, tendance logistique

§ gnuplot plot 0.5 exp(0.05x), 77 − 77 exp(−0.04x), 1.0/((1.0/75) + 0.3 exp(−0.1x))

2.5 Estimation des variations saisonnières
Après la recherche du modèle et l’estimation de la tendance, nous étudions dans ce chapitre l’estimation des variations
saisonnières

2.5.1 Estimation des coefficient saisonniers

Il s’agit maintenant d’estimer les variations saisonnières en supposant connue la tendance. Nous utilisons un tableau
(Table de Buys Ballot 2.2). où p est la période, et n est le nombre de période considérée (bien sûr ici, on considère

1 2 ··· j ··· p
1 X1,1 X1,2 ··· X1,j ··· X1,p X 1·
2 X2,1 X2,2 ··· X2,j ··· X2,p X 2·
··· ··· ··· ··· ··· ··· ··· ···
i Xi,1 Xi,2 ··· Xi,j ··· Xi,p X i·
··· ··· ··· ··· ··· ··· ··· ···
n Xn,1 Xn,2 ··· Xn,j ··· Xn,p X n·
X ·1 X ·2 ··· X ·j ··· X ·p X

Table 2.2: Table de Buys Ballot

que p est connu). On a la relation suivante

t = (i − 1)p + j 0<j≤p i≥1

D’après le tableau on a aussi

p n
1X 1X
X i· = Xi,j X ·i = Xi,j
p j=1 n i=1
et
n p
1X 1X
X= Xi· = X·j
n i=1 p j=1

Alors le modèle
X t = T t + St + ε t
devient
Xi,j = Ti,j + Sj + εi,j

Estimation de coefficients saisonnièrs

On déduit de l’expression du modèle

∀j = 1, · · · , p, Sj = Xi,j − Ti,j − εi,j ∀i = 1, · · · n

Il est naturel de choisir pour estimateur Sj la moyenne

n
1X
sj = (Xi,j − Ti,j )
n i=1

Quelques modifications sont utilsées

· Moyenne en éliminant de la somme ci-dessus les valeurs extrêmes des différences.
· Principe de la conservation des aires: La somme des coefficients saisonnièrs sur une période doit être nulles, pour
cela, il sufffit de considérer les valeurs sj − s où s est la moyenne de sj sur j.
p
1X
s= sj
p j=1

A noter aussi que le filtre moyennes mobiles permet d’estimer la tendance sauf pour les premièrs et les derniers termes,
l’estimation des coefficients ne peut être effectuée que sur n − 1 périodes dans le tableau Buys Ballot.
Les coefficients saisonniers sont utlisés pour désaisonnalisér les séries. L’objectif est d’éliminer des observations
les variations saisonnières afin de pouvoir les comparer entre elles, et comparer la dernière observation connue à la
précédente. On peut ainsi suivre l’évolution mois par mois du nombre de chomeurs, de l’indice des prix sans que la
variation sainsonnière soit prise en considération.
Définition 2.5.1 Soit sj l’estimateur des coefficients saisonniers, on appelle série désaisonnalisé (ou série corrigée
des variations saisonnières) la série des observations après élimination des variations saisonnières
′
Xi,j = Xi,j − sj
′
On obtient ainsi une série Xt .
Example 2.5.2 (TP). On se donne dans un fichier obs une chronique Xt , t = 1, · · · , 100 de période p = 4. Calculer
les coefficients saisonnièrs. On en utlise pour cela la moyenne mobiles de longueur q = 4, pour estimer la tendance,
et moyenne centrée pour estimer les coefficients saisonniers en éliminant la plus grande et la plus petite valeur dans
chaque estimation.

2.5.2 Recherche simultanée de la tendance et des variations saisonnières

Dans les paragraphes précédents, nous avons proposé des méthodes d’estimation des variations saisonnières indépendantes
de la procédure choisie pour ajuster la tendance. Il est clair qu’il est plus intéressant de lier les deux démarches

Filtre de Buys Ballot

Le filtre de Buys Balot consiste à ajuster, par le critère des moindres carrés, à la série Xi,j un modèle additif intégrant
une tendance linéaire et des variations saisonnières constantes de période connue p. Il s’applique evidemment aux
séries qui suivent le premier modèle multiplicatif. Soit

Xt = at + b + St + εt

en adoptant les notations de Buys Ballot,

Xi,j = a(p(i − 1) + j) + b + Sj + εi,j

Les coefficients Sj (j = 1, · · · p) et les paramètres a et b de la tendance sont déterminés par mimnimisation de la somme
des carrés des résidu. XX XX
Q= e2i,j = (Xi,j − a(p(i − 1) + j) − bj )2
i j i j

△
où bj = b + Sj . Nous admettrons que cette somme est minimale
P lorsque les dérivées partielles en bj et a sont nulles. En
ajoutant en plus le principe de conservation des aires ( Sj = 0), on déduit pour les paramètres les valeurs suivante
" n #
12 X n(n + 1)
a = (iX i· ) − X
np(n2 − 1) i=1 2

a(np + 1)
b = X−
2
p+1
sj = X ·j − X − a(j − )
2
On montre en fait que cet ajustement revient à déterminer la droite des moindres carrés sur la tendance définie
par les moyennes annuelles affetées aux fins de périodes

(X i··· , ip), i = 1, · · · , n

et à estimer les coefficients saisonniers en écrivant que

∀j = 1, · · · , p sj = X ·j − aj − cte
P
La constante étant déterminée par le principe de conservation des aires: sj = 0.

Filtre de Buys Ballot

La méthode Census est employée surtout en économie pour estimer simultanément la tendance et les variations
saisonnières d’une chronique. Elle consiste à faire passer le filtre des moyennes mobiles, à estimer les coefficients
saisonniers, à désaisonnaliser la série avec ces premières estimations, et à lisser à nouveau la série désaisonnalisée pour
réestimer la tendance, ce qui permet enfin de réestimer définitivement les coefficients saisonniers.
La difficulté réside dans le choix des longueurs des moyennes mobiles. La première doit être évidemment de
longueur égale à la période des variations saisonnières, et la longueur de la seconde est choisie après examen de la série
désaisonnalisée.
2.5.3 Autocorrélation et variations saisonnières
On étudie l’influence de variations saisonnière sur les coefficients d’autocorrélation. Le but de ces études : mettre en
évidence la période p.

Théorème 2.5.3 Soit Xt une série chronique périodique de période p et de modèle Xt = Sr + εt , avec t = pq + r,
r > 0, r ≤ p. Alors les coefficients d’autocorrélations rk sont approximativement périodique de période p.

différence de rang p
On appelle différence de rang p de la série Xt les différence de deux observations séparées par p terme Xt+p − Xt .
On remarque que seul les différence premières sont égales à la différence de rang 1. Etudions l’effet de rang p sur une
chronique sans tendance de période p:

Xt = Sr + εt t = pq + r 0 < r ≤ p εt bruit blanc

Soit
Yt = Xt+p − Xt = εt+p − εt

Yt+1 = Xt+p+1 − Xt+1 = εt+p+1 − εt+1

La loi de Xt+p − Xt est donc indépendante de t; Sous l’hypothèse que la série εt est un bruit blanc , il est évident que
le coefficient d’autocorrélation de rang 1 est nul et qu’il en est de même pour tous ceux dont le rang est différent de
p puisque les termes Yt et Yt′ n’ont alors aucun point commun. Les coefficients d’autocorrélation de rang p sont par
contre differents de 0¶ :
Yt+p = Xt+2p − Xt+p = εt+2p − εt+p
On en déduit facilement rp = − 21
Théorème 2.5.4 Soit Xt une chronique sans tendance et périodique de période p, alors la série des différence de
rang p ne présente pas de variation saisonnière et tous les coefficients d’autocorrélation théorique sont nuls sauf le rp
est égale à − 12 .
Ce résultat, facile à montrer, est illustré dans la figure 2.17. La présence d’une tendance dans la chronique peut être
envisagée : si elle est linéaire, le résultat est identique puisque cela revient à ajouter une constante aux différences Yt .
Si elle est de degre 2, elle n’est pas éliminée par le filtre des différences de rang p, qu’il faut éventuellement compléter
par le filtre des différences premières.

0.8

0.6

0.4

0.2

0.0

-0.2

-0.4

-0.6

4 8 12

Figure 2.17: coefficients d’autocorrélation de Yt = Xt+4 − Xt pour une série de sans tendance de période p = 4

¶Y pour t′ > t, et t′ 6= t + p
t ⊥Yt′
Chapter 3

Modèles auto-régressifs et moyennes

mobiles

3.1 Introduction aux processus stationnaires

On sait qu’une variable aléatoire réelle (v.a.r) est une application d’un espace probabilité (Ω, A, P ) dans l’ensemble
de nombre réels (R, BR )
Définition 3.1.1 Processus aléatoires discrets : une famille de v.a.r. {Xt }t ∈ Z définies sur un même espace
probabilisé (Ω, A, P ) et à valeurs dans (R, BR )
L’indice t qui désigne le temps Z est un ensemble discret, souvent Z est N , l’ensemble des entiers.
Définition 3.1.2 Processus du 2 ordre : Si E Xt et Var Xt existent pour ∀t ∈ Z, alors on appelle {Xt } un
processus du 2 ordre
On ne considère que des processus du 2 ordre dans ce cours
Remarque 3.1.3 On remarque
· Dans le cas d’une chronique, on ne connait, pour un nombre fini de valeurs de t, qu’une seule réalisation des
v.a.r. Xt , on ne peut pas utiliser les résultats des probabilité ni ceux de la statistique inférentielle qui sont établis
à partir de la répétition de l’épreuve, on ne peut pas se payer le luxe de répéter ces expériences. La chronique
est appelée trajectoire du processus.
· La connaissance d’une seule réalisation d’une v.a. ne permet pas d’en déduire des propriétés intéressantes.
L’hypothèse de stationnarité, que nous définissons ci-dessous, permet de pallier cette difficulté, elle correspond
à la modélisation suivante : chaque observation d’une série chronique est considérée comme la somme d’une
valeur réelle déterminée par une fonction du temps et d’une valeur réelle résultant du hasard. Ces dernières
sont supposées être des réalisations d’une même v.a.r. X dont les paramètres ne dépendent pas du temps. Ce
qui permet d’appliquer les résultats de la statisitique, puisque on dispose alors d’un échantillon {Xt }t=1,··· de la
v.a.r. cet échantillon n’est pas nécessairement constitué des v.a.r. indépendantes

Définition 3.1.4 Processus stationnaire Un processus noté {Xt } est stationnaire si la loi du k−uple (Xt1 , Xt2 , · · · , Xtk )
est la même que celle du k−uple (Xt1 +r , Xt2 +r , · · · , Xtk +r ) pour toutes valeurs de t1 , t2 , · · · , tk , k, r

C’est une propriété très stricte, qui implique en particulier qu’une même loi régit chacun des variables {Xt }, t ∈ I,
et la loi de (X1 , X2 , · · · , Xn ) est la même que celle de (X1+k , X2+k , · · · , Xn+k ). Cette dernière implique qu’observer
un échantillon de mesure (X1 , X2 , · · · , Xn ) au temps 1, 2, · · · , n, plutôt que faire plus tard cette observation de n
mesures successives n’est pas important, que les lois mutuelles déterminées à partir d’observation du passé pourront
être légitimement utilisées dans le future puisque qu’elles sont invariantes pour une translation du temps. En pratique,
cette propriété de stationnarité sera supposée sans pouvoir être vérifiable ou testée.
Définition 3.1.5 Processus stationnaire du second ordre On appelle processus du 2 ordre un processus tel que
la loi de la v.a. Xt est indépendante de t, et que la covariance de Xt et Xt+k ne dépend que de k pour ∀t et k.
d’une façon équivalente, Xt est stationnaire du 2 ordre si la loi de Xt est indépendante de t et ∀t, s, (s, t) ∈ I, la
quantité Cov (Xt , Xs ) ne dépend que de l’écart t − s et non de s.

Remarque 3.1.6 En pratique

· Une série observée sera dite stationnaire quand un processus stationnaire modélise convenablement cette série.
 · Compte tenu des propriétés de la loi normale, un processus gaussien stationnaire du second ordre est stationnaire.
Nous nous limitons dans la suite de ce chapitre aux processus stationnaires du 2 ordre que nous appellerons
simplement processus stationnaires.
· On aura intérêt à avoir plutôt une série observée modélisable par un processus gaussien. Sinon, on envisage à
le transformer en un processus qui s’avérerait modélisable par un processus gaussien.

Définition 3.1.7 processus centré : Un processus {Xt } est centré si E Xt = 0, ∀t

Dans un premier temps, nous interessons aux processus centrés et stationnaires, ils vérifient aussi la propriété dite
d’homoscédasticité :
Var Xt = E Xt2 = σX 2
∀t
de plus pour tout k,
Cov (Xt , Xt−k ) = Cov (Xt , Xt+k ) = E Xt Xt−k = ck , ∀t
ck
et ρk = 2
σX

Définition 3.1.8 Processus ergodique Un processus est dit ergodique quand les limites des moments empiriques
convergent en moyenne quadratique vers les espérance correspondantes, lorsque la taille de la série s’accroı̂t.
Par example, la moyenne d’un processus stationnaire est définie par E Xt , mais pour t fixé, nous avons aucune répétition
d’échantillon pour estimer cette moyenne, l’ergodicité du processus stationnaire suppose que pendant un temps assez
long, la moyenne empirique temporelle va être approximativement indépendante et qui permettra une estimation
n
1X
E Xt = lim Xi
n→∞ n
i=1

Lorqu’on observe le processus suffisament longtemps, les estimations des paramètres des modèles considérés faites
à partir des observations sont alors valides et ces estimations ne seront pas influencées par des conditions initiales
particulières
Théorème 3.1.9 Tout processus stationnaire gaussien Xt tel que

Cov (Xt , Xs ) −→ 0 quand |t − s| → ∞

est ergodique
Pour les processus stationnaires considérés, ρk tens vers zéro quand k → ∞, alors ils sont ergodique s’ils sont gaussien.
Définition 3.1.10 Opérateur retard B Soit {Xt } un processus quelconque, on appelle opérateur retard l’opérateur
B défini par
BXt = Xt−1 ∀t
Ce qui permet de définir un nouveau processus B{Xt }. On note 1 l’opérateur identité et B 2 l’opérateur B ◦ B.

1Xt = Xt , B 2 Xt = Xt−2 , ∀t

Les 1, B, B 2 font parties de ce qu’on appelle opérateurs linéaires, qui définissent des nouveaux processus à partir de
combinaisons linéaires de processus. Plus généralement, on note P (B) (comme on parle de polynôme) de degré n, en
confondant la notation a0 avec celle de a01:

P (B) = a0 + a1 B + a2 B 2 + · · · + an B n

et
P (B)Xt = (a0 + a1 B + · · · + an B n )Xt = a0 Xt + a1 Xt−1 + · · · + an Xt−n
Les opérations habituelles sont valides, par example
△
(1 − B)2 Xt = (1 − 2B + B 2 )Xt

On fera usage de fractions rationnelles, de décomposition en éléments simples de ces fractions, de développements
polynomiaux. Par example, soit λ un nombre complexe de module inférieur à 1, on sait que

(1 − λx)−1 = 1 + λx + λ2 x2 + · · · + λn xn + · · ·

de sorte que avec |λ| < 1, on note

∞
X
(1 − λB)−1 = λi B i
0

ce qui permet de définir un nouveau processus par (1 − λB)−1 Xt .

Théorème 3.1.11 Des processus précédents tels que P (B)Xt et (1 − λB)−1 Xt sont stationnaires si le processus {Xt }
est stationnaire.

plus généralement, un opérateur f est un filtre linéaire quand c’est un opérateur invariant par changement de l’origine
du temps. Le processus filtré {Yt } défini par
Yt = f (Xt )
vérifie ainsi
Yt−τ = f (Xt−τ ) ∀t, ∀τ
est alors stationnaire si {Xt } est stationnaire.

3.2 Fonction d’autocorrelation partielle

Dans le chapitre 2.4, les coefficients de corrélation ρk est estimé par les fonctions d’autocorrélation rk . Sa représentation
graphique est le corrélogramme.

Remarque 3.2.1 On remarque que

1. Si ρk est proche de 1, c’est qu’une grande(resp. petite) valeur observée au temps t − k est suivie le plus souvent
d’une grande (resp. petite) valeur au temps t et ∀t.

2. Si ρk est proche de −1, c’est qu’une grande (resp. petite) valeur observée au temps t − k est suivie le plus souvent
d’une petite (resp. grande) valeur au temps t, et cela pour ∀t

3. Si ρk est proche de 0, qu’une grande valeur ou une petite valeur observée au temps t − k est suivie le plus souvent
d’une valeur quelconque au temps t

Soit {Xt } un processus tel que∗

Xt = φ11 Xt−1 + Rt,1 , ∀t
où φ11 ∈ R, la variable résidu Rt est non corrélée à la variable explicative Xt−1 . Soit en générale

Xt = φ1,k Xt−1 + φ2,k Xt−2 + · · · + φk,k Xt−k + Rt,k

et Rt,k est non corrélée aux variables explicatives Xt−1 , · · · , Xt−k . Cela implique que les variables du passé proche,
associées au premières variables explicatives, pourraient suffire pour expliquer totalement le présent. Dans le cas où
k=1
Xt Xt−1 = φ11 Xt−1 Xt−1 + Rt,1 Xt−1 ∀t
Donc pour ∀t ≥ 1
E Xt Xt−1 = φ11 E Xt−1 Xt−1 + E Rt,1 Xt−1 = φ11 Var Xt−1
c’est à dire φ11 = ρ1 de même, dans le cas où k = 2 on a

Xt Xt−1 = φ12 Xt−1 Xt−1 + φ22 Xt−2 Xt−1 + Rt,2 Xt−1

Xt Xt−2 = φ12 Xt−1 Xt−2 + φ22 Xt−2 Xt−2 + Rt,2 Xt−2

l’espérance mathématique
2
E Xt Xt−1 = φ12 E Xt−1 + φ22 E Xt−2 Xt−1

2
E Xt Xt−2 = φ12 E Xt−1 Xt−2 + φ22 E Xt−2

On divise tout par E Xt = σx2 et en utilisant le fait que Xt est stationnaire centrée, on a

 ρ1 = φ12 + ρ1 φ22 ρ2 − ρ21
=⇒ φ22 =

ρ2 = ρ1 φ12 + φ22 1 − ρ21

En générale la relation entre φkk et ρ dans le régressions

Xt = φ1,k Xt−1 + φ2,k Xt−2 + · · · + φk,k Xt−k + Rt,k

∗ régression linéaire de Xt en Xt−1
 1 ρ1 ρ2 · · · ρk−2 ρ1

ρ1 1 ρ1 · · · ρk−3 ρ2

··· ··· ··· ··· ··· ···

ρk−1 ρk−2 ρk−3 ··· ρ1 ρk k×k

φkk = (Yules-Walker)
1 ρ1 ρ2 · · · ρk−2 ρk−1

ρ1 1 ρ1 · · · ρk−3 ρk−2

··· ··· ··· ··· ··· ···

ρk−1 ρk−2 ρk−3 ··· ρ1 1 k×k

Par example si k = 2
1 ρ1

ρ1 ρ2 ρ2 − ρ21
φ22 = =
1 ρ1 1 − ρ21

ρ1 1
△
Définition 3.2.2 On appelle la fonction Φ(k) = φkk fonction autocorrélation partielle.
Exercice 3.2.3 Montrer le corrélogramme d’un bruit blanc est nul pour ∀k ≥ 1. son corrélogramme partiel Φ l’est
aussi.
nous dirons qu’un processus {Xt } est un processus centré engendré à partir d’un bruit blanc (noté {εt }), si et seulement
si on peut écrire
∞
X
Xt = ψj εj
0
P∞
avec ψ0 = 1, ψj ∈ R et 0 |ψj | < ∞ Tous les processus stationnaires considérés dans ce chapitre sont des processus
qui peuvent être considérés comme engendrés à partir d’un bruit blanc.
Définition 3.2.4 Le modèle auto-régressifs et moyenne mobile ARMA (autoregressive moving average), est
Xt = φ1 Xt−1 + · · · + φp Xt−p + εt + θ1 εt−1 + · · · + θq εt−q
où {εt } est un bruit blanc.
1. φ1 = φ2 = · · · = φp = 0 alors on appelle Xt moyenne mobile MA, si de plus θq 6= 0, alors il s’agit d’un processus
MA(q) d’ordre q.
2. θ1 = θ2 = · · · = θq = 0 alors on appelle Xt modèle autoregressif, si de plus φp 6= 0, alors il l’ordre de ce AR est
p, noté AR(p)
En générale, on note ce modèle par ARMA(p, q). On présente d’abord quelques modèles particuliers.

3.3 Modèle en moyenne mobile

3.3.1 modèles en moyenne mobile d’ordre 1
On considère dans ce paragraphe le modèle MA(1).
Xt = µ + εt + θεt−1
où θ est une constante et {εt } un bruit blanc. On peut alors calculer de manière explicite
E Xt = µ, Var Xt = (1 + θ2 )σ 2

Cov (Xt , Xt−1 ) = θσ 2 , Cov (Xt , Xt−j ) = 0, j > 0

Comme E Xt et Cov (Xt , Xt−j ) ne dépende pas de t, MA(1) est donc un processus stationnaire. Les coefficients
d’autocorrélation théorique sont (voir la figure 3.1)
θσ 2 θ
ρ1 = 2 2
= , ρ2 = ρ3 = · · · = 0
(1 + θ )σ (1 + θ2 )
 1

0.8

0.6

0.4

0.2

0.0

-0.2

-0.4

-0.6

2 4 6

Figure 3.1: Coefficients d’autocorrélation de MA(1).

3.3.2 modèles en moyenne mobile d’ordre q

On considère le processus
Xt = µ + εt + θ1 εt−1 + θ2 εt−2 + · · · + θq εt−q , θq 6= 0

Après calculs on a

E Xt = µ

2
Var Xt = E (εt + θ1 εt−1 + · · · + θq εt−q )

= (1 + θ12 + θ22 + · · · + θq2 )σ 2

Cov (Xt , Xt−j ) = E [(εt + θ1 εt−1 + · · · + θq εt−q ) × (εt−j + θ1 εt−j−1 + · · · + θq εt−j−q )]


 (θj + θj+1 θ1 + θj+1 θ2 + · · · + θq θq−j )2 j = 1, 2, · · · , q
=

0, j>q

Par example, pour MA(2) nous avons

γj = Cov (Xt , Xt−j ) Donc γ0 = (1 + θ12 + θ22 )σ 2

γ1 = (θ1 + θ2 θ1 )σ 2

γ2 = θ2 σ 2

γ3 = γ4 = · · · = 0

γj
On pourra donc calculer facilement les coefficients de corrélation ρj = γ0 , j = 1, 2, · · ·. le corrélogramme est alors
Question : quelle est la condition stationnarité ?

3.4 Modèle autorégressif

Nous commençons par

3.4.1 Modèle autorégressif du premier ordre

Soit
Xt = c + φXt−1 + εt
 1

0.8

0.6

0.4

0.2

0.0

-0.2

-0.4

-0.6

2 4 6

Figure 3.2: Corrélogramme MA(4), Xt = εt − 0.6εt−1 + 0.3εt−2 − 0.5εt−3 + 0.5εt−4

En fait si on prend Wt = c + ε, alors AR(1) est une équation différentielle Xt = φXt−1 + Wt

X1 = φX0 + W1 = φ2 X−1 + φW0 + W1

Xt = φt+1 X−1 + φt−1 W0 + · · · + φWt−1 + Wt et

Xt+j = φj+1 Xt−1 + φj Wt + φj−1 Wt+1 + · · · + φWt+j−1 + Wt+j

Supposons Xt−1 et Wt , Wt+1 , · · · , Wt+j connus alors le multiplicateur dynamique est

∂Xt+j
= φj
∂Wt
∂Xt+j
Remarque 3.4.1 Si Wt varie (et Wt+1 , · · · , Wt+j ne changent pas), l’effet† sur Xt+j est décrit par ∂Wt , d’où le
nom multiplicateur dynamique.
On distingue alors plusieurs cas possibles
∂Xt+j
1. 0 < φ < 1 alors ∂Wt décroit géométriquement vers 0.
∂Xt+j
2. −1 < φ < 0 alors ∂Wt change de signe de manière alternative et amplitude décroit géométriquement vers 0.
3. φ > 1 effet de Wt sur Xt+j est exponentiel par rapport à j, explosif.
4. φ < −1 oscillation explosive.
Donc si |φ| < 1 le système est stable, et si |φ| > 1 non stable. cela se traduit par
1. si |φ| < 1 alors Xt est stationnaire
2. si |φ| > 1 alors Xt est non stationnaire
3. si |φ| = 1 alors Xt est non stationnaire (à réflechir)
Si |φ| < 1 on a

Xt = (c + εt ) + φ(c + εt−1 ) + φ2 (c + εt−2 ) + φ3 (c + εt−3 ) + · · ·

c
= + εt + φεt−1 + φ2 εt−2 + · · · (représentation MA(∞) de AR)
1−φ
c
E Xt = 1−φ

Var Xt = E (Xt − E Xt )2 = E (εt + φεt−1 + φ2 εt−2 + · · ·)2

σ2
= (1 + φ2 + φ4 + φ6 + · · ·)σ 2 = 1−φ2

† il s’agit d’effet du présent sur la future

Attention, contrairement aux chapitres précedents, on a considéré la série comme

(· · · , X−2 , X−1 , X0 , X1 , X2 , · · · , Xt ).

Par ailleurs ∀j ≥ 0
γj = Cov (Xt , Xt−j )

= E (εt + φεt−1 + · · · + φj εt−j + φj+1 εt−j−1 + φj+2 εt−j−2 + · · ·)

×(εt−j + φεt−j−1 + φ2 εt−j−2 + · · ·) σ 2

  φj
= φj + φj+2 + φj+4 + · · · = σ2
1 − φ2
Les coefficients autocorrélations sont (voir la figure 3.3 par example)
γj
ρj = = φj
γ0

1 1

0.8 0.8

0.6 0.6

0.4 0.4

0.2 0.2

0.0 0.0

-0.2 -0.2

-0.4 -0.4

-0.6 -0.6

2 4 6 2 4 6

Figure 3.3: Corrélogramme AR(1) : Xt = 0.8Xt−1 + εt et AR(1) : Xt = −0.8Xt−1 + εt

3.4.2 Modèle autorégressif du second ordre

On considère dans cette paragraphe le processus AR(2)

Xt = c + φ1 Xt−1 + φ2 Xt−2 + εt

La condition pour que AR(2) soit stationnaire est que l’équation 1 − φ1 z − φ2 z 2 = 0 possède deux racines (réelles ou
complexes) de module > 1.‡ On calcule maintenant quelques statistiques pour AR(2).
c
µ = E Xt = c + φ1 µ + φ2 µ =⇒ µ =
1 − φ1 − φ2
Pour calculer E Xt2 (ou Var Xt ) on utilise
Xt = µ(1 − φ1 − φ2 ) + φ1 Xt−1 + φ2 Xt−2 + εt

Xt − µ = φ1 (Xt−1 − µ) + φ2 (Xt−2 − µ) + εt
multipli Xt−j − µ et calculer l’espérance on a

γj = φ1 γj−1 + φ2 γj−2 , j = 1, 2, · · ·

diviser tout par γ0 on a aussi

ρj = φ1 ρj−1 + φ2 ρj−2
en particulier, j = 1 et j = 2
φ1
ρ1 = φ1 + φ2 ρ1 =⇒ ρ1 =
1 − φ2

ρ2 = φ1 ρ1 + φ2
‡ On admet ce résultat, il faut étudier la condition de stabilité d’équation différentielle Xt = φ1 Xt−1 + φ2 Xt−2 + Wt .
On peut utiliser ces méthdoes pour trouver
(1 − φ1 )σ 2
γ0 = exercice
(1 + φ2 )[(1 − φ2 )2 − φ21 ]
En effet, il suffit de calculer
(Xt − µ)2 = φ21 (Xt−1 − µ)2 + φ22 (Xt−2 − µ)2 + ε2t
+2φ1 φ2 (Xt−1 − µ)(Xt−2 − µ) + 2φ1 (Xt−1 − µ)εt + 2φ2 (Xt−2 − µ)εt
Xt − µ étant stationnaire et moyenne nulles, εt un bruit blanc indépendant de Xt−1 et Xt−2 , prenons espérance
mathématique nous avons
γ0 = φ21 γ0 + φ22 γ0 + σ 2 + 2φ1 φ2 γ1 = φ21 γ0 + φ22 γ0 + σ 2 + 2φ1 φ2 ρ1 γ0
Donc  
1 − φ21 − φ22 − 2φ1 φ2 ρ1 γ0 = σ 2
φ1
Pour montrer le résultat, il suffit alors de remplacer ρ1 par 1−φ 2
. On peut aussi montrer la relation entre (φ1 , φ2 ) et
(γ1 , γ2 ) (exercice)
φ1 = ρ11−ρ(1−ρ2 )
2
1

ρ2 −ρ21
φ2 = 1−ρ21
(Yule-Walker).
la condition stationnarité pour AR(2) est présentée par la figure 3.4.

φ2
1

-1 1 φ1

-2 2

-1
Figure 3.4: Condition de stationnarité pour un modèle de AR(2)

3.4.3 Modèle autorégressif d’ordre p

AR(p) : Xt = c + φ1 Xt−1 + φ2 Xt−2 + · · · + φp Xt−p + εt

p
X
Soit le polynôme caractéristique Φ(z) = 1 − φi z k z ∈ C ne possède pas de racines réelles ou complexes de module
1
inférieur ou égale à 1 (condition stationnarité). On peut effectuer quelques calculs simples
E Xt = µ = c + φ1 µ + φ1 µ + · · · φp µ
c
E Xt = µ=
1 − φ1 − φ2 − · · · − φp
En écrivant AR(p) par
Xt − µ = φ1 (Xt−1 − µ) + φ2 (Xt−2 − µ) + · · · + φp (Xt−p − µ) + εt
Multipli par Xt−j − µ et calculer l’espérance, on obtient

 φ1 γj−1 + φ2 γj−2 + · · · + φp γj−p j = 1, 2, · · ·
γj =

φ1 γ1 + φ2 γ2 + · · · + φp γp j=0
diviser par γ0 on obtient l’équation de Yule-Walker:

ρj = φ1 ρj−1 + φ2 ρj−2 + · · · + φp ρj−p , j = 1, 2, · · ·

3.4.4 Inversibilité
On présente dans la suite la notion de inversibilité pour certain ARMA, elle est nécessaire pour identifier un ARMA.
Pour comprendre cette notion, on introduit d’abord la notion représentation AR(∞) d’un modèle ARMA. Soit
Φ(B)Xt = Θ(B)εt .

Représentation MA(∞) et AR(∞) d’un ARMA

Pour un modèle ARMA(p, q),

Xt = φ1 Xt−1 + φ2 Xt−2 + · · · + φp Xt−p + εt + θ1 εt−1 + · · · + θq εt−q

on peut obtenir une représentation MA(∞)

∞
X
Xt = ψj εt−j ou Xt = Ψ(B)εt
j=0

Par example, pour ARMA(1, 1) on peut

Xt = φXt−1 + εt + θεt−1 = φ(φXt−2 + εt−1 + θXt−2 ) + ε2 + θεt−1

= φ2 Xt−2 + εt + (θ + φ)εt−1 + θεt−2

= φ2 (φXt−3 + εt−2 + θεt−3 ) + εt + (θ + φ)εt−1 + θεt−2 = etc = etc

Il y a une deuxième méthode (plus efficace). Soit Ψ(B) = ψ0 +ψ1 B+ψ2 B 2 +· · · on identifie les coefficients ψ0 , ψ1 , ψ2 , · · ·
par équation
Φ(B)Ψ(B)εt = Θ(B)εt
Dans l’example précédent, on obtient (1 − φB)Ψ(B) = 1 + θB on obtient donc

ψ0 = 1 = coefficient devant B 0

ψ1 = θ + φ = coefficient devant B

ψ2 = φψ1 = φ(θ + φ) = coefficient devant B 2

etc etc

donc
Xt = εt + (θ + φ)εt−1 + φ(θ + φ)εt−2 + · · ·
La troisième méthode consiste à inverser AR polynôme.

(1 − φB)Xt = (1 + θB)εt alors Xt = (1 − φB)−1 (1 + θB)εt

Xt = (1 + φB + φ2 B 2 + · · ·)(1 + θB)

Si p ≥ 2, on doit alors inverser un polynôme de degré ≥ 2 i.e.

Xt = (1 − φ1 B − φ2 B 2 − · · · − φp B p )−1 (1 + θB)εt

Or l’inversion d’un polynôme est souvent difficile, on suppose d’abord ∃ un factorisation de Φ(B) = (1 − φ1 B − φ2 B 2 −
· · · − φp B p ) :
    
1 1 1
Φ(B) = 1 − B 1 − B ··· 1 − B
r1 r2 rp
où r1 , r2 , · · · , rp sont les racines de Φ(B), alors
a1 a2 ap
Φ−1 (B) = + + ···+
1 − r11 B 1 − r12 B 1 − r1p B
  2
1 1 1
et on écrit pour chaque terme de Φ−1 (B) par 1 = 1+ B+ B + · · · pour obtenir la representation
1− r1 B
r1 r2
inifinie MA(∞).
Le modèle ARMA a une représentation AR(∞) également :
∞
X
εt = ξj Xt−j
j=0

Par example pour ARMA(1, 1)

εt = Xt − φXt−1 − θεt−1 = Xt − φXt−1 − θ(Xt−1 − φXt−2 − θεt−2 )

= Xt − (φ + θ)Xt−1 + θ(φ + θ)Xt−2 − θ2 (φ + θ)Xt−2 + · · · (−θ)k−1 (φ + θ)Xt−k + · · ·

On peut biensur utiliser l’identification du polynôme Θ(B). Soit εt = Ξ(B)Xt où Ξ(B) = ξ0 + ξ1 B + ξ2 B 2 + · · · est à
identifier. L’équation Φ(B) = Θ(B)Ξ(B) on obtient pour ARMA(1, 1)

ξ0 = 1, ξ1 = −(θ + φ), ξ2 = θξ1 , ξ3 = −θξ2 , · · ·

La troisième méthode consiste aussi d’inverser un polynôme εt = Θ−1 (B)Φ(B)Xt pour ARMA(p, 1) on utilise (1 +
θB)−1 = 1 − θB + θ2 − θ3 B 3 + · · ·. Le cas où q ≥ 2 est difficile.

Inversibilité d’un modèle ARMA

∞
X
L’existance d’une représentation AR(∞) : εt = ξj Xt−j est importante pour l’estimation et la prédiction (c’est
j=0
même essentielle si on veut faire de calculs numériques). Par example : pour un modèle MA(1) on a

Xt = (1 + θB)εt = εt + θεt−1

εt = (1 + θB)−1 Xt = (1 − θB + θ2 B 2 − θ3 B 3 + · · ·)Xt

· Cas 1, si θ > 1, le poid de Xt−k est (−θ)k croissant, on a desoin de tous les passés pour calculer εt .
P∞
· Cas 2, si θ = 1, on a besoin de tous les passés Xt−1 , Xt−2 , · · · pour calculer εt et j=0 |ξj | diverge
· Cas 3, |θ| < 1, alors le poids de Xt−k décroit vers 0, on appelle ce cas inversible.

Théorème 3.4.2 La condition pour que un modèle ARMA(p, q) (Φ(B)Xt = Θ(B)εt ) soit stationnaire est que tous
les racines du polynôme Φ(z) (réelles ou complexes) ont des modules supérieurs à 1. La condition pour que un modèle
ARMA soit inversible est que tous les racines du polynôme Θ(B) (réelles ou complexes) ont des modulès supérieurs à
1.

Preuve: Admis. En effet, cela vient du fait que il faut

· Stationnarité : avoir une représentation convergente de MA(∞).
· Inversibilité : avoir une représentation convergente de AR(∞).
On remarque que le calcul des corrélogramme d’un modèle de ARMA mixte ne permet pas de identifier l’ordre (p, q)
d’un ARMA(p, q).
Chapter 4

Estimation, Identification et Prévision

des modèles ARIMA

4.1 Introduction
L’identification de l’ordre (p, q) du modèle mixte est délicate dans la mesure où les méthodes à mettre en oeuvre sont
tous imprécises ou deviennent si coûteuses quand la précision est recherchée, qu’elle sont pas toujours portés sur les
logiciels. Il existe des méthodes heuristiques comme par example la méthode dite du coin. Une autre approche consiste
à ne s’interesser aux modèles mixtes qu’après avoir constaté qu’il ne peut s’agiter ni d’un AR ni d’un MA :

· l’ordre p d’un AR(p) peut se lire sur corrélogramme Φ des autocorrélations partielles, ce corrélogramme étant
nulle à partir de k = p + 1.

· l’ordre q d’un MA(q) se lit sur corrélogramme, se comportant comme un corrélogramme nul à partir des décalages
superieurs ou égaux à q + 1 et non nul en deça.

avec p, q données, il reste à estimer les paramètres de ce modèle en utilisant les observations (X1 , · · · , Xt ): la variation
du bruit blanc, les coefficient φ1 , φ2 , · · · , φp et θ1 , θ2 , · · · , θq . En pratique, la valeur de p et q sont rarement ≥ 3. Alors
que signifie une bonne estimation ? Les critères sont

• sans biais, • petite erreur (au sens moyenne quadratique)

Un des objectifs principaux est de prédire l’avenir :

E (Xt+h |X1 , X2 , · · · , Xt )

4.2 Estimation
On présente dans cette section deux méthodes d’estimation différentes, on considère tout d’abord des examples simples.

4.2.1 Méthode de moindre carré

Soit
Xt = φ1 Xt−1 + φ2 Xt−2 + εt
Les paramètres à estimer sont σ 2 , φ1 et φ2 , on dispose des observations sous forme d’un vecteur de taille N :
(X1 , · · · , XN ). Il s’agit donc de chercher φ1 et φ2 t. q.

N
X
(Xt − φ1 Xt−1 − φ2 Xt−2 )2 soit minimale (4.1)
t=1

Les variables X−1 et X0 ne sont pas disponibles, on cherche

· φˆ1 , φˆ2 t. q.
N
X
(Xt − φˆ1 Xt−1 − φˆ2 Xt−2 )2 soit minimale
t=3

On suppose en fait que les 2 premiers termes dans 4.1 sont nuls
 · La valeur de σ 2 sont estimée par
N
1 X
σˆ2 = (Xt − φˆ1 Xt−1 − φˆ2 Xt−2 )2
N − 2 t=3

On peut facilement appliquer le principe de ces méthodes pour des problèmes plus généraux, la méthode de moindre
carré est facile à mettre en oeuvre, elle est valide quand la taille N est importante.

4.2.2 Méthode de maximum de vraisemblance

La méthode de maximum de vraisemblance est une des méthodes la plus utilisée en statistique, elle est très intéressant
car sous des conditions assez générales, l’estimateur est asymptotiquement sans biais et de variance minimale. On
suppose au départ que les {Xt }, t = 1, · · · , N sont indépendantes de même distribution f (x, Θ), alors la fonction de
vraisemblance est
N
Y
f (X1 , X2 , · · · , XN , Θ) = f (Xt , Θ)
t=1
L’estimateur EMV est obtenu en maximisant cette fonction. Le problème dans les séries temporelles est que les
{Xt } ne sont pas indépendantes, On peut néanmoins écrire la loi conjointe f (X1 , X2 , · · · , XN , Θ) à l’aide de densité
conditionnelle.
f (X1 , · · · , XN ) = f (XN |XN −1 )f (X1 , X2 , · · · , XN −1 )

=f (XN |XN −1 )f (XN −1 |XN −2 )f (X1 , X2 , · · · , XN −2 )

Q
où X = (X1 , X2 , · · · , Xn ), on obtient f (X1 , · · · , XN ) = N
t=2 f (Xt |Xt−1 )f (X1 ), produit des densités conditionnelles.
la méthode de MV est applicable si le modèle considéré nous permet de calculer ces densités.
En effet, c’est le cas des modèles ARMA sous certains conditions. Par example dans le cas où AR(p) est gaussien.
En particulier, soit εt ∈ N (0, σ 2 ), et que X0 est gaussienne, alors la densité conditionnelle d’une v. a. gaussienne
sanchant une autre variable gaussienne est une densité gaussienne∗ .

Examples
On étudie maintenant quelques examples.
Example 4.2.1 Soit un modèle AR(1) gaussien

Xt = c + φXt−1 + εt εt i.i.d. et εt ∈ N (0, σ 2 )

Les paramètres à estimer sont Θ = (c, φ, σ 2 ). Calculer son EMV.

Sol La méthode de MV consiste donc à
· Evaluer la fonction de vraisemblance

fX1 ,X2 ,···,XN (x1 , x2 , · · · , xN , Θ)

· Trouver l’estimateur qui maximise cette fonction.

Pour commencer, nous considérons la distribution de X1 , d’après les résultats du chapitre 3 nous avons
c σ2
µ = E Xt = , Var Xt =
1−φ 1 − φ2
 
c σ2
donc X1 est gaussienne, X1 ∼ N , :
1 − φ 1 − φ2
c 2
!
1 −(x1 − 1−φ )
fX1 (x1 , Θ) = √ p exp
2π σ/1 − φ2 2σ 2 /1 − φ2

Considérons maintenant la distribution de X2 sachant X1 = x1 , X2 = c + φX1 + ε2 . le conditionnement X1 = x1

est traité comme si X1 est une constante déterministe donc la distribution conditionnelle est (X2 |X1 = x1 ) ∼ N (c +
φx1 , σ 2 ), c’est à dire
 
1 (x2 − c − φx1 )2
fX2 |X1 (x2 |x1 , Θ) = f (x1 , x2 , Θ|x1 ) = √ exp −
2πσ 2σ 2
∗ C’est une propriété très importante pour les vecteurs gaussiens, c’est en partie grâce à cette propriété, on a trouvé le filtre de Kalman,

voir le chapitre 7.
d’une manière récurrente, la distribution de XN sachant X1 = x1 , X2 = x2 , · · · , XN −1 = xN −1 est
(XN |X1 = x1 , X2 = x2 , · · · , XN −1 = xN −1 , Θ) ∼ N (c + φxN −1 , σ 2 )
 2

fXN |X1 ,···,XN −1 (xN |x1 , · · · , xN −1 , Θ) = √1
2πσ
exp − (xN −c−φx
2σ 2
N −1 )

car l’influence de (X1 = x1 , X2 = x2 , · · · , XN −1 = xN −1 ) sur XN est faite seulement via XN −1 = xN −1 , et

XN = c + φXN −1 + εN
N
Y
Finalement, la formule f (X1 , X2 , · · · , XN , Θ) = f (Xt |Xt−1 ) · f (X1 ) se traduit par
t=1

N
Y
f (x1 , x2 , · · · , xN , Θ) = fX1 (x1 , Θ) fXt |Xt−1 (xt |xt−1 , Θ)
t=2

N
Y
L(x1 , x2 , · · · , xN , Θ) = fX1 (x1 , Θ) fXt |Xt−1 (xt |xt−1 , Θ)
t=2

1 1 [x1 − (c/1 − φ2 )]2

ln L(x1 , x2 , · · · , xN , Θ) = − ln 2π − ln(σ 2 /1 − φ2 ) −
2 2 2σ 2 /1 − σ 2

N  
N −1 N −1 2
X (xt − c − φxt−1 )2
− ln 2π − ln σ −
2 2 t=2
2σ 2

La distribution f (x1 , Θ) prend une place importante dans l’expression L(xN , Θ). Souvent quand N est grand, on
peut négliger X1 en fixant x1 : c’est à dire on considère X1 = x1 comme connu et fX1 (x1 , Θ) est simplifiée par
fX1 (x1 , Θ) = 1, alors la formule ln L est beaucoup plus facile à évaluer. Maximiser la fonction ln L(x1 , x2 , · · · , XN , Θ)
par rapport à (c, φ) en considérant σ 2 comme connu est équivalent à minimiser la fonction
N
X
(xt − c − φxt−1 )2
t=2

Donc l’EMV pour AR(1) est (exercice)

   P  P 
ĉ N −1 Xt−1 Xt X N
X
 =   où =
P P 2
P
φ̂ Xt−1 Xt−1 Xt−1 Xt t=2

∂
En faisant ∂σ2 ln L(x1 , · · · , xn , σ 2 ) = 0 on obtient
N  
T −1 X (xt − c − φxt−1 )2
− + =0
2σ 2 t=2
2σ 4
N
!
X (xt − ĉ − φ̂xt−1 )2
donc σˆ2 = (exercice)
t=2
N −1

Example 4.2.2 On considère le modèle Xt = φ1 Xt−1 + φ2 Xt−2 + εt , où εt est un bruit blanc gaussien, alors la densité
conditionnelle de Xt sachant X1 , · · · , Xt−1 est normale, la moyenne E (Xt |Xt−1 ) est égale à φ1 Xt−1 + φ2 Xt−2 et la
variance est
E (Xt − (φ1 Xt−1 + φ2 Xt−2 ))2 = E ε2t = σ 2
c’est à dire  
1 (Xt − φ1 Xt−1 − φ2 Xt−2 )2
f (Xt |X1 , X2 , · · · , Xt−1 , Θ) = √ exp −
2πσ 2σ 2
Montrer que l’EMV de Θ = (φ1 , φ2 , σ 2 ) est
 N
 X
 (φˆ1 , φˆ2 ) = (Xt − φ1 Xt−1 − φ2 Xt−2 )2


 arg min

 t=3


 N
 1 X
σˆ2 (Xt − φˆ1 Xt−1 − φˆ2 Xt−2 )2


 =
 N − 3 t=3
Sol Etant donné que X−1 , X−2 et X0 sont inconnus, et que les densités f (X1 ), f (X2 |X1 ) sont fixées à 1, on considère
N
Y
L(XN , Θ) = f (Xt |Xt−1 )
t=3

N
N −2 1 X
ln L(XN , φ1 , φ2 ) = − ln(2πσ 2 ) − 2 (Xt − φ1 Xt−1 − φ2 Xt−2 )2
2 2σ t=3

En utilisant la même méthode que l’example 4.2.1, on obtient le résultat. 2

Example 4.2.3 Modèle AR(p) : Xt = c − φ1 Xt−1 − · · · − φp Xt−p + εt . Les paramètres sont Θ = (c, φ1 , · · · , φp , σ 2 ).

Sol Maximiser la fonction de vraisemblance pour un AR(p) peut être réalisé par des méthodes numériques, on cherche
alors à maximiser la fonction
ln fXN ,XN −1 ,···,Xp+1 |Xp ,Xp−1 ,···,X1 ( xN , xN −1 , · · · , xp+1 | xp , · · · , x1 , Θ)

N
N −p X (xt − c − φ1 xt−1 − φ2 xt−2 − · · · − φp xt−p )2
=− ln(2πσ 2 ) −
2 t=p+1
2σ 2

d’où  N  2
1 X
σˆ2



 = Xt − ĉ − φˆ1 Xt−1 − φˆ2 Xt−2 − · · · − φˆp Xt−p


 N − p t=p+1


 N
 X 2
 ˆ ˆ ˆ
 (φ1 , φ2 , · · · , φp , ĉ) = arg min (Xt − c − φ1 Xt−1 − φ2 Xt−2 − · · · − φp Xt−p )


t=p+1
2
Example 4.2.4 Calculer la fonction de vraisemblance pour un MA(1).
Sol Soit
Xt = µ + εt + θεt−1 où εt ∼ N (0, σ 2 ) i.i.d.
Les paramètres à estimer sont Θ = (µ, θ, σ 2 ). Dans le cas gaussien
 
1 (xt − µ − θεt−1 )2
Xt |εt−1 ∼ N (µ + θεt−1 , σ 2 ) = √ exp −
2πσ 2σ 2
On suppose que ε0 = 0 alors
 
2 1 (x1 − µ)2
X1 |ε0 = 0 ∼ N (µ, σ ) = √ exp − et ε1 = x1 − µ
2πσ 2σ 2
 
1 (x2 − µ − θε1 )2
X2 |X1 , ε0 = 0 ∼ fX2 |X1 ,ε0 =0 (x2 |x1 , ε0 , Θ) = √ exp −
2πσ 2σ 2

et ε2 = x2 − µ − θε1

··· ··· ···

Xt |Xt−1 , · · · , X1 , ε0 = 0 ∼ fXt |Xt−1 ,···,X1 ,ε0 =0 (xt |xt−1 , · · · , x1 , ε0 , Θ) = fXt |εt−1 (xt |εt−1 , Θ)
   
1 (xt − µ − θεt−1 )2 1 ε2t
= √ exp − = √ exp − 2
2πσ 2σ 2 2πσ 2σ
Donc la fonction de vraisemblance est
 N  N (N )
1 1 X ε2
t
f (XN , XN −1 , · · · , X1 |ε0 = 0) = √ √ exp et εt = xt − µ − θεt−1
2π σ2 t=1
2σ 2

Les εt sont calculées de manière itérative. 2

On remarque que contrairement aux problèmes 4.2.1, 4.2.2 et 4.2.3 où on peut calculer les EMV Θ̂ = arg min ln L(xN , Θ)
par des méthodes directes, dans l’example 4.2.4, on ne pourra pas le faire. Pour cela on présentera dans la section
suivante une méthode implicite itérative.
Example 4.2.5 Processus MA(q), Θ = (µ, θ1 , θ2 , · · · , θq , σ 2 ).

Sol ln L(Θ) est

N
N X ε2t
ln L(XN , Θ) = ln f (XN , XN −1 , · · · , X1 |~ε0 = 0, Θ) = − ln(2πσ 2 ) −
2 t=1
2σ 2

Où le q-vecteur ~ε0 = (ε0 , ε−1 , · · · .ε−q+1 ) et les {εt } sont obtenues par

εt = Xt − µ − θ1 εt−1 − · · · − θq εt−q

Example 4.2.6 Processus ARMA(p, q). Soit

Xt = c + φ1 Xt−1 + φ2 Xt−2 + · · · + φp Xt−p + εt + θ1 εt−1 + · · · + εq εt−q

Θ = (c, φ1 , φ2 , · · · , φp , θ1 , θ2 , · · · , θq , σ 2 )

Sol Soit 
~ε0 = (ε0 , ε−1 , · · · , ε−q+1 ) 
condition initiale

~x0 = (x0 , x−1 , · · · , x−p+1 )
On peut calculer successivement pour t = 1, 2, · · · , N

εt = xt − c − φ1 xt−1 − φ2 xt−2 − · · · − φp xt−p − θ1 εt−1 − θ2 εt−2 − · · · − θq εt−q

Les conditions initiales sont choisies :

p
z }| {
~ε0 = (0, 0, · · · , 0)
q
z }| {
c c c
~x0 = ( , ,···, )
1 − φ1 − φ2 − · · · − φp 1 − φ1 − φ2 − · · · − φp 1 − φ1 − φ2 − · · · − φp

Elles représentent espérances mathématiques de ~ε0 et ~x0 . On en déduit alors

N
N X ε2t
ln L(XN , Θ) = f (XN , XN −1 , · · · , X1 |~x0 , ~ε0 , Θ) = − ln(2πσ 2 ) −
2 t=1
2σ 2

2
Box et Jenkin suggèrent une formulation légèrement différente : les εt ne sont calculées qu’à partir de t = p + 1, p +
2, · · · , N avec εp−q+1 = εp−q = · · · = εp−1 = εp = 0 et (x1 , x2 , · · · , xp ) (p valeurs observée) comme conditions initiales,
pour t = p + 1, · · · , N
εt = xt − c − φ1 xt−1 − · · · − φp xt−p − θ1 εt−1 − · · · − θq εt−q
La fonction ln L(Θ) devient

ln L(XN , Θ) = f (XN , XN −1 , · · · , X1 |Xp , Xp−1 , · · · , X1 , εp = εp−1 = · · · = εp−q+1 = 0, Θ)

N
N −p X ε2t
= − ln(2πσ 2 ) −
2 t=p+1
2σ 2

4.3 Méthode numérique

On présente dans cette section une méthode itérative pour trouver l’EMV. Supposons qu’on dispose comme dans la
figure 4.6 d’une boı̂te noire subroutine), capable de calculer à chaque valeur de Θ et les observations (x1 , x2 , · · · , xN )
la fonction ln L(xN , Θ) correspondante. L’objectif est de trouver une suite Θm qui converge vers Θ̂ maximisant
△
L(Θ) = ln L(Θ).
 xN Calculer ln L(xN , Θ )

Input ln L(xN , Θ ) Output

Figure 4.1: Boı̂te noire pour calculer ln L(xN , Θ)

4.3.1 Méthode de recherche par grille

C’est la méthode la plus simple, on établi d’abord une grille (qui discrétise le domaine admissible de Θ), on teste
ensuite point par point de ce domaine, cela permet d’identifier arg max L(Θ). En pratique

· Il faut que la grille soit assez fine, pour capter le comportement de L(Θ) par rapport à Θ.
· Raffiner sur la zone où arg max se trouve.

· Fixer un critère d’arrêt pour le raffinement.

Example 4.3.1 En dim = 1, AR(1) : Xt = c + φXt−1 + εt , on suppose que c = 0, σ 2 = 1 est connu. Alors selon
l’example 4.2.1
N
N 1 1X
ln L(φ) = − ln(1 − φ2 ) − (1 − φ2 )x21 − (xt − φxt−1 )2
2 2 2 t=2

Soit xN = (x1 = 0.8, x2 = 0.2, x3 = −1.2, x4 = −0.4, x5 = 0.0), La fonction L(φ) est illustrée par la figure 4.3

xN Calculer ln L(xN , Θ )

Input ln L(xN , Θ ) Output

Figure 4.2: Méthode de recherche par grille

Cette méthode devient très couteuse lorsque dim ≥ 2.

4.3.2 méthode de descente gradiant

Pour maximiser la fonction L(Θ), on définit
∂L(Θ)
g(Θ) = (g1 (Θ), g2 (Θ), · · · , ga (Θ) = ∂Θ = DΘ L(Θ)

= gradΘ L(Θ)

La méthode consiste donc à calculer successivement la suite

Θm+1 = Θm + s · g(Θm ). s : paramètre de contrôle

Si on ne dispose pas de g(Θ) théorique, on se contentera d’avoir une approximation† :

1 
gi (Θm ) = L(Θm m m m m m
1 , Θ2 , · · · , Θi−1 , Θi + ∆, Θi+1 , · · · , Θa )
∆

− L(Θm m m m m m
1 , Θ2 , · · · , Θi−1 , Θi , Θi+1 , · · · , Θa )

† gradiant numérique
 1 1 1 1
, , , , 1, 2, 4, 8 etc (selon problème),
Pour i = 1, · · · , a. La valeur de s est choisie de la façon suivante : on test s =
m m
16
m
8 4 2
on prendra la veleur pour la quelle L(Θ + sg(Θ )) est maximale, quand Θ est proche de la solution finale, il faut
alors raffiner encore plus. Enfin, pour éviter de tomber dans des maximums locaux, on effectue ces calculs avec
plusieurs valeurs initiales distincts.

g(θ1 )
g(θ2 )

θ3

Figure 4.3: Méthode de descente gradiant

4.3.3 Méthode de Newton

Soit H(Θm ) une matrice définie par
∂ 2 L(Θ)
H(Θm ) = −
∂Θ∂Θ Θ=Θm
La méthode de Newton consiste à itérer
Θm+1 = Θm + [H(Θm )]−1 g(Θm )
On remarque que si la fonction L(Θ) est une forme quadratique, alors une seule itération suffise : Θ1 = EMV. En
générale, on utilise la formule
Θm+1 = Θm + s[H(Θm )]−1 g(Θm ) s : paramètre de contrôle

4.4 Identification
4.4.1 Choix de p et q dans un processus AR(p) et MA(q)
Comme signalé dans l’introduction, les indicateurs pour le choix de p et q sont essentiellement la forme des fonctions
d’autocorrélation partielle Φ(k) et les corrélogrammes rk . On peut penser par example qu’une chute de Φ(k) vers de
valeurs proches de zéro pour k > p indique un AR(p), ou qu’une chute rk vers des valeurs proches de zéro pour k > q
indique un MA(q). Pour juger si les Φ(k) ou rk sont significativement différents de zéro, il est bon d’avoir une idée de
leur écart-type. On peut utiliser deux résultats
1. Pour un MA(q), la variance de rk vérifie
" q
#
1 X
∀k > q, Var rk ∼
= 1−2 2
rk
T
i=1

2. Pour un AR(p), la variance de Φ(k) vérifie

1
∀k > p, Var Φ(k) ∼
=
T
L’application de ces résultats sont assez simple, par example pour identifier l’ordre p d’un processus autorégressif,
on peut calculer
 la suite
 des autocorrélations partielles estimées Φ(k) et examiner à partir quelle valeur elle reste dans
1.96 1.96
la bande − √ , √ .
T T
Ainsi dans l’example ci-dessus, on peut considérer p = 3. Cette démarche consiste donc à tester la nullité des
autocorrélations partielle Φ(k) à partir d’un certain rang.
 1

1.96
T
k

1.96
T
2 4 6

Figure 4.4: Corrélogramme des autocorrélations partielles

4.4.2 Choix de p et q dans le cas d’un processus ARMA(p, q)

Nous présentons la méthode de coin en admetttant le résultat suivant
Théorème 4.4.1 Un processus ARMA(p, q) stationnaire satisfait la propriété suivant (propriété du coin).

 ∆(i, j) = 0


i ≥ p, j ≥ q


∆(i, q − 1) 6= 0 i≥p−1





∆(p − 1, j) 6= 0 j ≥q−1
avec
γ(j + 1) γ(j + 2) · · · γ(j + i + 1)

γ(j) γ(j + 1) · · · γ(j + i)

∆(i, j) = det
··· ··· ··· ···

γ(j + 1 − i) γ(j − i) ··· γ(j + 1) (i+1)×(i+1)

Les valeurs des ∆(i, j) sont inconnues, elles peuvent être estimées en remplaçant dans le déterminant les corrélations
γ(k) par leurs estimateurs rk . Il est vraisemblable que dans le tableau dont les éléments sont de ∆(i, j) on verra
appaı̂tre une rupture entre les lignes i = p − 1 et i = p, et entre les colonnes j = q − 1 et j = q.
En fait, ces nullités sont vérifiées par des tests statistiques (test student ou test chi-2). En pratique, à l’issue
de ces tests, un ensemble de modèles correspondants aux différents choix de (p, q) seront acceptés. On peut par
example accepter simultanément plusieurs modèles: MA(3), ARMA(1, 2), ARMA(3, 1) et AR(5). Après cet étape de
identification (dite Identification à priori), on procède alors aux estimations des paramètres de chacun de ces modèles,
une deuxième étape d’identification (Identification a posteriori) est alors effectuée afin de choisir le ou les meilleurs
modèles. Le schéma de cette procédure peut être résumé par l’organigramme 4.5 ci-dessous

4.5 Prévision
4.5.1 Principe de prévision
∗
On s’interesse à la prévision de la v.a. Xt+1 basée sur le vecteur Xt = (Xt , Xt−1 , · · · , Xt−m+1 ). On désigne Xt+1|t la
prévision de Xt+1 basée sur l’observation Xt . Il faut d’abord mesurer la qualité d’une prévision.
Définition 4.5.1 Erreur quadratique moyenne (MSE) On appelle MSE ( erreur quadratique moyenne) d’une
∗
prévision Xt+1|t l’espérance
∗ △ ∗
MSE Xt+1|t = E (Xt+1 − Xt+1|t )2
L’objectif est alors de minimiser cette erreur.
∗
Proposition 4.5.2 La meilleure prévision Xt+1|t au sense de MSE est l’espérance conditionnelle de Xt+1 sachant
Xt .
∗
Xt+1|t = E (Xt+1 |Xt )
 Identification
a priori

Estimation des
modèles retenus

vérification
Identification
a priori

Choix de modèle
(s’il reste plusieurs modèles)

Prévisions

Figure 4.5: Schéma : Identification, estimation, prévision

Preuve‡ Soit g(Xt ) une prévision quelconque Xt+1|t

∗
= g(Xt ) alors MSE correspondante est

2
E (Xt+1 − g(Xt ))2 = E [Xt+1 − E (Xt+1 |Xt ) + E (Xt+1 |Xt ) − g(Xt )]

2
= E [Xt+1 − E (Xt+1 |Xt )]

+2E {[Xt+1 − E (Xt+1 |Xt )][E (Xt+1 |Xt ) − g(Xt )]}

2
+E [E (Xt+1 |Xt ) − g(Xt )]

= E I12 + E I2 + E I32

or E I2 = E (E (I2 |Xt )) et

E (I2 |Xt ) = [E (Xt+1 |Xt ) − g(Xt )] E {[Xt+1 − E (Xt+1 |Xt )]|Xt }

car E (Xt+1 |Xt ) − g(Xt ) est Xt mesurable§ , donc E I2 = [E (Xt+1 |Xt ) − g(Xt )] × 0 = 0.
E I12 ≥ 0 ne dépend pas du choix de g, et E I32 ≥ 0 atteint le minimum avec g(Xt ) = E (Xt+1 |Xt ) donc la meilleure
prévision est E (Xt+1 |Xt ). 2

4.5.2 Espace L2 , régression affine

R 1
Définition 4.5.3 X ∈ L2 (Ω, A, P ) (processus carrée intégrable) si ||Xt ||2 = ( |x|2 dP ) 2 < ∞.

Il s’agit d’un espace Hilbert muni de produit scalaire < ·, · >, et la norme || · ||2 où < X, Y >= E XY et ||X||2 =
E X 2 = Var X(E X)2 .

Théorème 4.5.4 théorème de projection Si H est un sous espace fermé de L2 pour tout y ∈ L2 , il existe une
unique v.a. ŷ ∈ H t.q.
||y − ŷ|| = min ||Y − H||
H∈H

caractérisé par Ŷ ∈ H et y − ŷ ∈ H⊥ (i.e.: < Y − Ŷ , H >= 0, ∀H ∈ H).

‡ On trouve un résultat similaire Proposition 7.3.3 dans le chapitre 7

§ E (XY |ξ) = XE (Y |ξ) si X est ξ mesurable
 Y

Y − Ŷ

Ŷ
H

Régression affine dans L2

Soit X1 , X2 , · · · , Xn n v.a. Y une v.a. dans L2 , alors la régression de Y sur X1 , · · · , Xn est la projection orthogonale
dans L2 de Y sur H∗ = Vect(X1 , · · · , Xn ) et la régression affine de Y sur X1 , · · · , Xn est la projection orthogonale
dans L2 de Y sur H = Vect(1, X1 , X2 , · · · , Xn ) et on note Ŷ = E L(Y |1, X1 , X2 , · · ·).

Exercice 4.5.5 Soit X1 , · · · , Xn un vecteur, Xi ∈ L2 et H = Vect(1, X1 , · · · , Xn ) espace combinaison linéaire affine

de X1 , · · · , Xn . Soit X = X0 ∈ L2 ,
X
γ = Cov (X, Xi )i=1,···,n , = matrice Cov (Xi , Xj )i,j=1,2···,n

Alors
X̂ = E L(X|1, X1, · · · , Xn ) = a0 + a1 + · · · + an
P−1
avec a = (a1 , · · · , an ) vérifiant a = γ

Solution : On suppose X̂ = a0 + a1 + · · · + an , ai à déterminer. On sait que X̂ ∈ H , et X − X̂ ∈ H⊥ cad

< X − X̂, 1 >= 0 (4.2)

< X − X̂, Xi >= 0, (4.3)

avec 4.2 on a E X = E X̂ = a0 + a1 E X1 + · · · + αn E Xn avec 4.3 on a E XXi = E X̂Xi = a0 E Xi + a1 E X1 Xi + · · · +

αn E Xn Xi donc

a0 = E X − a1 E X 1 − · · · − an E X n
E XXi = E XE Xi − (a1 E X1 + · · · + αn E Xn )E Xi +a1 E X1 Xi + · · · + αn E Xn Xi

∀i = 1, · · · , n nous avons Cov (X, Xi ) = a1 Cov (X1 , Xi ) + α2 Cov (X2 , Xi ) + · · · + an Cov (Xn , Xi )


 Cov (X, X1 ) = a1 Cov (X1 , X1 ) + α2 Cov (X2 , X1 ) + · · · + an Cov (Xn , X1 )

Cov (X, X2 ) = a1 Cov (X1 , X2 ) + α2 Cov (X2 , X2 ) + · · · + an Cov (Xn , X2 )

 ···

Cov (X, Xn ) = a1 Cov (X1 , Xn ) + α2 Cov (X2 , Xn ) + · · · + an Cov (Xn , Xn )
P
γ= a sous forme matricielle, d’où le résultat.

Prévision linéaire
Définition 4.5.6 Soit (Xt )t∈Z un processus L2 , on appelle meilleur prévision linéaire de Xt sachant son passé la
régression linéaire de Xt sur son passé H = Vect(1, Xt−1 , Xt−2 , · · ·)

Définition 4.5.7 processus d’innovation est le processus

ηt = Xt − E L(Xt |1, Xt−1 , Xt−2 , · · ·)

Exemple : Soit Yt ∼ BBG(0, σ 2 ) le processus Xt est définit par

Xt = Yt − ρYt−1 , |ρ| < 1

∞
X
Alors E L(Xt |1, Xt−1 , Xt−2 , · · ·) = − ρi Xt−i est la régression linéaire de Xt sur son passé. et en plus, le processus
i=1
d’innovation est Yt .
 Pn Pn
Solution : Il est évident que Xt est stationnaire. Soit St,n = i=1 ρi Xt−i alors Xt = Yt −ρYt−1 , ajouer St,n = i=1
à deux côtés et passer à la limite

Xt + St,n = Yt − ρYt−1 + ρXt−1 + ρ2 (Xt−2 + · · · + ρn Xt−n

= Yt − ρYt−1 + ρ(Yt−1 − ρYt−2 ) + ρ2 (Yt−2 − ρYt−3 ) + · · · + ρn (Yt−n − ρYt−n )
= Yt + ρn+1 Yt−n−1
∞
X
Xt + ρi Xt−i = Yt car |ρ| < 1
i=1
Pn Pn
or la suite St,n est convergante : ||St,n − St,m || = || i=m+1 ρi Xt−i || ≤ ||Xt || i=m+1 |ρi | −→ 0 quand n → ∞ et m →
∞, c’est une suite de cauchy, donc convergante. Soit St,n → St on a Xt + St = Yt . Maintenant < Xs , Yt >= 0, ∀s < t
et < 1, Yt >= 0 (Yt est un BB!) et St = E L(Xt |1, Xt−1 , · · ·) donc tous les conditions sont réunies pour que Yt soit un
processus d’innovation de (X1 , X2 , · · · , Xn ).

4.6 Prévision des modèles ARIMA, méthode de Box-Jenkins

Soit (X1 , X2 , · · · , XT ) une série stationnaire, on cherche à faire prévision à horizon h. Dans un premier temps, les
modèles AR, MA et ARMA sont supposés mis sous forme canonique, et n’avoir aucune racine unité, et tous les racines
de Φ(z) et Θ(z) sont de modules strictement supérieur à 1. On va considérer succesivement le cas AR(p), MA(q) et
ARMA(p,q)

4.6.1 Prévision de AR(p)

C’est cas le facile à faire. Xt = φ1 Xt−1 + · · · + φp Xt−p + εt , Φ(B)Xt = εt . Donc la prévision optimale à horizon 1,
faite à la date T est
XT∗ +1|T = E L(XT +1 |XT ) = φ1 XT + · · · + φp XT −p
car εt est l’innovation. De la même façon, XT +h = φ1 XT +h−1 + φ2 XT +h−2 + · · · + φn XT +h−n alors XT∗ +h|T =
E L(XT +h |XT ) est donné de façon récursive par
(
∗ φ1 XT∗ +h−1|T + · · · + φh−1 XT∗ +1|T + φh XT + · · · + φp XT +h−p pour h ≤ p
XT +h|T = (4.4)
φ1 XT∗ +h−1|T + · · · + φp XT∗ +h−p|T pour h > p

Exercice 4.6.1 AR(1) moyenne non nulle Calculer la prévision à horizon h du modèle AR(1) : Xt = µ+φXt−1 +εt
(µ considéré connu).
Solution :
XT∗ +1|T = µ + φXT
, XT∗ +2|T = µ + φXT∗ +1|T = µ + φ(µ + φXT ) = µ(1 + φ) + φ2 XT
de façon itérative
XT∗ +h|T = µ + φXT∗ +h−1|T = µ(1 + φ + φ2 + · · · + φh−1 ) + φh XT .
µ
On remarque que lorsque h → ∞, XT∗ +h|T tend vers 1−φ , la vrai moyenne du processus Xt . L’erreur de prévision à
l’horizon h est
εh = XT∗ +h|T − XT +h = XT∗ +h|T − (φXT +h−1 + µ + εT +h ) = · · ·
= XT∗ +h|T − (φh XT + (φh−1 + · · · + φ + 1)µ + εT +h + φεT +h−1 + · · · + φh−1 εT +1 )
= εT +h + φεT +h−1 + · · · + φh−1 εT +1 )

Son MSE est alors

MSE XT∗ +1|h = (1 + φ2 + φ4 + · · · + φ2h−2 )σ 2
où σ 2 est la variance de εt .
Exercice 4.6.2 On consideère le processus suivant,

Xt = 4.0 + 0.8Xt−1 + εt , ε ∼ N (0, 1)

dont les dernières observations ont été

21.287168, 19.002169, 18.664791, 17.550591, 21.724832, 21.8215, 22.59378, 23.025447, 22.023941,
22.795715, 20.833277, 20.608492, 19.533135, 20.927484, 21.572549, 20.819312, 20.63309, 20.998086,
22.245658, 20.709627, 21.132428, 20.935427, 20.701293, 20.16264, 18.862392, 20.378879, 21.843372,
21.366562, 22.460745, 22.286604,

Calculer la prévision à l’horizon 1, 2, · · · , 7. Calculer intervalle de coniance de votre prévision.

 Solution.
La prévision XT∗ +1|T est

XT∗ +1|T = µ + φXT = 4.0 + 0.8 ∗ 22.286604 = 21.829283

XT∗ +2|T = µ + φXT∗ +1|T = 4.0 ∗ 0.8 ∗ 21.829283 = 21.463426
De meme on peut trouver
XT∗ +3|T = 21.170741, XT∗ +4|T = 20.936593, XT∗ +5|T = 20.749274
XT∗ +6|T = 20.599419 XT∗ +7|T = 20.479536

On peut biensur calculer MSE XT∗ +h|T facilement. On suppose ensuite que la variable aléatoire XT∗ +h|T suit une loi
normale de moyenne XT∗ +h|T , et de variance MSE XT∗ +h|T , pour calculer ensuite intervalle de confiance de XT∗ +h|T :
h p p i
XT∗ +h|T − µ1−α/2 σ 1 + φ2 + φ4 + · · · + φ2h−2 , XT∗ +h|T + µ1−α/2 σ 1 + φ2 + φ4 + · · · + φ2h−2

Par exemple, h = 2 nous avons ici σ = 1 XT∗ +2|T = 21.463, avec intervalle de confiance de niveau 0.95
h p p i
21.829 − 1.96σ 1 + 0.82 21.829 + 1.96σ 1 + 0.82 = [21.829 − 2.5121.829 + 2.51]] = [19.319, 24.339]

où α = 0.05 et le niveau de signification est 1 − α = 0.95, donc la valeur critique est µ0.975 = 1.96 Voici le graphe
explicatif

24

23

22

21

20

19

18

17
0 5 10 15 20 25 30 35 40

Figure 4.6: Xt et XT∗ +h|T

µ
Exercice 4.6.3 suite Soit Yt = Xt − 1−φ calculer les prévisions YT∗+h|T

Solution : Facile, après cette transformation, Yt est un AR(1) de moyenne nulle :

µ µ
Yt = Xt − 1−φ = µ + φXt−1 + εt − 1−φ
φµ φµ µ
= µ + φXt−1 − 1−φ + 1−φ + εt + 1−φ
µ φµ µ µ
= φ(Xt−1 − 1−φ ) + ε + µ + 1−φ − 1−φ = φ(Xt−1 − 1−φ ) +ε

D’après le résultat 4.4, nous avons YT∗+h|T = φYT∗+h−1|T = · · · = φh YT . Revenons sur Xt nous avons
µ µ
XT∗ +h|T = YT∗+h|T + 1−φ = φh YT + 1−φ
n
h µ
= φ (XT − 1−φ ) + 1−φ = φ XT + µ( 1−φ
µ h
1−φ )

Comparer avec exercice précédent, on a le même résultat.

Exercice 4.6.4 AR(2) moyenne non nulle On considère le processus suivant :
Xt = 0.5 + 0.8εt−1 − 0.2εt−2 , εt ∼ N (0, σ 2 )
Est-ce qu’il s’agit d’un processus stationnaire? Calculer la moyenne E Xt de ce processus, et les prévisions à l’horizon
1, 2, · · · 5. Sachant les 9 premières observations sont
13.787168, 10.472434, 10.08357, 9.5911267, 14.340547, 14.995847, 15.265148, 15.163372, 13.681251,
Solution. Le calcul est presque identique que le cas, AR(1).
4.6.2 Prévision de MA(q)
Xt = εt + θ1 + · · · + θ1 εt−q = Θ(B)εt . La prévision optimale à horizon 1, faite à la date T est

XT∗ +1|T = E (XT +1 |XT , XT −1 , · · ·) = E L(XT +1 |εT εT −1 , · · ·)

car σ(XT , XT −1 , · · ·) = σ(εT , εT −1 , · · ·), par conséquence

XT∗ +1|T = θ1 εT + · · · + θq εT +1−q + 0

de même ∀h ≥ 1, nous avons le résultat suivant

XT∗ +h|T = θh εT + · · · + θq εT +h−q pour h ≤ q
(4.5)
0 pourh > q

Inconvenient de cette méthode : Il faut connaitre (estimer) les résidus passés, a priori non observables. E (XT +1 |XT , XT −1 , · · ·) =
E L(XT +1 |εT εT −1 , · · ·) . Revenons sur Xt nous avons

utilisation de AR(∞)
On peut écrire Xt = Θ(B)εt par Θ−1 (B)Xt = εt ce qui nous donne
∞
X ∞
X
Xt = ak Xt−k + εt donc Xt+h = ak Xt+h−k + εt+h
k=1 k=1

On peut donc calculer XT∗ +h|T de façon itérative

h−1
X ∞
X
XT∗ +h|T = ak XT∗ +h−k|T + ak XT +h−k|T
k=1 k=h

cette formule n’est pas utilisable, car les (Xt ) ne sont pas observés en pratique pour t < 0. On utlise une autre écriture
∞
X h−1
X ∞
X
XT +h = ak XT +h−k + εT +h = ak XT +h−k + ak XT +h−k + εT +h
k=1 k=1 k=h
P∞
k=1 ak Xt+h−k est une série convergante, le reste de cette série est négligéable. On peut considérer quand T suffisa-
ment grand que
h−1
X T
X +h X∞
XT∗ +h|T = ak XT∗ +h−k|T + ak XT +h−k + ak XT +h−k
k=1 k=h k=T +h+1

le troisième terme étant négligeable, on approxime XT∗ +h|T par X̂T∗ +h|T

h−1
X T
X +h
X̂T∗ +h|T = ak X̂T∗ +h−k|T + ak XT +h−k
k=1 k=h

Cette formule devient calculable

Formule de mise à jour des résultats

Dans les calculs en temps réel, les observations Xt arrive successivement, idée est de comparer l’estimation XT∗ +1|T
faite en T de XT +1 avec la vraie valeur, XT +1 − XT∗ +1|T = εT +1 . L’équation 4.5 peut être modifiée. En se plaçant à
la date T + 1, on peut noter que XT∗ +h|T +1 = X(T∗
+1)+(h−1)|T +1 donc

θh−1 εT +1 + XT∗ +h|T pour h ≤ q + 1
XT∗ +h|T +1
0 pour h > q + 1
donc pour h ≤ q + 1
XT∗ +h|T +1 = XT∗ +h|T + θh−1 (XT +1 − XT∗ +1|T )− (4.6)
Cette formule s’appelle formule de mise à jour.
Exercice 4.6.5 MA(1) moyenne non nulle On considère le processus

Xt = 10.0 + εt − 0.5εt−1 , où ε ∼ N (0, 1)

Les observations ont été 11.066351, 7.3288504, 10.476839, 8.8872297, 14.374981 Calculer la prévision à l’horizon
h = 1, 2, · · · de ce processus.
Solution. On défini d’abort un processus centré Yt = Xt −10 = (1.0663511, −2.6711496, 0.4768388, −1.1127703, 4.3749811)
Pour calculer la prévision d’un modèle MA(q) on cherche d’abord sa forme AR(∞).

εt = Yt + 0.5εt−1 = Yt + 0.5(Yt−1 + 0.5εt−2 ) = Yt + 0.5(Yt−1 + 0.5(Yt−2 + 0.5εt−3 )) = · · ·

X∞ ∞
X
= 0.5i Yt−i = 0.5i (Xt−i − 10.0)
i=0 i=0

On obtient donc une forme AR(∞)

∞
X
Yt = εt + ak Xt−k = εt − 0.5Yt−1 − 0.52 Yt−2 − 0.53 Yt−3 − · · ·
k=1

La formule qui permet de calculer la prévision à l’horizon h est

h−1
X T
X +h
YT∗+h|T = ak YT∗+h−k|T + ak YT +h−k
k=1 k=h

Ici ak = −0.5k , et T = 4, lorsque h = 1 nous avons

4+1
X
∗
Y4+1|4 = ak Y4+1−k = −0.5Y4 − 0.52 Y3 − 0.53 Y3 − 0.54 Y1 − 0.55 Y0 = −1.8352794
k=1

Lorsque h=2, nous avons

Ph−1 P +h P
∗
Y4+2|4 = ak YT∗+h−k|T + Tk=h
k=1
∗
ak YT +h−k = a1 Y5|4 + 6k=2 ak Y6−k
∗
= −0.5Y5|4 − 0.52 Y4 − · · · − 0.56 Y0 = 0.0

En fait, il est tres facile de montrer que avec cette formule nous avons toujours

YT +h|T = 0.0, si h ≥ 2

Pour calculer XT∗ +h|h , il suffit alors de utiliser XT∗ +h|T = YT∗+h|T + 10, ce qui donne
∗ ∗ ∗
X4+1|4 = 8.1642, X4+2|4 = X4+3|4 = 10.0

Attention, ici le vecteur d’observation est X4 = (X0 , X1 , · · · , X4 ).

4.6.3 Prévision de ARMA(p,q)

Le processus Φ(B)Xt = Θ(B)εt :
p
X q
X
Xt = φi Xt−i + εt + θj εt−j
i=1 j=1
p
X q
X
Xt+h = φi Xt+h−i + εt+h + θj εt+h−j
i=1 j=1

On sait XT∗ +h|T = E L(XT +h |XT , XT −1 , · · ·) = E L(XT +h |εT , εT −1 , · · ·). Quand h > q la prévision est similaire que le
cas AR(p) (les termes qui contiennent ε sont tous indépedant de moyenne nulle par rapport à (XT , XT −1 , · · ·) quand
h > q). (
∗ φ1 XT∗ +h−1|T + · · · + φh−1 XT∗ +1|T + φh XT + · · · + φp XT +h−p pour h ≤ p
XT +h|T =
φ1 XT∗ +h−1|T + · · · + φp XT∗ +h−p|T pour h > p

mais quand t ≤ q, on ne peut pas enlever ε, et elle sont non observable....

utilisation de AR(∞) pour ARMA(p,q)

Tout comme le cas de MA(q). on se serve de la forme AR(∞) pour contruire un algorithme calculable. L’équation
Φ(B)Xt = Θ(B)εt peut toujours être écrite sous la forme (on suppose le processus est stationnaire et inversible)
εt = Θ−1 (B)Φ(B)Xt . On obtient donc
∞
X ∞
X
Xt = ak Xt−k + εt =⇒ Xt+h = ak Xt+h−k + εt+h
k=1 k=1
On utilise la même technique que dans le cas MA(q), en écrivant
∞
X h−1
X T
X +h ∞
X
XT +h = ak XT +h−k + εt+h = ak XT +h−k + ak XT +h−k + ak XT +h−k + εt+h
k=1 k=1 k=h k=T +h+1
P∞
le terme k=T +h+1 ak XT +h−k est négligeable quand T est grande. D’où la formule itérative :

h−1
X T
X +h
X̂T∗ +h|T = ak X̂T∗ +h−k|T + ak XT +h−k
k=1 k=h

Ici, XT∗ +h|T est remplacée par X̂T∗ +h|T , car c’est une approximation.