\( \newcommand\R{\mathbb{R}} \newcommand\Z{\mathbb{Z}} \newcommand\LL{\mathbb{L}} \newcommand{\E}[1]{\mathbb{E}\left[#1\right]} \newcommand{\V}[1]{\mathbb{V}\left[#1\right]} \newcommand{\ps}[2]{\left\langle #1 \,,\, #2 \right\rangle} \newcommand\1{\mathbb{1}} \newcommand\N{\mathbb{N}} \newcommand\Norm{\mathcal{N}} \newcommand{\transp}[1]{{}^t\!#1} \newcommand\ud{\,\mathrm{d}} \DeclareMathOperator*{\argmax}{argmax} \DeclareMathOperator*{\argmin}{argmin} \DeclareMathOperator{\e}{e} \DeclareMathOperator{\Cov}{Cov} \DeclareMathOperator{\Determinant}{det} \newcommand{\determinant}[1]{\Determinant\left(#1\right)} \)

Chapitre 5 Régression non paramétrique et régression polynomiale locale

Comme notamment montré par Loader (1999), la régression locale est un cas particulier de la régression non paramétrique. Supposons que l’on ait un ensemble de points \((x_i,y_i)_{1\leq i\leq n}\). La régression non paramétrique consiste à supposer qu’il existe une fonction \(\mu\), à estimer, telle que \(y_i=\mu(x_i)+\varepsilon_i\) avec \(\varepsilon_i\) un terme d’erreur. D’après le théorème de Taylor, pour tout point \(x_0\), si \(\mu\) est différentiable \(d\) fois, alors : \[ \forall x \::\:\mu(x) = \mu(x_0) + \mu'(x_0)(x-x_0)+\dots + \frac{\mu^{(d)}(x_0)}{d!}(x-x_0)^d+R_d(x), \] où \(R_d\) est un terme résiduel négligeable au voisinage de \(x_0\). Dans un voisinage \(h(x_0)\) autour de \(x_0\), \(\mu\) peut être approchée par un polynôme de degré \(d\). La quantité \(h(x_0)\) est appelée fenêtre (bandwidth). Si \(\varepsilon_i\) est un bruit blanc, on peut donc estimer par les moindres carrés \(\mu(x_0)\) en utilisant les observations qui sont dans \(\left[x_0-h(x_0),x_0+h(x_0)\right]\).

5.1 Régression polynomiale : approche de Proietti et Luati

5.1.1 Filtres symétriques

Reprenons maintenant les notations de Proietti et Luati (2008) : supposons que notre série temporelle \(y_t\) peut être décomposée en \[ y_t=\mu_t+\varepsilon_t, \] où \(\mu_t\) est la tendance et \(\varepsilon_{t}\overset{i.i.d}{\sim}\mathcal{N}(0,\sigma^{2})\) est le bruit¹⁵. La tendance \(\mu_t\) est localement approchée par un polynôme de degré \(d\), de sorte que dans un voisinage \(h\) de \(t\) \(\mu_t\simeq m_{t}\) avec : \[ \forall j\in\left\llbracket -h,h\right\rrbracket :\: y_{t+j}=m_{t+j}+\varepsilon_{t+j},\quad m_{t+j}=\sum_{i=0}^{d}\beta_{i}j^{i}. \] Le problème d’extraction de la tendance est équivalent à l’estimation de \(m_t=\beta_0\). En notation matricielle : \[ \underbrace{\begin{pmatrix}y_{t-h}\\ y_{t-(h-1)}\\ \vdots\\ y_{t}\\ \vdots\\ y_{t+(h-1)}\\ y_{t+h} \end{pmatrix}}_{y}=\underbrace{\begin{pmatrix}1 & -h & h^{2} & \cdots & (-h)^{d}\\ 1 & -(h-1) & (h-1)^{2} & \cdots & (-(h-1))^{d}\\ \vdots & \vdots & \vdots & \cdots & \vdots\\ 1 & 0 & 0 & \cdots & 0\\ \vdots & \vdots & \vdots & \cdots & \vdots\\ 1 & h-1 & (h-1)^{2} & \cdots & (h-1)^{d}\\ 1 & h & h^{2} & \cdots & h^{d} \end{pmatrix}}_{X}\underbrace{\begin{pmatrix}\beta_{0}\\ \beta_{1}\\ \vdots\\ \vdots\\ \vdots\\ \vdots\\ \beta_{d} \end{pmatrix}}_{\beta}+\underbrace{\begin{pmatrix}\varepsilon_{t-h}\\ \varepsilon_{t-(h-1)}\\ \vdots\\ \varepsilon_{t}\\ \vdots\\ \varepsilon_{t+(h-1)}\\ \varepsilon_{t+h} \end{pmatrix}}_{\varepsilon} \]

Pour estimer \(\beta\) il faut \(H\geq d+1\) et l’estimation est faite par moindres carrés pondérés — weighted least squares (WLS) —, ce qui revient à minimiser la fonction objectif suivante : \[ S(\hat{\beta}_{0},\dots,\hat{\beta}_{d})=\sum_{j=-h}^{h}\kappa_{j}(y_{t+j}-\hat{\beta}_{0}-\hat{\beta}_{1}j-\dots-\hat{\beta}_{d}j^{d})^{2} \] où \(\kappa_j\) est un ensemble de poids appelés noyaux (kernel). On a \(\kappa_j\geq 0:\kappa_{-j}=\kappa_j\), et en notant \(K=diag(\kappa_{-h},\dots,\kappa_{h})\), l’estimateur \(\beta\) peut s’écrire \(\hat{\beta}=(X'KX)^{1}X'Ky\). Avec \(e_{1}=\begin{pmatrix}1&0&\cdots&0\end{pmatrix}'\), l’estimateur de la tendance peut donc s’écrire : \[ \hat{m}_{t}=e_{1}\hat{\beta}=\theta'y=\sum_{j=-h}^{h}\theta_{j}y_{t-j}\text{ avec }\theta=KX(X'KX)^{-1}e_{1} \] En somme, l’estimation de la tendance \(\hat{m}_{t}\) est obtenue en appliquant une moyenne mobile symétrique \(\theta\) à \(y_t\)¹⁶. De plus, \(X'\theta=e_{1}\) donc : \[ \sum_{j=-h}^{h}\theta_{j}=1,\quad\forall r\in\left\llbracket 1,d\right\rrbracket :\sum_{j=-h}^{h}j^{r}\theta_{j}=0. \] Ainsi, la moyenne mobile \(\theta\) préserve les polynômes de degré \(d\).

Concernant le choix des paramètres, l’idée générale qui prévaut est que le choix entre ces différents noyaux est secondaire¹⁷ et qu’il vaut mieux se concentrer sur deux autres paramètres :

• le degré du polynôme \(d\) : s’il est trop petit on risque d’avoir des estimations biaisées de la tendance-cycle et s’il est trop grand on risque d’avoir une trop grande variance dans les estimations (du fait d’un sur-ajustement) ;

• le nombre de voisins \(H=2h+1\) (ou la fenêtre \(h\)) : s’il est trop petit alors trop peu de données seront utilisées pour les estimations (ce qui conduira à une grande variance dans les estimations) et s’il est trop grand alors l’approximation polynomiale sera vraisemblablement fausse ce qui conduira à avoir des estimations biaisées.

5.1.2 Les différents noyaux

Dans les problèmes d’extraction du signal, les observations sont généralement pondérées par rapport à leur distance à la date \(t\) : pour estimer la tendance-cycle à la date \(t\), on accorde généralement plus d’importance aux observations qui sont proches de \(t\).

Dans le cas continu, un noyau \(K\) est une fonction positive, paire et intégrable telle que \(\int_{-\infty}^{+\infty}\kappa(u) \ud u=1\) et \(\kappa(u)=\kappa(-u)\). Dans le cas discret, un noyau est un ensemble de poids \(\kappa_j\), \(j=0,\pm1,\dots,\pm h\) avec \(\kappa_j \geq0\) et \(\kappa_j=\kappa_{-j}\).

Une classe importante de noyaux est celle des noyaux Beta. Dans le cas discret, à un facteur multiplicatif près (de sorte que \(\sum_{j=-h}^h\kappa_j=1\)) : \[ \kappa_j = \left( 1- \left\lvert \frac j {h+1} \right\lvert^r \right)^s,\quad\text{avec }r>0,s\geq 0 \] Cette classe englobe la majorité des noyaux présentés dans cette étude, à l’exception des noyaux d’Henderson, trapézoïdal et gaussien. Les principaux noyaux (qui sont également implémentés dans rjdfilters) sont :

• \(r=1,s=0\) noyau uniforme : \[\kappa_j^U=1\]

• \(r=s=1\) noyau triangulaire : \[\kappa_j^T=\left( 1- \left\lvert \frac j {h+1} \right\lvert \right)\]

• \(r=2,s=1\) noyau d’Epanechnikov (ou parabolique) : \[\kappa_j^E=\left( 1- \left\lvert \frac j {h+1} \right\lvert^2 \right)\]

• \(r=s=2\) noyau quadratique (biweight) : \[\kappa_j^{BW}=\left( 1- \left\lvert \frac j {h+1} \right\lvert^2 \right)^2\]

• \(r = 2, s = 3\) noyau cubique (triweight) : \[\kappa_j^{TW}=\left( 1- \left\lvert \frac j {h+1} \right\lvert^2 \right)^3\]

• \(r = s = 3\) noyau tricube : \[\kappa_j^{TC}=\left( 1- \left\lvert \frac j {h+1} \right\lvert^3 \right)^3\]

• noyau d’Henderson (voir partie 5.1.3 pour plus de détails) : \[ \kappa_{j}=\left[1-\frac{j^2}{(h+1)^2}\right] \left[1-\frac{j^2}{(h+2)^2}\right] \left[1-\frac{j^2}{(h+3)^2}\right] \]

• noyau trapézoïdal : \[ \kappa_j^{TP}= \begin{cases} \frac{1}{3(2h-1)} & \text{ if }j=\pm h \\ \frac{2}{3(2h-1)} & \text{ if }j=\pm (h-1)\\ \frac{1}{2h-1}& \text{ otherwise} \end{cases} \]

• noyau gaussien¹⁸: \[ \kappa_j^G=\exp\left( -\frac{ j^2 }{ 2\sigma^2h^2 }\right) \]

Les noyaux d’Henderson, trapézoïdal et gaussien sont particuliers :

• Les fonctions noyau d’Henderson et trapézoïdal changent avec la fenêtre (les autres dépendent uniquement du rapport \(j/h+1\)).

• Pour les noyaux trapézoïdal et gaussien, d’autres définitions pourraient être utilisées et sont donc définis arbitrairement.
  Le noyau trapézoïdal est implémenté dans rjdfilters car il permet de calculer les moyennes mobiles utilisées dans l’algorithme X-13ARIMA pour l’extraction des composantes saisonnières. Il n’est pas adapté dans le cas de l’extraction de la tendance-cycle.

5.1.3 Quelques filtres symétriques particuliers

Lorsque \(p=0\) (ajustement local par une constante) on obtient l’estimateur de Nadaraya-Watson (ou l’estimateur par noyaux).

Avec le noyau uniforme on obtient le filtre de Macaulay et al. (1931). Lorsque \(p=0\) ou \(p=1\), on retrouve la moyenne arithmétique : \(w_j=w=\frac{1}{2h+1}\).

Le noyau d’Epanechnikov est souvent recommandé comme le noyau optimal car il minimise l’erreur quadratique moyenne de l’estimation par polynômes locaux.

Le Loess, locally estimated scatterplot smoothing (utilisé dans la méthode STL), est une régression locale pondérée qui utilise le noyau tricube.

Le filtre d’Henderson est un cas particulier de l’approximation locale cubique (\(p=3\)), couramment utilisée pour l’extraction de la tendance-cycle (c’est par exemple le filtre utilisé dans le logiciel de désaisonnalisation X-13ARIMA). Pour une fenêtre fixée, Henderson a trouvé le noyau qui donnait l’estimation la plus lisse de la tendance. Il montre l’équivalence entre les trois problèmes suivants :

1. minimiser la variance de la différence d’ordre trois de la série lissée par l’application d’une moyenne mobile ;
   
2. minimiser la somme du carré de la différence d’ordre trois des coefficients du filtre, c’est le critère de lissage (smoothness) : \(S=\sum_j(\nabla^{3}\theta_{j})^{2}\) ;
   
3. estimer une tendance localement cubique par les moindres carrés pondérés, où les poids sont choisis de sorte à minimiser la smoothness (cela conduit au noyau présenté dans la section 5.1.2).

Le filtre d’Henderson étant couramment utilisé pour l’extraction de la tendance-cycle, nous nous intéresserons uniquement aux filtres issus du noyau d’Henderson.

5.1.4 Filtres asymétriques

Comme mentionné dans la partie 2.3.1, pour l’estimation en temps réel, plusieurs approches peuvent être utilisées :

1. Construire un filtre asymétrique par approximation polynomiale locale sur les observations disponibles (\(y_{t}\) pour \(t\in\left\llbracket n-h,n\right\rrbracket\)).

2. Appliquer les filtres symétriques sur les séries prolongées par prévision \(\hat{y}_{n+l\mid n},l\in\left\llbracket 1,h\right\rrbracket\).

3. Construire des filtres asymétriques qui minimisent l’erreur quadratique moyenne de révision sous des contraintes de reproduction de tendances polynomiales.

Proietti et Luati (2008) montrent que les deux premières approches sont équivalentes lorsque les prévisions sont faites par extrapolation polynomiale de degré \(d\). Elles sont également équivalentes à la troisième approche sous les mêmes contraintes que celles du filtre symétrique. Cette méthode est appelée direct asymmetric filter (DAF). C’est cette méthode qui est utilisée pour l’estimation en temps réel dans la méthode de désaisonnalisation STL (Seasonal-Trend decomposition based on Loess, voir R. B. Cleveland et al. (1990)). Même si les estimations des filtres DAF sont sans biais, c’est au coût d’une plus grande variance dans les estimations.

Pour résoudre le problème de la variance des estimations des filtres temps réel, Proietti et Luati (2008) proposent une méthode générale pour construire les filtres asymétriques qui permet de faire un compromis biais-variance. Il s’agit d’une généralisation des filtres asymétriques de Musgrave (1964) (utilisés dans l’algorithme de désaisonnalisation X-13ARIMA).

On modélise ici la série en entrée par : \[\begin{equation} y=U\gamma+Z\delta+\varepsilon,\quad \varepsilon\sim\mathcal{N}(0,D) \tag{5.1} \end{equation}\] où \([U,Z]\) est de rang plein et forme un sous-ensemble des colonnes de \(X\). L’objectif est de trouver un filtre \(v\) qui minimisent l’erreur quadratique moyenne de révision (au filtre symétrique \(\theta\)) sous certaines contraintes. Ces contraintes sont représentées par la matrice \(U=\begin{pmatrix}U_{p}'&U_{f}'\end{pmatrix}'\) : \(U_p'v=U'\theta\) (avec \(U_p\) la matrice \((h+q+1)\times (d+1)\) qui contient les observations de la matrice \(U\) connues lors de l’estimation par le filtre asymétrique). Le problème est équivalent à trouver \(v\) qui minimise : \[\begin{equation} \varphi(v)= \underbrace{ \underbrace{(v-\theta_{p})'D_{p}(v-\theta_{p})+ \theta_{f}'D_{f}\theta_{f}}_\text{variance de l'erreur de révision}+ \underbrace{[\delta'(Z_{p}'v-Z'\theta)]^{2}}_{biais^2} }_\text{Erreur quadratique moyenne de révision}+ \underbrace{2l'(U_{p}'v-U'\theta)}_{\text{contraintes}} \tag{5.2} \end{equation}\] où \(l\) est le vecteur des multiplicateurs de Lagrange.

Lorsque \(U=X\), la contrainte équivaut à préserver les polynômes de degré \(d\) : on retrouve les filtres directs asymétriques (DAF) lorsque \(D=K^{-1}\).

Lorsque \(U=\begin{pmatrix}1&\cdots&1\end{pmatrix}'\), \(Z=\begin{pmatrix}-h&\cdots&+h\end{pmatrix}'\), \(\delta=\delta_1\), \(D=\sigma^2I\) et lorsque le filtre symétrique est le filtre d’Henderson, on retrouve les filtres asymétriques de Musgrave. Ce filtre suppose, que pour l’estimation en temps réel, les données sont générées par un processus linéaire et que les filtres asymétriques préservent les constantes (\(\sum v_i=\sum \theta_i=1\)). Ces filtres asymétriques dépendent du rapport \(\delta_1/\sigma\), qui est lié à l’I-C ratio \(R=\frac{\bar{I}}{\bar{C}}=\frac{\sum\lvert I_t-I_{t-1}\rvert}{\sum\lvert C_t-C_{t-1}\rvert}\) (\(\delta_1/\sigma=2/(R\sqrt{\pi})\)), qui est notamment utilisé dans X-13ARIMA pour déterminer la longueur du filtre d’Henderson¹⁹.

Lorsque \(U\) correspond aux \(d^*+1\) premières colonnes de \(X\), \(d^*<d\), la contrainte consiste à reproduire des tendances polynomiales de degré \(d^*\). Cela introduit du bais mais réduit la variance. Ainsi, Proietti et Luati (2008) proposent trois classes de filtres asymétriques :

1. Linear-Constant (LC) : \(y_t\) linéaire (\(d=1\)) et \(v\) préserve les constantes (\(d^*=0\)). On obtient le filtre de Musgrave avec le filtre d’Henderson comme filtre symétrique.

2. Quadratic-Linear (QL) : \(y_t\) quadratique (\(d=2\)) et \(v\) préserve les tendances linéaires (\(d^*=1\)).

3. Cubic-Quadratic (CQ) : \(y_t\) cubic (\(d=3\)) et \(v\) préserve les tendances quadratiques (\(d^*=2\)).

Le tableau 5.1 compare les critères de qualité des différentes méthodes en utilisant le filtre d’Henderson et \(h=6\) (filtre symétrique de 13 termes). Pour les filtres en temps réel (\(q=0\)), plus le filtre asymétrique est complexe (en termes de préservation polynomiale), moins la timeliness est élevée et plus la fidelity/smoothness est grande : la réduction du déphasage se fait au détriment d’une augmentation de la variance. Ce résultat varie lorsque \(q\) augmente : pour \(q=2\) le filtre QL a une plus grande timeliness que le filtre LC. Ce résultat étonnant souligne le fait que le déphasage n’est pas contrôlé par l’approche de Proietti et Luati (2008).

En termes de révision, (\(A_w+S_w+T_w+R_w\)), les filtres LC et QL donnent toujours de meilleurs résultats que les filtres CQ et DAF.

Table 5.1 : Critères de qualité des filters asymétriques (\(q=0,1,2\)) calculés par polynômes locaux en utilisant le noyau d’Henderson avec \(h=6\) et \(R=3,5\).

Method	\(b_c\)	\(b_l\)	\(b_q\)	\(F_g\)	\(S_g\)	\(T_g \times 10^{-3}\)	\(A_w\)	\(S_w\)	\(T_w\)	\(R_w\)
\(q=0\)
LC	0	-0,41	-2,16	0,39	1,27	30,34	0,10	0,49	0,41	0,55
QL	0	0,00	-0,47	0,71	5,15	0,05	0,07	1,89	0,00	0,11
CQ	0	0,00	0,00	0,91	11,94	0,01	0,02	2,23	0,00	0,10
DAF	0	0,00	0,00	0,94	14,20	0,00	0,01	2,18	0,00	0,10
\(q=1\)
LC	0	-0,12	-0,52	0,27	0,43	4,80	0,01	0,12	0,06	0,11
QL	0	0,00	-0,06	0,29	0,71	0,69	0,00	0,19	0,01	0,04
CQ	0	0,00	0,00	0,37	0,57	0,16	0,02	0,58	0,00	0,06
DAF	0	0,00	0,00	0,41	0,37	0,06	0,02	0,76	0,00	0,06
\(q=2\)
LC	0	0,00	1,08	0,20	0,08	0,35	0,01	0,01	0,00	0,01
QL	0	0,00	0,03	0,22	0,05	2,08	0,00	0,01	0,02	0,07
CQ	0	0,00	0,00	0,37	0,66	0,13	0,02	0,56	0,00	0,06
DAF	0	0,00	0,00	0,40	0,77	0,02	0,02	0,68	0,00	0,05

Une application en ligne, disponible à l’adresse https://aqlt.shinyapps.io/FiltersProperties/, permet de comparer les coefficients, les fonctions de gain et de déphasage entre les différentes méthodes et les différents noyaux.

Filtres locaux polynomiaux (Proietti et Luati (2008))

Avantages :

• Modèles avec une interprétation simple.

• Le filtre asymétrique est indépendant de la date d’estimation. Toutefois, il dépend indirectement des données si le filtre est calibré sur l’I-C ratio.

Inconvénients :

• La timeliness n’est pas contrôlée (mais peut être introduite dans le programme de minimisation).

5.2 Extension avec le critère de timeliness

Un inconvénient de la méthode précédente est que le déphasage n’est pas contrôlé. Il est en revanche possible de généraliser d’avantage la modélisation en ajoutant le critère de timeliness définit par Grun-Rehomme, Guggemos, et Ladiray (2018) dans l’équation (5.2). C’est ce qui a été proposé par Jean Palate, puis codé en Java²⁰ et intégré dans rjdfilters.

En utilisant les mêmes notations que dans 5.1.4, \(\theta\) le filtre symétrique et \(v\) le filtre asymétrique. Notons également \(\theta=\begin{pmatrix}\theta_p\\\theta_f\end{pmatrix}\) avec \(\theta_p\) de même longueur que \(v\), et \(g=v-\theta_p\). Le critère de timeliness s’écrit : \[ T_g(v)=v'Tv=g'Tg+2\theta_p'Tg+\theta_p'T\theta_p \quad(T\text{ étant symétrique)}. \] De plus, la fonction objectif \(\varphi\) de l’équation (5.2) peut se réécrire : \[\begin{align*} \varphi(v)&=(v-\theta_p)'D_{p}(v-\theta_p)+ \theta_f'D_{f}\theta_f+ [\delta'(Z_{p}'v-Z'\theta)]^{2}+ 2l'(U_{p}'v-U'\theta)\\ &=g'Qg-2Pg+2l'(U_{p}'v-U'\theta)+c\quad\text{avec } \begin{cases} Q=D_p+Z_p\delta\delta'Z'_p \\ P=\theta_fZ_f\delta\delta'Z_p'\\ c\text{ une constante indépendante de }v \end{cases}. \end{align*}\]

En ajoutant le critère de timeliness, on obtient : \[ \widetilde\varphi(v)=g'\widetilde Qg- 2\widetilde Pg+2l'(U_{p}'v-U'\theta)+ \widetilde c\quad\text{avec } \begin{cases} \widetilde Q=D_p+Z_p\delta\delta'Z'_p +\alpha_TT\\ \widetilde P=\theta_fZ_f\delta\delta'Z_p'-\alpha_T\theta_pT\\ \widetilde c\text{ une constante indépendante de }v \end{cases} \] où \(\alpha_T\) est le poids associé au critère de timeliness. Avec \(\alpha_T=0\) on retrouve \(\varphi(v)\). Cette extension permet donc de retrouver tous les filtres symétriques et asymétriques présentés dans la section précédente mais généralise également l’approche de Gray et Thomson (1996) présentée dans la section 5.3. Elle est proche de l’approche FRST.

Cette extension s’inscrit dans le cadre de la théorie générale définie dans la section 3.1. Cela revient en effet à minimiser une somme pondérée de l’erreur quadratique de révision : \[ \E{\left( \sum_{i=-h}^h\theta^s_{i}y_{t+s}-\sum_{i=-h}^qv_iy_{t+s} \right)^2} = I(v,\,0,\,y_t,\,M_{\theta^s} y_t) \] et du critère de timeliness : \[ T_g(\theta) = J(f\colon(\rho,\varphi)\mapsto\rho^2\sin(\varphi)^2,\,\omega_1, \,\omega_2) \] sous une contrainte linéaire.

5.3 Régression polynomiale : Gray et Thomson

5.3.1 Filtres symétriques

L’approche de Gray et Thomson (1996) est proche de celles de Proietti et Luati (2008) et de Grun-Rehomme, Guggemos, et Ladiray (2018). De la même façon que pour les autres méthodes, ils considèrent que la série initiale \(y_t\) peut se décomposer entre une somme entre la tendance-cycle \(g_t\) et d’un bruit blanc \(\varepsilon_t\) de variance \(\sigma^2\) : \[y_t = g_t+\varepsilon_t.\] Toutefois, plutôt que de directement remplacer \(g_t\) par un polynôme local de degré \(d\), ils prennent en compte l’erreur d’approximation de la tendance : \[ g_t=\sum_{j=0}^{d}\beta_{j}t^{j}+\xi_{t}, \] où \(\xi_t\) est un processus stochastique de moyenne nulle, autocorrélé mais non corrélé à \(\varepsilon_t\).

La tendance \(g_t\) est estimée par une moyenne mobile : \[ \hat{g}_{t}=\sum_{s=-r}^{r}\theta_{s}y_{t+s}. \]

Pour le filtre central, les auteurs cherchent à avoir un estimateur \(\hat g_t\) qui soit sans biais (ce qui implique que \(\theta\) conserve les tendances de degré \(d\)) et qui minimise une somme pondérée d’un critère de fidelity et d’un critère de smoothness : \[\begin{equation} Q=\alpha\underbrace{\E{(\hat{g}_{t}-g_{t})^{2}}}_{=F_{GT}}+ (1-\alpha)\underbrace{\E{ (\Delta^{p+1}\hat{g}_{t})^{2}} }_{=S_{GT}} \tag{5.3} \end{equation}\] La solution est un filtre symétrique qui peut s’écrire sous la forme \[ \theta=E_{\alpha}^{-1}X\left[X'E_{\alpha}^{-1}X\right]^{-1}e_{1}\text{ avec }E_{\alpha}=\alpha\left(\sigma^{2}I+\Omega\right)+(1-\alpha)\left(\sigma^{2}B_{p+1}+\Gamma\right) \] où : \[ \begin{cases} \Omega_{jk} & =cov\left(\xi_{t+j}-\xi_{t},\xi_{t+k}-\xi_{t}\right)\\ \Gamma_{jk} & =cov\left(\Delta^{p+1}\xi_{t+j},\Delta^{p+1}\xi_{t+k}\right)\\ \sigma^{2}\left(B_{p+1}\right)_{jk} & =cov\left(\Delta^{p+1}\varepsilon_{t+j},\Delta^{p+1}\varepsilon_{t+k}\right) \end{cases}. \] Les deux critères utilisés dans le programme de minimisation (5.3) sont des cas particuliers du critère \(I\) défini dans l’équation (3.1) : \[\begin{align*} F_{GT}(\theta)&=I(\theta,0,y_t,M_\theta y_t)\\ S_{GT}(\theta)&=I(\theta,d+1,y_t,0). \end{align*}\] La théorie générale définie dans la section 3.1 permet donc de retrouver les filtres de Gray et Thomson (1996).

En ne minimisant que la smoothness et avec \(\xi_t=0\) on retrouve le filtre d’Henderson. En ne minimisant que la fidelity, cette méthode est équivalente à l’estimation de polynômes locaux par moindres carrés généralisés : on retrouve donc les filtres de Proietti et Luati (2008) avec \(\sigma^2=0\) et \(\Omega =K^{-1}\), ainsi que le filtre de Macaulay.

L’avantage de la modélisation de Gray et Thomson est que le paramètre \(\xi_t\) permet une spécification plus précise du modèle en prenant notamment en compte la corrélation entre les observations. Par exemple, McLaren et Steel (2001) ont étudié le lien entre le plan de sondage et l’estimation de la composante tendance-cycle et de la composante saisonnière. Cette modélisation leur permet de prendre en compte, dans l’estimation de la tendance-cycle, la structure de corrélation induite par le plan de sondage de l’enquête emploi mensuelle de l’Australie (groupe de rotations avec une période de recouvrement). Cependant, les auteurs avertissent que dans leur simulations (et dans la modélisation de Gray et Thomson) la structure d’autocorrélation de la variable aléatoire \(\xi_t\) est supposée connue. Ce n’est généralement pas le cas en pratique, où cette structure doit être estimée, ce qui rajoute de l’incertitude dans les estimations.

5.3.2 Filtres asymétriques

L’approche retenue par Gray et Thomson (1996) est une approche de minimisation des révisions sous contraintes. Étant donné un filtre symétrique \(\theta^s\) utilisé pour estimer la tendance au centre de la série, l’objectif est de chercher un filtre asymétrique \(v=(v_{-h},\dots,v_q)\) de sorte à minimiser l’erreur quadratique moyenne de révision : \[ \E{\left(Y-\hat Y\right)^2} = \E{\left( \sum_{i=-h}^h\theta^s_iy_{t+s}-\sum_{i=-h}^qv_iy_{t+s} \right)^2}. \] Les auteurs étudient deux cas :

1. Dans le premier cas, ils cherchent un estimateur sans biais : cela implique que \(v\) conserve les mêmes tendances polynomiales que \(\theta^s\). \(\hat Y\) est alors le meilleur prédicteur linéaire sans biais — best linear unbiased predictor (BLUP) — de \(Y\).

2. Dans le second cas, ils autorisent l’estimateur à être biaisé mais imposent que ce biais soit constant dans le temps : si l’on modélise localement la tendance par un polynôme de degré \(d\), cela implique que \(v\) conserve les tendances polynomiales de degré \(d-1\). \(\hat Y\) est alors le meilleur prédicteur linéaire à biais constant — best linear time invariant predictor (BLIP) — de \(Y\). Cela permet notamment de reproduire les filtres asymétriques de Musgrave.

La méthode utilisée est donc très proche de celle de Proietti et Luati (2008) : on retrouve d’ailleurs le filtre DAF avec \(\sigma^2=0\) et \(\Omega =K^{-1}\) et en utilisant la première méthode (estimation du BLUP) et les méthodes LC (filtre de Musgrave), QL et CQ avec la seconde méthode en utilisant respectivement \(d=1\), \(d=2\) et \(d=3\).

La théorie générale définie dans la section 3.1 permet également de retrouver les filtres asymétriques puisqu’ils sont construits en minimisant l’erreur quadratique moyenne des révisions sous contraintes linéaires (préservation d’un polynôme de degré \(p\)).

Pour la construction des filtres asymétriques, une approche alternative pourrait être d’utiliser la même méthode que celle utilisée pour construire les filtres symétriques. C’est-à-dire minimiser \(Q\) (équation (5.3)) sous contrainte que le filtre asymétrique fournisse un estimateur dans biais de la tendance. Comme discuté dans Gray et Thomson (1996), les auteurs ne retiennent pas cette méthode pour deux raisons :

• Il n’est pas évident qu’il faudrait chercher à maintenir le même équilibre entre smoothness et fidelity en fin de série et au centre de la série. Le problème rencontré en fin de série est transitoire et disparaît au fur et à mesure que l’on a de nouvelles observations. Minimiser des critères de révision serait donc préférable puisque cela reviendrait à minimiser le coût de la transition (mais dans le cas où l’on ne minimise que la fidelity les deux méthodes sont équivalentes).

• Les valeurs de la fidelity et de la smoothness ne dépendent pas du temps au centre de la série mais en dépendent en fin de série. Ainsi, même si au centre de la série le choix des poids entre les deux critères contrôle indirectement le niveau des indicateurs, ce n’est plus le cas en fin de série. De plus, en fin de série, cela pourrait introduire des déphasages plus importants car \(S_{GT}\) dépend du temps et des valeurs passées (du fait du l’utilisation de l’opérateur différence).

Inversement, Grun-Rehomme, Guggemos, et Ladiray (2018) justifie de ne pas intégrer le critère de révision dans leur problème car ce critère est fortement corrélé à une combinaison fixée, donc non ajustable par l’utilisateur, des critères fidelity et timeliness.

Filtres locaux polynomiaux (Gray et Thomson (1996))

Avantages :

• Modèles généraux qui permettent de prendre en compte l’autocorrélation entre les observations.

• Interprétation statistique des différentes méthodes.

• Le filtre asymétrique est indépendant de la date d’estimation. Toutefois, il dépend indirectement des données si le filtre est calibré sur l’I-C ratio.

Inconvénients :

• La timeliness n’est pas contrôlée.

• La spécification du modèle (i.e., du paramètre \(\xi_t\)) peut être compliquée : si la structure d’autocorrélation est estimée à partir des données, cela rajoute de l’incertitude dans les estimations, ce qui peut avoir des effets indésirables.

5.4 Reproducing Kernel Hilbert Space (RKHS) : approche de Dagum et Bianconcini

La théorie des Reproducing Kernel Hilbert Space (RKHS) — espaces de Hilbert à noyau reproduisant — est une théorie générale dans l’apprentissage statistique non-paramétrique qui permet d’englober un grand nombre de méthodes. C’est par exemple le cas des méthodes de régression par moindres carrés pénalisés, des Support Vector Machine (SVM), du filtre d’Hodrick-Prescott (utilisé pour décomposer tendance et cycle) ou encore des moyennes mobiles telles que celle d’Henderson. Ainsi, Dagum et Bianconcini (2008) utilise la théorie des RKHS pour approcher le filtre d’Henderson et en dériver des filtres asymétriques associés.

Un RKHS \(\mathbb{L}^{2}(f_{0})\) est un espace de Hilbert caractérisé par un noyau qui permet de reproduire toutes les fonctions de cet espace. Il est caractérisé par une fonction de densité \(f_0\) et un produit scalaire \(\ps{\cdot}{\cdot}\) définit par : \[ \left\langle U(t),V(t)\right\rangle =\E{U(t)V(t)}=\int_{\R}U(t)V(t)f_{0}(t)\ud t\quad \forall U,V\in\mathbb{L}^{2}(f_{0}). \] La fonction \(f_0\) pondère donc chaque valeur en fonction de sa position temporelle : il s’agit de la version continue des noyaux définis dans la partie 5.1.2.

Dans notre cas, on suppose que notre série initiale \(y_t\) est désaisonnalisée et peut s’écrire comme la somme d’une tendance-cycle, \(TC_t\), et d’une composante irrégulière, \(I_t\) (qui peut être un bruit blanc ou suivre un modèle ARIMA) : \(y_t=TC_t+I_t.\) La tendance-cycle peut être déterministe ou stochastique. On suppose que c’est une fonction régulière du temps, elle peut être localement approchée par un polynôme de degré \(d\) : \[ TC_{t+j}=TC_t(j)=a_0+a_1j+\dots+a_dj^d+\varepsilon_{t+j},\quad j\in\llbracket-h,h\rrbracket, \] où \(\varepsilon_t\) est un bruit blanc non corrélé à \(I_t\).

Les coefficients \(a_0,\dots,a_d\) peuvent être estimés par projection des observations au voisinage de \(y_t\) sur le sous-espace \(\mathbb P_d\) des polynômes de degré \(d\), ou, de manière équivalente, par minimisation de la distance entre \(y_t\) et \(TC_t(j)\) : \[\begin{equation} \underset{TC\in\mathbb P_d}{\min}\lVert y -TC \rVert^2 = \underset{TC\in\mathbb P_d}{\min}\int_\R (y(t+s)-TC_t(s))^2f_0(s)\ud s. \tag{5.4} \end{equation}\] L’espace \(\mathbb P_d\) étant un espace de Hilbert à dimension finie, il admet un noyau reproduisant (voir, par exemple, Berlinet et Thomas-Agnan (2004)). Il existe ainsi une fonction \(R_d(\cdot,\cdot)\) telle que : \[ \forall P\in \mathbb P_d: \forall t: R_d(t,\cdot)\in\mathbb P_d\quad\text{et}\quad P(t)=\ps{R_d(t,\cdot)}{P(\cdot)}. \]

Le problème (5.4) admet une solution unique qui dépend d’une fonction \(K_{d+1}\), appelée fonction de noyau (kernel function). Cette fonction est dite d’ordre \(d+1\) car elle conserve les polynômes de degré \(d\)²¹. Cette solution s’écrit : \[\begin{equation} \widehat{TC}(t)=\int_\R y(t-s)K_{d+1}(s) \ud s. \tag{5.5} \end{equation}\] Généralement \(f_0(t) = 0\) pour \(\lvert t \rvert>1\). Cette solution s’écrit alors : \[\begin{equation} \widehat{TC}(t)=\int_{[-1,1]} y(t-s)K_{d+1}(s) \ud s. \tag{5.6} \end{equation}\] On peut, par ailleurs, montrer que \(K_{d+1}\) s’écrit en fonction \(f_0\) et du noyau reproduisant \(K(\cdot,\cdot)\) et que ce dernier peut s’écrire en fonction de polynômes \((P_i)_{i\in \llbracket 0, d-1 \rrbracket}\) qui forme une base orthonormée de \(\mathbb L^2(f_0)\) (voir par exemple Berlinet (1993)) : \[ K_{d+1}(t) = R_d(t,0)f_0(t) = \sum_{i=0}^dP_i(t)P_i(0)f_0(t). \]

De plus, dans le cas discret, la solution (5.6) s’écrit comme une somme pondérée au voisinage de \(y_t\) : \[\begin{equation} \widehat{TC}_t=\sum_{j=-h}^h w_j y_{t+j}\quad \text{où} \quad w_j=\frac{K_{d+1}(j/b)}{\sum_{i=-h}^{^h}K_{d+1}(i/b)}. \tag{5.7} \end{equation}\] Le paramètre \(b\) est choisi de sorte que les \(2h+1\) points autour de \(y_t\) soient utilisés avec un poids non nul.

Pour les filtres asymétriques, la formule (5.7) est simplement adaptée au nombre d’observations connues : \[\begin{equation} \forall j\in\left\llbracket -h,q\right\rrbracket\::\: w_{a,j}=\frac{K_{d+1}(j/b)}{\sum_{i=-h}^{^q}K_{d+1}(i/b)}. \tag{5.8} \end{equation}\] En utilisant \(b=h+1\) on retrouve les filtres symétriques obtenues par polynômes locaux.

Comme notamment montré par Dagum et Bianconcini (2016), \(K_{d+1}\) peut s’exprimer simplement à partir des moments de \(f_0\)²². Ainsi, notons \(H_{d+1}\) la matrice de Hankel associée aux moments de \(f_0\) : \[ \forall i,j\in \llbracket 0, d\rrbracket: \left(H_{d+1}\right)_{i,j}=\ps{X^i}{X^j}=\int s^{i+j}f_0(s)\ud s. \] Notons également \(H_{d+1}[1,t]\) la matrice obtenue en remplaçant la première ligne de \(H_{d+1}\) par \(\begin{pmatrix} 1 & t & t^2 & \dots & t^d\end{pmatrix}\). On a : \[\begin{equation} K_{d+1}(t)=\frac{\det{H_{d+1}[1,t]}}{\det{H_{d+1}}}f_0(t). \tag{5.9} \end{equation}\] C’est cette formule qui est utilisée dans le package rjdfilters pour calculer les différentes moyennes mobiles.

Comme discuté dans la partie 5.1, le noyau d’Henderson dépend de la fenêtre utilisée. Ainsi, tous les moments de l’équation (5.9) doivent être recalculés pour chaque valeur de \(h\). Pour éviter cela, Dagum et Bianconcini (2008) suggèrent d’utiliser le noyau quadratique (biweight) pour approcher le noyau d’Henderson lorsque \(h\) est petit (\(h< 24\)) et le noyau cubique (triweight) lorsque \(h\) est grand \(h\geq 24\).

Dans Dagum et Bianconcini (2015), les auteures suggèrent de faire une sélection optimale du paramètre \(b\), par exemple en minimisant l’erreur quadratique moyenne (option "frequencyresponse" dans rjdfilters::rkhs_filter) : \[ b_{q,\Gamma}=\underset{b_q\in[h; 3h]}{\min} 2\int_{0}^{\pi} \lvert \Gamma_s(\omega)-\Gamma_\theta(\omega)\rvert^2\ud \omega. \] Cela suppose en fait que la série entrée \(y_t\) est un bruit blanc. En supposant \(y_t\) stationnaire, les critères définis dans l’article originel peuvent donc être étendus en multipliant les quantités sous les intégrales par la densité spectrale de \(y_t\) notée \(h\) : \[ b_{q,\Gamma}=\underset{b_q\in[h; 3h]}{\min} 2\int_{0}^{\pi} \lvert \Gamma_s(\omega)-\Gamma_\theta(\omega)\rvert^2h(\omega)\ud \omega. \] Cette erreur quadratique moyenne peut également se décomposer en plusieurs termes (voir équation (4.1) de la section 4) qui peuvent également être minimisés :

• l’accuracy qui correspond à la part de la révision liée aux différences de fonction de gain dans les fréquences liées à la tendance-cycle \[ b_{q,G}=\underset{b_q\in[h; 3h]}{\min} 2\int_{0}^{\omega_1} \left(\rho_s(\omega)-\rho_\theta(\omega)\right)^{2} h(\omega)\ud \omega \]

• la smoothness qui correspond à la part de la révision liée aux différences de fonction de gain dans les fréquences liées aux résidus \[ b_{q,s}=\underset{b_q\in[h; 3h]}{\min} 2\int_{\omega_1}^{\pi} \left(\rho_s(\omega)-\rho_\theta(\omega)\right)^{2} h(\omega)\ud \omega \]

• la timeliness qui correspond à la part de la révision liée au déphasage \[ b_{q,\varphi}=\underset{b_q\in[h; 3h]}{\min} 8\int_{0}^{\omega_1} \rho_s(\lambda)\rho_\theta(\lambda)\sin^{2}\left(\frac{\varphi_\theta(\omega)}{2}\right)h(\omega)\ud \omega \] Dans rjdfilters \(h\) peut être fixée à la densité spectrale d’un bruit blanc (\(h_{WN}(x)=1\), comme c’est le cas dans Dagum et Bianconcini (2015)) ou d’une marche aléatoire (\(h_{RW}(x)=\frac{1}{2(1-\cos(x))}\)).

Pour assurer une cohérence dans les définitions entre les différentes sections, les définitions de \(b_{q,G}\), \(b_{q,\Gamma}\) et \(b_{q,\varphi}\) ont été légèrement modifiées par rapport à celles définies dans Dagum et Bianconcini (2015) où :

• dans \(b_{q,G}\) le terme à minimiser est sous une racine carrée (sans impact sur le minimum) : \[ b_{q,G}=\underset{b_q}{\min}\sqrt{ 2\int_{0}^{\pi} \left(\rho_s(\omega)-\rho_\theta(\omega)\right)^{2}\ud \omega} \]

• \(b_{q,s}\) n’est pas considéré, \(b_{q,\Gamma}\) et \(b_{q,G}\) sont définis avec \(\omega_1=\pi\)

• \(b_{q,\varphi}\) est défini par : \[ b_{q,\varphi}=\underset{b_q}{\min} \sqrt{2\int_{\Omega_S} \rho_s(\lambda)\rho_\theta(\lambda)\sin^{2}\left(\frac{\varphi_\theta(\omega)}{2}\right)\ud \omega} \] où \(\Omega_S=[0,2\pi/32]\) est l’intervalle de fréquences associées aux cycles d’au moins 16 mois.

• une formule différente est utilisée pour la fonction de réponse (\(\Gamma_\theta(\omega)=\sum_{k=-p}^{+f} \theta_k e^{2\pi i \omega k}\)), ce qui conduit à des bornes d’intégrales légèrement différentes, sans effet sur le résultat.

Un des inconvénients de cette méthode est qu’il n’y a pas unicité de la solution et donc qu’il y a parfois plusieurs extremum (uniquement pour le calcul de \(b_{q,\varphi}\). Ainsi, la valeur optimal retenue par défaut par rjdfilters produit des discontinuités dans l’estimation de la tendance-cycle.
Par ailleurs, les valeurs de \(b_{q,G}\) varient fortement en fonction de si l’on retient \(\omega_1=2\pi/12\) ou \(\omega_1=\pi\) (tableau 5.2). Par cohérence, simplicité, nous utiliserons dans cet article les valeurs optimales présentées dans Dagum et Bianconcini (2015).

Table 5.2 : Fenêtres optimales pour les filtres asymétriques associés à un filtre symétrique de 13 termes (\(h=6\)) avec le noyau biweight.

	\(q=0\)	\(q=1\)	\(q=2\)	\(q=3\)	\(q=4\)	\(q=5\)
\(b_{q,\Gamma}\)
\(\omega_1 = \pi\)	9,54	7,88	7,07	6,88	6,87	6,94
\(\omega_1 = 2\pi/12\)	9,54	7,88	7,07	6,88	6,87	6,94
\(b_{q,G}\)
\(\omega_1 = \pi\)	11,78	9,24	7,34	6,85	6,84	6,95
\(\omega_1 = 2\pi/12\)	8,61	7,64	6,01	6,01	6,01	6,59
\(b_{q,\varphi}\)
Valeurs de l’article (avec \(\omega_1 = 2\pi/36\))	6,01	6,01	7,12	8,44	9,46	10,39
\(\omega_1 = 2\pi/36\)	6,01	6,01	7,21	8,47	9,46	6,01
\(\omega_1 = 2\pi/12\)	6,01	6,01	6,38	8,15	9,35	6,01

RKHS filters - Dagum et Bianconcini (2008)

Avantages :

• Le filtre asymétrique est indépendant des données et de la date d’estimation.

• La méthode est généralisable à des séries avec des fréquences irrégulières (par exemple avec beaucoup de valeurs manquantes).

Inconvénient :

• Il peut y avoir des problèmes de minimisation (notamment en minimisant la timeliness).

Tous les critères utilisés de la sélection optimale du paramètre \(b\) pouvant s’écrire comme cas particuliers du critère \(J\) défini dans l’équation (3.1) (voir section 4), cette méthode s’inscrit dans le cadre de la théorie générale définie dans la section 3.1 en imposant comme contrainte linéaire que les coefficients soient sous la forme \(w_j=\frac{K_{d+1}(j/b)}{\sum_{i=-h}^{^p}K_{d+1}(i/b)}\).

5.5 Liens entre les différentes méthodes

5.5.1 Critères de Gray et Thomson et ceux de Grun-Rehomme et alii

Les critères \(F_g\) et \(S_g\) peuvent se déduire de \(F_{GT}\) et \(S_{GT}\). Les approches de Gray et Thomson (1996) et Grun-Rehomme, Guggemos, et Ladiray (2018) sont donc équivalentes pour la construction de filtres symétriques.

Notons \(x_{t}=\begin{pmatrix}1 & t & t^{2} & \cdots & t^{d}\end{pmatrix}\), \(\beta_{t}=\begin{pmatrix}\beta_{0} & \cdots & \beta^{d}\end{pmatrix}'\).

Pour le critère de fidelity : \[ \hat{g}_{t}-g_{t}=\left(\sum_{j=-h}^{+h}\theta_{j}x_{t+j}-x_{t}\right)\beta+\sum_{j=-h}^{+h}\theta_{j}\varepsilon_{t+j}+\sum_{j=-h}^{+h}\theta_{j}(\xi_{t+j}-\xi_{t}), \] Si \(\theta\) préserve les polynômes de degré \(d\) alors \(\sum_{j=-h}^{+h}\theta_{j}x_{t+j}=x_{t}\). où \(x_{t}=\begin{pmatrix}1 & t & t^{2} & \cdots & t^{d}\end{pmatrix}\). Puis, comme \(\xi_{t}\) et \(\varepsilon_{t}\) sont de moyenne nulle et sont non corrélés : \[ F_{GT}(\theta)=\E{(\hat{g}_{t}-g_{t})^{2}}=\theta^{'}\left(\sigma^{2}I+\Omega\right)\theta. \] Si \(\xi_t=0\) alors \(\Omega=0\) et \(F_{GT}(\theta)=F_g(\theta)\).

Pour la smoothness on a : \[ \nabla^{q}\hat{g}_{t}=\sum_{j=h}^{h}\theta_{j}\underbrace{\nabla^{q}\left(\left(x_{j}-x_{0}\right)\beta\right)}_{=0\text{ si }q\geq d+1}+\sum_{j=h}^{h}\theta_{j}\nabla^{q}\varepsilon_{t+j}+\sum_{j=h}^{h}\theta_{j}\nabla^{q}\xi_{t+j}. \] D’où pour \(q=d+1\) : \[ S_{GT}(\theta)=\E{(\nabla^{q}\hat{g}_{t})^{2}}=\theta^{'}\left(\sigma^{2}B_{q}+\Gamma_{q}\right)\theta. \] On peut par ailleurs montrer que pour toute série temporelle \(X_t\), \[ \nabla^{q}(M_{\theta}X_{t})=\left(-1\right)^{q}\sum_{k\in\Z}\left(\nabla^{q}\theta_{k}\right)X_{t+k-q} \] avec \(\theta_k=0\) pour \(|k|\geq h+1\). Avec \(\xi_t=0\) on trouve donc que \(S_{GT}(\theta)=\sigma^2S_g(\theta)\).

5.5.2 Équivalence avec les moindres carrés pondérés

Du fait de la forme des filtres obtenus par la méthode de Grun-Rehomme, Guggemos, et Ladiray (2018), lorsque les contraintes imposées sont la préservation des tendances de degré \(d\), celle-ci est équivalente à une estimation locale d’une tendance polynomiale de degré \(d\) par moindres carrés généralisés. En effet, dans ce cas, la solution est \(\hat \theta = \Sigma^{-1}X_p'\left(X_p\Sigma^{-1}X_p'\right)^{-1}e_1\) avec \(\Sigma=\alpha F+\beta S+ \gamma T\), et c’est l’estimation de la constante obtenue par moindres carrés généralisés lorsque la variance des résidus est \(\Sigma\). L’équivalence entre les deux méthodes peut donc se voir comme un cas particulier de l’équivalence entre les moindres carrés pondérés et les moindres carrés généralisés. C’est par exemple le cas des filtres symétriques d’Henderson qui peuvent s’obtenir par les deux méthodes.

Dans ce sens, Henderson (1916) a montré que les poids \(w=(w_{-p},\dots w_{f})\) associés à une moyenne mobile issue de la régression polynomiale locale par moindres carrés pondérés pouvaient s’écrire sous la forme : \[ w_i = \kappa_i P\left(\frac{i}{p+f+1}\right)\text{ où }P\text{ est un polynôme de degré }d. \] Il a également montré l’inverse : toute moyenne mobile \(\theta=(\theta_{-p},\dots, \theta_{f})\) qui préserve les tendances de degré \(d\) et dont le diagramme des coefficients change au plus \(d\) fois de signes peut être obtenue par une régression polynomiale locale de degré \(p\) estimée par moindres carrés pondérés. Pour cela il suffit de trouver un polynôme \(P\left(\frac{X}{p+f+1}\right)\) de degré inférieur ou égal à \(d\) et dont les changements de signes coïncident avec les changements de signes de \(\theta\). Le noyau associé est alors \(\kappa_i=\frac{ \theta_i}{P\left(\frac{i}{p+f+1}\right)}\). C’est le cas de tous les filtres symétriques issues de l’approche FST et de la majorité des filtres asymétriques. L’annexe C présente les quelques poids pour lesquels il n’y a pas équivalence.

Plus récemment, Luati et Proietti (2011) se sont intéressés aux cas d’équivalences entre les moindres carrés pondérés et les moindres carrés généralisés pour déterminer des noyaux optimaux (au sens de Gauss-Markov). Ils montrent que le noyau d’Epanechnikov est le noyau optimal associé à la régression polynomiale locale où le résidu, \(\varepsilon_t\), est un processus moyenne mobile (MA) non inversible d’ordre 1 (i.e., \(\varepsilon_t=(1-B)\xi_t\), avec \(\xi_t\) un bruit blanc). Dans ce cas, la matrice \(\Sigma\) de variance-covariance correspond à la matrice obtenue par le critère de smoothness avec le paramètre \(q=2\) (\(\sum_{j}(\nabla^{2}\theta_{j})^{2} = \theta'\Sigma\theta\)) : il y a donc équivalence avec l’approche FST. De même, le noyau d’Henderson est le noyau optimal associé à la régression polynomiale locale où le résidu est un processus moyenne mobile (MA) non inversible d’ordre 2 (i.e., \(\varepsilon_t=(1-B)^2\xi_t\), avec \(\xi_t\) un bruit blanc).

5.5.3 RKHS et polynômes locaux

Comme montré dans la section précédente, la théorie des espaces de Hilbert à noyau reproduisant permet de reproduire les filtres symétriques par approximation polynomiale locale. Comme le montrent Luati et Proietti (2011), cette théorie permet donc également de reproduire les filtres directs asymétriques (DAF), qui sont équivalents à l’approximation polynomiale locale mais en utilisant une fenêtre d’estimation asymétrique. Cependant, ils ne peuvent pas être obtenus par la formalisation de Dagum et Bianconcini (2008) mais par une discrétisation différente de la formule (5.9) : \[ K_{d+1}(t)=\frac{\det{H_{d+1}[1,t]}}{\det{H_{d+1}}}f_0(t). \] Dans le cas discret, \(f_0(t)\) est remplacé par \(\kappa_j\) et en remplaçant les moments théoriques par les moments empiriques \(H_{d+1}\) devient \(X'_pK_pX_p\) et les coefficients du filtre asymétrique sont obtenus en utilisant la formule : \[ w_{a,j}=\frac{\det{X'_pK_pX_p[1,j]} }{ \det{X'_pK_pX_p} }\kappa_j. \] En effet, la règle de Cramer permet de trouver une solution explicite à l’équation des moindres carrés \((X'_pK_pX_p)\hat \beta=X'_pK_py_p\) où \(\hat \beta_0=\hat m_t\) : \[ \hat \beta_0 = \frac{\det{X'_pK_pX_p[1,b]}}{\det{X'_pK_pX_p}}f_0(t) \quad\text{où}\quad b=X'_pK_py_p. \] Comme \(b=\sum_{j=-h}^qx_j\kappa_jy_{t+j}\) il vient : \[ \det{X'_pK_pX_p[1,b]} = \sum_{j=-h}^q\det{X'_pK_pX_p[1,x_j]}\kappa_jy_{t+j}. \] Et enfin : \[ \hat \beta_0 = \hat m_t= \sum_{j=-h}^q\frac{\det{X'_pK_pX_p[1,j]} }{ \det{X'_pK_pX_p} }\kappa_j y_{t+j}. \]

Références

Berlinet, Alain. 1993. « Hierarchies of higher order kernels ». Probability Theory and Related Fields 94: 489‑504.

Berlinet, Alain, et Christine Thomas-Agnan. 2004. Reproducing Kernel Hilbert Spaces in Probability and Statistics. Springer.

Cleveland, Robert B., William S. Cleveland, Jean E. McRae, et Irma Terpenning. 1990. « STL: A Seasonal-Trend Decomposition Procedure Based on Loess (with Discussion) ». Journal of Official Statistics 6: 3‑73.

Cleveland, William S., et Clive Loader. 1996. « Smoothing by local regression: Principles and methods ». In Statistical theory and computational aspects of smoothing, 10‑49. Springer.

Dagum, Estela Bee, et Silvia Bianconcini. 2008. « The Henderson Smoother in Reproducing Kernel Hilbert Space ». Journal of Business & Economic Statistics 26: 536‑45. https://ideas.repec.org/a/bes/jnlbes/v26y2008p536-545.html.

———. 2015. « A new set of asymmetric filters for tracking the short-term trend in real-time ». The Annals of Applied Statistics 9 (3): 1433‑58. https://doi.org/"10.1214/15-AOAS856".

———. 2016. Seasonal adjustment methods and real time trend-cycle estimation. Springer.

Gray, Alistair, et Peter Thomson. 1996. « Design of Moving-Average Trend Filters using Fidelity and Smoothness Criteria ». In Athens Conference on Applied Probability and Time Series Analysis, édité par P. M. Robinson et Murray Rosenblatt, 205‑19. New York, NY: Springer New York.

Grun-Rehomme, Michel, Fabien Guggemos, et Dominique Ladiray. 2018. « Asymmetric Moving Averages Minimizing Phase Shift ». Handbook on Seasonal Adjustment. ec.europa.eu/eurostat/web/products-manuals-and-guidelines/-/KS-GQ-18-001.

Henderson, Robert. 1916. « Note on graduation by adjusted average ». Transactions of the actuarial society of America 17: 43‑48.

Loader, Clive. 1999. Local regression and likelihood. New York: Springer-Verlag.

Luati, Alessandra, et Tommaso Proietti. 2011. « On the equivalence of the weighted least squares and the generalised least squares estimators, with applications to kernel smoothing ». Annals of the Institute of Statistical Mathematics 63 (4): 851‑71. https://doi.org/10.1007/s10463-009-0267-8.

Macaulay, Frederick R et al. 1931. « The smoothing of time series ». NBER Books.

McLaren, Craig H, et David G Steel. 2001. « Rotation patterns and trend estimation for repeated surveys using rotation group estimates ». Statistica Neerlandica 55 (2): 221‑38. https://documents.uow.edu.au/~craigmc/sn_2001.pdf.

Musgrave, John C. 1964. « A set of end weights to end all end weights ». US Census Bureau [custodian]. https://www.census.gov/ts/papers/Musgrave1964a.pdf.

Proietti, Tommaso, et Alessandra Luati. 2008. « Real time estimation in local polynomial regression, with application to trend-cycle analysis ». Ann. Appl. Stat. 2 (4): 1523‑53. https://doi.org/10.1214/08-AOAS195.

---

15. La série est donc désaisonnalisée.↩︎

16. \(\theta\) est symétrique du fait de la symétrie des noyaux \(\kappa_j\).↩︎

17. Voir par exemple W. S. Cleveland et Loader (1996) ou Loader (1999). Les seules contraintes souhaitées sur le noyau est qu’il accorde un poids plus important à l’estimation centrale (\(\kappa_0\)) et qu’il décroit vers 0 lorsque l’on s’éloigne de l’estimation centrale. Le noyau uniforme est donc à éviter.↩︎

18. Dans rjdfilters \(\sigma^2\) est fixé arbitrairement à \(\sigma^2=0,25\).↩︎

19. Dans la majorité des cas un filtre de 13 termes est utilisé. Si le ratio est grand alors un filtre de 23 termes est utilisé (pour supprimer d’avantage de bruit) et si le ratio est petit un filtre de 9 termes est utilisé.↩︎

20. https://github.com/palatej/jdemetra-core.↩︎

21. C’est-à-dire \(\int_\R K_{d+1}(s)\ud s = 1\) et \(\int_\R K_{d+1}(s) s^i\ud s = 1\) pour \(i\in \llbracket 1, d\rrbracket\).↩︎

22. Cela vient en fait du procédé d’orthonomalisation de Gram-Schmidt.↩︎