Séries temporelles multivariées : quand une série en regarde une autre
Neuf chapitres à considérer les séries « une à la fois », puis la question inévitable : que se passe-t-il quand elles réagissent les unes aux autres ? Corrélations croisées, processus VAR, condition de stationnarité, et le piège du préblanchiment. Avec un laboratoire VAR(1) interactif et du code Python exécutable dans la page.
- Pourquoi c'est important
- Le vocabulaire vectoriel
- Stationnarité et corrélations croisées
- Le bruit blanc multivarié
- Le processus VAR(p)
- Laboratoire : le VAR(1) de l'exemple 10.2
- Un mot sur le spectre croisé
- Estimation et choix de l'ordre
- Tester une association : la procédure de Haugh
- Le cas non stationnaire
- Limites et pièges
- Références
Pourquoi c'est important
Voici deux séries temporelles. Leur corrélation vaut \(-0{,}02\) : rigoureusement rien. Vous les superposez, vous calculez un coefficient de Pearson, vous concluez qu'elles n'ont aucun rapport, et vous passez à autre chose.
Vous venez de manquer une dépendance de \(-0{,}76\). Elle était là, entière, mais décalée de trois pas de temps. La première série ne réagit pas à la seconde au même instant — elle y réagit trois périodes plus tard. Aucun coefficient de corrélation ne peut voir ça, parce qu'un coefficient de corrélation ne regarde qu'un seul alignement : celui de maintenant.
C'est ce que le cadre vectoriel rend visible. Là où l'analyse univariée traite chaque série « une à la fois », le modèle VAR laisse chacune écouter le passé des autres. Trois taux d'intérêt qui réagissent ensemble aux annonces de la banque centrale, deux capteurs qui subissent le même choc thermique, une variable économique qui en entraîne une autre avec un délai : dans tous ces cas, l'information utile est précisément celle que l'approche série-par-série jette à la poubelle.
Et il y a un piège symétrique, plus vicieux encore. Deux séries lisses et sans aucun lien paraissent corrélées, à presque tous les décalages. Vous verrez ici pourquoi — et le remède, le préblanchiment, tient en quatre lignes. Deux laboratoires vous laissent déformer un modèle VAR à la main, regarder ses racines franchir le cercle unité, et faire s'effondrer des corrélations fantômes d'un clic.
Le vocabulaire vectoriel
Un processus stochastique multivarié \(\{\mathbf{X}_t;\, t \in T\}\) est simplement une collection de variables aléatoires vectorielles. Chaque composante \(\{X_{ti}\}\) est une série temporelle univariée ordinaire, et on les empile :
Pour les données de taux d'intérêt américains mensuels de janvier 1991 à avril 2010, \(X_{t1}\) désigne le taux des certificats de dépôt à six mois, \(X_{t2}\) celui des obligations Aaa de Moody's, et \(X_{t3}\) celui des hypothèques conventionnelles à trente ans. Les trois réalisations univariées forment une seule réalisation d'un processus trivarié.
La moyenne devient un vecteur, la covariance une matrice :
L'élément \(\gamma_{ij}(t_1,t_2)\) de cette matrice est la covariance entre \(X_{t_1 i}\) et \(X_{t_2 j}\) : deux séries différentes, à deux instants différents. C'est exactement l'information qu'une analyse univariée jette.
Stationnarité et corrélations croisées
Comme dans le cas univarié, on impose que les moments ne dépendent pas de l'origine du temps. Le processus est stationnaire en covariance si \(E[\mathbf{X}_t] = \boldsymbol\mu\) est constant, si \(\Gamma(0)\) existe et est constante, et si \(\Gamma(t_1,t_2)\) ne dépend que de l'écart \(t_2 - t_1\). On écrit alors :
Quand \(i \neq j\), \(\gamma_{ij}(k)\) s'appelle la covariance croisée au retard \(k\). Les éléments diagonaux de \(\Gamma(0)\) sont les variances des composantes. La corrélation croisée normalise :
Deux propriétés méritent d'être mémorisées, parce qu'elles ne ressemblent pas à leurs analogues univariées. Dans le cas univarié, \(\gamma(k) = \gamma(-k)\) : l'autocovariance est paire. Ici, ce n'est pas vrai pour les termes croisés. On a seulement :
Notez le changement d'ordre des indices. C'est précisément ce qui permet de distinguer « \(X_1\) précède \(X_2\) » de « \(X_2\) précède \(X_1\) » — l'asymétrie porte toute l'information directionnelle. Les valeurs \(\rho_{ij}(0)\) mesurent quant à elles la relation contemporaine, et \(\rho(0)\) porte des 1 sur sa diagonale.
ccf(xi, xj) calcule la seconde ; pour obtenir la première, il faut appeler
ccf(xj, xi). Vérifiez toujours quelle convention emploie la source avant
d'interpréter un pic à un retard donné — sinon vous conclurez que la cause suit l'effet.
Le bruit blanc multivarié
La série \(\{\mathbf{a}_t\}\) est un bruit blanc multivarié si \(E[\mathbf{a}_t] = \mathbf{0}\) et \(\boldsymbol\rho_a(k) = \mathbf{0}\) pour tout \(k \neq 0\). Rien n'est imposé sur \(\boldsymbol\rho_a(0)\). Autrement dit, les composantes peuvent être corrélées de façon contemporaine : \(\rho_{ij}(0) \neq 0\) est parfaitement permis. On note \(\mathbf{a}_t \sim \text{MWN}(\mathbf{0}, \Gamma_a)\) avec \(\Gamma(0) = \Gamma_a\).
Cela paraît contradictoire au premier abord — du bruit blanc corrélé ? Mais les deux notions portent sur des axes différents. Le bruit blanc interdit toute structure le long de l'axe du temps ; il n'interdit rien en travers, entre composantes au même instant. Deux capteurs qui subissent le même choc thermique à chaque instant produisent exactement cela. L'exemple canonique du chapitre est :
Chaque réalisation prise isolément semble parfaitement non corrélée le long de l'axe du temps ; superposées, elles montent et descendent ensemble. Cette matrice \(\Gamma_a\) est celle qui alimente le VAR(1) du laboratoire ci-dessous.
Le processus VAR(p)
Le modèle central du chapitre. Chaque composante est régressée sur les \(p\) valeurs passées de toutes les composantes, y compris elle-même. En dimension deux et à l'ordre un :
Ce qui se condense en \(\mathbf{X}_t = \Phi_1 \mathbf{X}_{t-1} + \mathbf{a}_t\). Le cas général :
où \(\mathbf{X}_t\) est \(m \times 1\), chaque \(\Phi_k\) est une matrice réelle \(m \times m\), \(\mathbf{a}_t\) est un bruit blanc multivarié de covariance \(\Gamma_a\), et \(\boldsymbol\alpha = (I - \Phi_1 - \cdots - \Phi_p)\boldsymbol\mu\).
En notation opérateur, \(\Phi(B)(\mathbf{X}_t - \boldsymbol\mu) = \mathbf{a}_t\) avec \(\Phi(B) = I - \Phi_1 B - \cdots - \Phi_p B^p\). Le premier vrai résultat concerne la stationnarité, et c'est une condition purement algébrique :
Le processus est stationnaire et causal si et seulement si les racines de l'équation déterminantale
\[|\Phi(z)| = |I - \Phi_1 z - \cdots - \Phi_p z^p| = 0\]se situent toutes à l'extérieur du cercle unité.
C'est l'analogue exact de la condition univariée sur \(\phi(z) = 0\). Pour un VAR(1), \(|\Phi(z)| = |I - \Phi_1 z|\), et les racines sont les inverses des valeurs propres de \(\Phi_1\) : la condition revient à exiger que toutes les valeurs propres de \(\Phi_1\) soient de module strictement inférieur à 1. Le laboratoire ci-dessous rend cette condition visible.
Sous stationnarité, le processus admet une représentation linéaire générale \(\mathbf{X}_t - \boldsymbol\mu = \Psi(B)\mathbf{a}_t = \sum_{j \geq 0}\Psi_j \mathbf{a}_{t-j}\) avec \(\Psi(B) = \Phi^{-1}(B)\), et les autocovariances satisfont les équations de Yule–Walker multivariées :
Et pour \(k = 0\), la relation \(\Gamma(0) - \Phi_1\Gamma(0)\Phi_1' = \Gamma_a\) — une équation de Lyapunov discrète — détermine la covariance stationnaire. Elle se résout en vectorisant : \(\operatorname{vec}\Gamma(0) = (I - \Phi_1 \otimes \Phi_1)^{-1}\operatorname{vec}\Gamma_a\).
Laboratoire : le VAR(1) de l'exemple 10.2
L'exemple traité dans le chapitre est un VAR(1) bivarié dont les paramètres se reconstituent entièrement à partir des covariances publiées. En appliquant \(\Phi_1' = \Gamma(0)^{-1}\Gamma(1)\) aux matrices du livre, on retrouve :
Ces valeurs redonnent exactement l'équation caractéristique \(1 - 1{,}6z + 0{,}81z^2 = 0\) mentionnée dans le texte, dont les racines sont complexes de module \(1/\sqrt{0{,}81} = 1{,}11\) : à l'extérieur du cercle unité, donc stationnaire, et complexes, d'où le comportement pseudo-cyclique de période environ 13 visible sur les réalisations. Elles redonnent aussi, par l'équation de Lyapunov, la covariance stationnaire du livre :
Voilà le fait le plus instructif de tout le chapitre. La corrélation contemporaine entre les deux séries vaut \(-0{,}02\) : rigoureusement rien. Un analyste qui superpose les deux courbes et calcule un coefficient de Pearson conclut à l'indépendance. Pourtant \(\rho_{12}(3) = -0{,}76\) : décalez l'une de trois périodes, et la relation apparaît, massive. Une relation forte peut être totalement invisible au retard zéro.
Un mot sur le spectre croisé
Le domaine fréquentiel se généralise de la même manière. Le spectre devient une matrice \(m \times m\) pour chaque fréquence :
Les éléments diagonaux \(P_{jj}(f)\) sont les spectres univariés habituels. Les éléments hors diagonale \(P_{jh}(f)\), appelés spectres croisés, sont complexes. On les décompose en \(P_{jh}(f) = c_{jh}(f) - i\,q_{jh}(f)\), où \(c_{jh}\) est le cospectre et \(q_{jh}\) le spectre de quadrature. En forme polaire \(P_{jh}(f) = \alpha_{jh}(f)e^{i\varphi_{jh}(f)}\), le module est le spectre d'amplitude et \(\varphi_{jh}(f) = \tan^{-1}(-q_{jh}(f)/c_{jh}(f))\) le spectre de phase — qui indique le déphasage entre les deux séries à chaque fréquence. Enfin la cohérence carrée :
Elle s'interprète comme un \(R^2\) fréquence par fréquence : la part de la variabilité de \(X_j\) à la fréquence \(f\) explicable linéairement par \(X_h\) à cette même fréquence. Deux séries peuvent être fortement cohérentes aux basses fréquences et pas du tout aux hautes — un diagnostic que la corrélation croisée seule ne donne pas.
Estimation et choix de l'ordre
Les estimateurs sont les analogues directs du cas univarié. La moyenne est estimée par le vecteur des moyennes empiriques \(\bar{\mathbf{X}}\), et les covariances par :
Notez le \(1/n\) et non \(1/(n-k)\) : estimateur biaisé mais de variance moindre, et qui garantit une matrice définie non négative. Pour \(k = 0\), \(\hat\rho_{ij}(0)\) est simplement le coefficient de Pearson usuel calculé sur les \(n\) observations.
Pour les paramètres du VAR, deux voies : les équations de Yule–Walker multivariées, ou les moindres carrés équation par équation. Elles donnent des résultats très proches en pratique. L'ordre \(p\) se choisit par AIC — indispensable, puisqu'un VAR(\(p\)) en dimension \(m\) compte \(m^2 p\) coefficients : un VAR(3) trivarié en compte déjà 27.
# ─── R : base + package vars ─── # la fonction ar() de base calcule les estimations de Yule-Walker pour un VAR ; # attention, les methodes mle et burg sont univariees uniquement donnees <- cbind(x1, x2) ar(donnees, order.max = 8, aic = TRUE, method = 'yw') ar(donnees, order.max = 8, aic = TRUE, method = 'ols') library(vars) ajust <- VAR(donnees, p = 1, type = c("const"), ic = c("AIC")) serial.test(ajust, lags.pt = 24, type = c("PT.asymptotic")) # Ljung-Box multivarie predict(ajust, n.ahead = 25) # prevision + limites a 95 % par defaut
Le code ci-dessous s'exécute réellement dans votre navigateur (Python compilé en WebAssembly). Modifiez les valeurs et relancez : le premier lancement prend quelques secondes le temps de charger Python.
Tester une association : la procédure de Haugh
Question fréquente : ces deux séries sont-elles liées ? On aimerait tester \(H_0 : \rho_{12}(k) = 0\) en comparant \(\hat\rho_{12}(k)\) à des bandes du type \(\pm 2/\sqrt{n}\). Mais la variance asymptotique de \(\hat\rho_{12}(k)\) vaut :
Regardez bien cette formule. Elle fait intervenir les autocorrélations des deux séries. Si les deux sont fortement autocorrélées, la somme est grande, et la variance réelle dépasse de loin \(1/n\). Les bandes \(\pm 2/\sqrt{n}\) sont alors beaucoup trop étroites : on voit des relations partout. La variance ne se réduit à \(1/n\) que si au moins une des deux séries est un bruit blanc.
D'où la procédure recommandée par Haugh (1976), en quatre étapes :
| Étape | Opération |
|---|---|
| 1 | Ajuster un ARMA(\(p_i, q_i\)) à chacune des deux séries — les ordres peuvent différer |
| 2 | Calculer les séries de résidus \(\hat a_{t1}\) et \(\hat a_{t2}\) (c'est le préblanchiment) |
| 3 | Calculer les corrélations croisées \(\hat\rho^{(a)}_{12}(k)\) entre les deux séries de résidus |
| 4 | Déclarer \(\hat\rho^{(a)}_{12}(k)\) significative au seuil \(\alpha = 0{,}05\) si \(|\hat\rho^{(a)}_{12}(k)| > 2/\sqrt{n}\) |
Les résidus sont du bruit blanc par construction, donc la variance asymptotique redevient \(1/n\) et les bandes usuelles retrouvent leur validité. Comme toujours avec des tests répétés, le seuil s'applique séparément à chaque \(k\) : sur 41 retards testés, deux dépassements fortuits sont attendus sous \(H_0\).
Le cas non stationnaire
Tout ce qui précède suppose la stationnarité. Quand elle manque, on généralise le modèle ARUMA au cas multivarié en insérant un opérateur diagonal \(\lambda(B)\) :
où chaque \(\lambda_j(z) = 0\) a toutes ses racines sur le cercle unité, tandis que celles de \(|\Phi(z)| = 0\) et \(|\Theta(z)| = 0\) restent à l'extérieur. Le cas qui a reçu le plus d'attention est celui où \(\lambda(B) = \operatorname{diag}\big((1-B)^{d_1}, \dots, (1-B)^{d_m}\big)\) — le pendant multivarié exact du modèle ARIMA(\(p,d,q\)). Si la composante \(j\) est déjà stationnaire, \(d_j = 0\) et aucune différenciation ne lui est appliquée.
Limites et pièges
| Hypothèse ou piège | Conséquence → parade |
|---|---|
| Nombre de paramètres en \(m^2p\) | Surajustement rapide dès \(m \geq 4\) → AIC/BIC, VAR bayésien, ou restrictions de nullité |
| Convention de la covariance croisée | Signe du retard inversé selon la source → vérifier, et préférer ccf(xj, xi) en R pour la convention de ce texte |
| Lecture d'une CCF brute | Corrélations fantômes entre séries autocorrélées → préblanchir systématiquement |
| Corrélation contemporaine nulle | Ne signifie pas indépendance : voir \(\rho_{12}(0) = -0{,}02\) et \(\rho_{12}(3) = -0{,}76\) dans le même modèle |
| Composantes non stationnaires | Régression fallacieuse → tester les racines unitaires, puis cointégration avant de différencier |
| Bruit blanc à corrélation contemporaine | \(\Gamma_a\) non diagonale complique l'interprétation causale → décomposition de Cholesky et ordre des variables, qui n'est pas neutre |
| VARMA au-delà du VAR | Problèmes d'identifiabilité sérieux → en pratique, on reste presque toujours sur du VAR pur |
Références
- Woodward, W. A., Gray, H. L. & Elliott, A. C. (2017). Applied Time Series Analysis with R, 2e éd., CRC Press — chapitre 10, la source directe de cet article (définitions, exemple 10.2, procédure de test).
- Haugh, L. D. (1976). Checking the independence of two covariance-stationary time series: a univariate residual cross-correlation approach, JASA 71 — la procédure de préblanchiment.
- Brockwell, P. J. & Davis, R. A. (1991, 2002). Time Series: Theory and Methods — la théorie asymptotique des estimateurs multivariés.
- Reinsel, G. C. (1997). Elements of Multivariate Time Series Analysis, 2e éd., Springer.
- Box, G. E. P., Jenkins, G. M. & Reinsel, G. C. (2008). Time Series Analysis: Forecasting and Control, 4e éd. — la convention de covariance croisée retenue ici.
- Lütkepohl, H. (2005). New Introduction to Multiple Time Series Analysis, Springer — la référence exhaustive sur les VAR.
- Pfaff, B. (2008). VAR, SVAR and SVEC models: implementation within R package vars, Journal of Statistical Software 27(4).
Les simulateurs de cette page sont écrits en JavaScript/Canvas sans dépendance externe ; les cellules Python s'exécutent via Pyodide. Le code est disponible sur GitHub.