% stephane.adjemian@ens.fr
% florian.pelgrin@unil.ch
\documentclass[10pt,subeqn]{article}

 \usepackage[latin1]{inputenc}
 \usepackage[T1]{fontenc}
 \usepackage{palatino}
 \usepackage{natbib}
 \bibliographystyle{mystyle}

%\usepackage{float}
%\usepackage{graphicx}
%\usepackage{dcolumn}
\usepackage{hyperref}
\usepackage[francais]{babel}
\usepackage{setspace}
\usepackage[active]{srcltx}
\usepackage{amssymb}
\usepackage[centertags]{amsmath}
\usepackage{amstext}
\usepackage{amsbsy}
\usepackage{amsopn}
\usepackage{amscd}
\usepackage{amsxtra}
\usepackage{amsthm}
%\usepackage{ae,aeguill}
%\usepackage{endnotes}

\newtheorem{theorem}{Theorem}%[section]
\newtheorem{algorithm}[theorem]{Algorithme}
\newtheorem{axiom}[theorem]{Axiome}
\newtheorem{case}[theorem]{Cas}
\newtheorem{conclusion}[theorem]{Conclusion}
\newtheorem{condition}[theorem]{Condition}
\newtheorem{conjecture}[theorem]{Conjecture}
\newtheorem{corollary}[theorem]{Corollaire}
\newtheorem{criterion}[theorem]{Crit\`{e}re}
\newtheorem{definition}[theorem]{D\'{e}finition}
\newtheorem{example}[theorem]{Exemple}
\newtheorem{lemma}[theorem]{Lemme}
\newtheorem{notation}[theorem]{Notation}
\newtheorem{problem}[theorem]{Probl\`{e}me}
\newtheorem{proposition}[theorem]{Proposition}
\newtheorem{remark}[theorem]{Remarque}
\newtheorem{summary}[theorem]{R\'{e}sum\'{e}}

%\newcolumntype{d}{D{,}{,}{-1}}
\newcommand{\norm}[1]{\left\Vert#1\right\Vert}
\newcommand{\abs}[1]{\left\vert#1\right\vert}
\newcommand{\set}[1]{\left\{#1\right\}}
\newcommand{\croc}[1]{\left[#1\right]}
\newcommand{\parent}[1]{\left(#1\right)}
\newcommand{\Real}{\mathbb R}
\newcommand{\sample}{\mathcal Y_T^*}
\newcommand{\samplet}[1]{\mathcal Y_{#1}^*}
\newcommand{\eps}{\varepsilon}
\newcommand{\trace}{\mathrm{tr}}
\newcommand{\To}{\longrightarrow}
\newcommand{\BX}{\mathbf{B}(X)}
\newcommand{\normal}[2]{ \mathcal{N}\left(#1,#2\right) }
\newcommand{\iid}[2]{ \text{iid}\left(#1,#2\right) }
\newcommand{\boxref}[1]{[\ref{#1}]}
\newcommand{\equaref}[1]{(\ref{#1})}
\newcommand{\fullref}[1]{[\ref{#1}, en page \pageref{#1}]}
\newcommand{\VEC}{\mathrm{vec}}

\setlength{\hoffset}{-18pt}
  \setlength{\oddsidemargin}{10pt}   % Marge gauche sur pages impaires
  \setlength{\evensidemargin}{9pt}   % Marge gauche sur pages paires
  \setlength{\marginparwidth}{54pt}   % Largeur de note dans la marge
  \setlength{\textwidth}{481pt}   % Largeur de la zone de texte (17cm)
  \setlength{\voffset}{-18pt}   % Bon pour DOS
  \setlength{\marginparsep}{7pt}   % Séparation de la marge
  \setlength{\topmargin}{0pt}   % Pas de marge en haut
  \setlength{\headheight}{7pt}%13pt}   % Haut de page
  \setlength{\headsep}{10pt}   % Entre le haut de page et le texte
  \setlength{\footskip}{27pt}   % Bas de page + séparation
  \setlength{\textheight}{22cm}%708pt}   % Hauteur de la zone de texte (25cm)

%% À décocher pour éco et prev
%\renewcommand{\footnote}{\endnote}
%\renewcommand{\enotesize}{\normalsize}

\begin{document}



\title{Un regard Bayésien sur les Modèles Dynamiques de la Macroéconomie}
\author{\texttt{Stéphane Adjemian\footnote{stephane.adjemian@ens.fr}}\\\texttt{Université du Maine, \textsc{gains} et \textsc{cepremap}}  \and
\texttt{Florian Pelgrin\footnote{florian.pelgrin@unil.ch}}\\ \texttt{Université de Lausanne - \textsc{hec}, \textsc{iems} et \textsc{cirano}}}

\date{\today}

\maketitle

\begin{abstract}
L'objet de cette contribution est de présenter l'approche bayésienne
des modèles dynamiques les plus considérés en macroéconomie~: les
modèles DSGE (\textit{Dynamic Stochastic General Equilibrium}) et
les modèles VAR.
\medskip

\noindent \textbf{Classification JEL :} C3, C5, E3

\noindent \textbf{Mots-Clefs :} \'Économétrie Bayésienne, Modèles
VAR, Modèles DSGE.

\end{abstract}



\bigskip
\bigskip


\begin{spacing}{1.3}

\section{Introduction}\label{intro}

\par{Ces dernières années, l'analyse des fluctuations économiques s'est
développée autour des Modèles d'\'{E}quilibre Général Intertemporels
Stochastiques (DSGE). Pour autant, jusqu'à très récemment, l'engouement 
pour l'approche DSGE comme outil d'analyse de la politique économique est demeuré relativement
faible, et l'approche des modèles Vectoriels Autorégressif (VAR) a
été (est) souvent privilégiée. Plusieurs raisons expliquent cette
préférence. D'une part, la modélisation VAR de la dynamique des
variables macro-économiques impose un nombre très restreint de
contraintes et offre une qualité d'ajustement aux données (et des
prévisions) relativement bonne. Au contraire, en augmentant le
nombre de contraintes sur les données, encourant ainsi le risque
d'une mauvaise spécification, les modèles DSGE de la première
génération (les modèles de la théorie des cycles réels) se sont
traduits par des performances d'ajustement et de prévisions très
pauvres. D'autre part, l'émergence d'une approche plus structurelle
des modèles VAR (par rapport à l'approche a-théorique, \cite{Sims80}
)--- autorisant des procédures d'identification des chocs à partir
de restrictions contemporaines, de court terme (\cite{Sims86};
\cite{Bernanke86}) ou de long terme (\cite{BlanchardQuah89})---ont
conduit à exiger que tout modèle théorique puisse reproduire les
fonctions de réponse des variables macro-économiques à des chocs
structurels identifiés dans les modèles VAR
(\cite{RotembergWoodford97}, \cite{ChristianoEichenbaumEvans03}).
Finalement, l'absence d'un traitement économétrique convaincant n'a
fait que renforcer la recommandation de
\cite{KydlandPrescott82}---l'étalonnage est préférable.}\newline

\par{Cependant, on a constaté un regain d'intérêt des modèles DSGE et cela essentiellement pour deux
raisons~: (\textit{i}) les avancées théoriques et notamment la prise
en compte de fondements micro-économiques des rigidités nominales
et/ou réelles (\textit{ii}) les progrès dans l'estimation et
l'évaluation des modèles sur la base de méthodes statistiques
formelles\footnote{Voir les travaux de \cite{Smith93},
\cite{Canova94}, \cite{DejongIngramWhiteman96}, \cite{Geweke99},
\cite{DridiGuayRenault07} et \cite{Bierens05}.}. Dans cette
perspective, l'idée suivant laquelle de tels modèles sont utiles
pour la prévision et l'analyse de la politique économique s'est
répandue dans le milieu académique ainsi qu'auprès des institutions
internationales et des banques centrales. Parmi toutes ces approches
économétriques, la littérature privilégie, pour de <<  bonnes >>
et <<  mauvaises >> raisons, la statistique bayésienne. Parmi les
<<  bonnes >> raisons, nous pourrions souligner le fait que la
fonction de vraisemblance d'un modèle de dimension élevée (de
nombreux paramètres à estimer) est souvent <<  plate >> dans
certaines directions. En d'autres termes, les données peuvent ne pas être
suffisamment informatives pour identifier (avec précision) les
paramètres structurels. En déformant la fonction de vraisemblance à l'aide
d'informations \textit{a priori} sur les paramètres, c'est-à-dire en
privilégiant une approche bayésienne, l'identification devient
possible. Il est néanmoins trop souvent ignoré que la mise en oeuvre
et l'interprétation des résultats de l'estimation bayésienne
requièrent un certain nombre d'hypothèses et de conditions de
validité, ou que nombre de problèmes rencontrés en économétrie
classique ont leur contrepartie en économétrie bayésienne. Toujours
est-il que l'approche bayésienne a considérablement favorisé le
développement des modèles DSGE comme outil d'analyse et de prévision
de la politique monétaire. Dans le même temps, il n'en demeure pas
moins que les modèles DSGE et VAR continuent à être opposés et que
nombre de papiers cherchent généralement à légitimer leurs résultats
en comparant, par exemple, les prévisions (ou tout autre statistique
ou quantité d'intérêt) de leur(s) modèle(s) avec ceux d'un VAR
(\cite{SmetsWouters2004}).}\newline

\par{L'objet de ce papier est de présenter l'approche bayésienne des modèles VAR et DSGE en mettant en avant les
principaux concepts, leur mise en oeuvre pratique et les limites
sous-jacentes. Nous montrons en quoi les modèles DSGE et VAR sont
des outils complémentaires que l'on ne doit pas nécessairement
chercher à opposer. Nous n'abordons pas ici certains problèmes
importants, comme l'estimation non linéaire des modèles
DSGE\footnote{Le lecteur intéressé trouvera une introduction
intéressante pour l'estimation des modèles nonlinéaires dans
\cite{AndrieuDoucetRobert04}, \cite{Arulampalam02}, et
\cite{Andrieu04}, ainsi que dans les contributions de
\cite{Gordon93} et \cite{Kitagawa96}. Pour des applications en
économie, voir \cite{ChopinPelgrin04}, \cite{VillaverdeRamirez05}
et \cite{AnSchorfheide07}.}}.\newline

\par{L'article est organisé comme suit. Dans une première section, nous présentons les principaux
concepts de l'analyse bayésienne et montrons comment les appliquer dans le cadre des modèles VAR. Une
attention particulière est attachée à la nature (informative, non
informative, empirique) des croyances \textit{a priori}. Dans une deuxième
section, nous abordons les spécificités de l'approche bayésienne des
modèles DSGE. Contrairement aux modèles VAR, il n'est
plus possible d'obtenir une expression analytique de la distribution
\textit{a posteriori}. Pour remédier à cette difficulté,
il est nécessaire de recourir à des méthodes de Monte-Carlo et
notamment à la théorie des chaînes de Markov. Dans cette
perspective, après avoir dérivé de manière générale la densité a
posteriori d'un modèle DSGE, nous expliquons les principaux
algorithmes d'estimation (algorithme de Metropolis-Hasting, par
fonction d'importance). Dans une troisième section, nous illustrons comment peuvent se combiner les
approches VAR et DSGE.}

\section{L'approche Bayésienne}\label{sec:1}

\subsection{Généralités}\label{sec:1:1}

\par{L'approche bayésienne propose un cadre rigoureux pour (\emph{i})
formaliser nos croyances\footnote{Le mot <<  croyance >> suggère
une dimension subjective de l'information \textit{a priori}. Il
convient néanmoins de signaler que parmi les
économètres bayésiens il n'y a pas de consensus sur l'interprétation
subjective ou objective des probabilités. Par exemple, l'approche
bayésienne empirique utilise l'échantillon pour définir l'information
\textit{a priori} (voir par exemple le prior Minnesota de la section
\ref{sec:1:4:pratique}).} \textit{a priori} et (\emph{ii})
déterminer comment celles-ci doivent être mises à jour une fois
que les données sont observées. Les croyances, \textit{a priori} ou
\textit{a posteriori}, sont représentées à l'aide
d'une densité de probabilité jointe sur les paramètres d'un modèle.
Cette densité jointe caractérise l'incertitude quant au
processus générateur des données (DGP, pour \textit{Data Generating
 Process}), en décrivant une famille (un continuum) de modèles.}\newline

\par{Imaginons que nous souhaitions caractériser nos croyances
sur le paramètre de Calvo d'une courbe de Phillips. Ce paramètre, $\xi_p$, est la probabilité
pour une firme, en concurrence monopolistique, de ne pas pouvoir
ajuster son prix de façon optimale à une date quelconque. Ainsi,
nous savons déjà que ce paramètre doit appartenir à l'intervalle
[0,1]. Nous pourrions donc utiliser une distribution
bêta\footnote{Pour une présentation des principales distributions
utilisées dans cette littérature, le lecteur peut se reporter aux
annexes de \cite{ZellnerBook}. La distribution uniforme est un cas
particulier de la bêta.} définie sur cet intervalle. A partir de la
probabilité $\xi_p$, nous pouvons définir le temps moyen pendant
lequel une firme ne pourra pas ajuster son prix de façon optimale~:
$\zeta_p\equiv\frac{1}{1-\xi_p}$. Si par ailleurs, à l'aide
d'enquêtes microéconomiques, nous savons que le temps moyen
durant lequel une firme ne réajuste pas son prix de façon
<<  optimale >> est de quatre trimestres, nous pouvons déduire qu'une
valeur pertinente de la probabilité $\xi_p$ est trois quarts.
L'économiste bayésien pourra donc formaliser son \textit{a priori}
sur le paramètre $\xi_p$ en sélectionnant une distribution bêta
ayant pour mode $3/4$ et en spécifiant une variance mesurant son
incertitude sur le paramètre d'intérêt. Il choisira une variance
d'autant plus grande qu'il est incertain des évaluations
microéconomiques dont il dispose\footnote{Cette incertitude peut
s'expliquer par une adéquation imparfaite entre le concept théorique
et l'enquête microéconomique.}. Notons qu'il pourrait directement
poser son \textit{a priori} sur le délai moyen d'attente avant qu'il
ne puisse optimalement ajuster son prix, $\zeta_p$~;
ceci résultera en une distribution différente pour le paramètre
$\xi_p$. Si $\xi_p$ est le seul paramètre du modèle pour lequel nous sommes incertain,
\textit{ie} si les autres paramètres ont des variances \textit{a priori} nulles, la densité
a priori sur ce paramètre décrit une famille de DGP, indexée par $\xi_p$~:
chaque valeur possible de $\xi_p$ correspond à un DGP.}\newline

\par{Plus généralement, nous noterons l'\textit{a priori} sur un vecteur
de paramètres $\theta_{\mathcal M}$ associé à un modèle paramétrique
$\mathcal M$, $\theta_{\mathcal M} \equiv \left(\theta_1^{\mathcal
M},\dots,\theta_{q_{\mathcal M}}^{\mathcal M}\right)$, de la façon
suivante~:
\begin{equation}\label{equ:bayes:prior}
p_0(\theta_{\mathcal M}|\mathcal M)
\end{equation}
Cette densité jointe définit notre incertitude quant aux paramètres $\theta_{\mathcal M}$
avant que nous ayons porté attention aux données. Il convient de noter
que nous raisonnons conditionnellement à un modèle. En toute
généralité l'incertitude pourrait aussi porter sur la forme du modèle
paramétrique $\mathcal M$. Plus loin nous omettrons généralement le conditionnement (ainsi que l'indexation)
par le modèle pour simplifier les notations.}\newline

\par{Nous observons un échantillon $\sample = \left\{y_t^*
\right\}_{t=1}^T$ où $y_t^*$ est un vecteur de $m$ variables. Nous
nous limiterons au cas où l'indice $t$ représente le temps. La
vraisemblance est la densité de l'échantillon conditionnellement au
modèle et ses paramètres~; on notera~:
\begin{equation}\label{equ:bayes:lik}
\mathcal L \left(\theta_{\mathcal M};\sample,\mathcal M\right)
\equiv p(\sample| \theta_{\mathcal M},\mathcal M)
\end{equation}
L'estimateur du maximum de vraisemblance (MV) des paramètres
$\theta_{\mathcal M}$ d'un modèle $\mathcal M$ est la valeur des
paramètres qui rend le plus <<  probable>> l'occurrence de
l'échantillon à notre disposition. Autrement dit, l'estimateur du
MV sélectionne le paramètre $\theta_{\mathcal M}$ définissant le
DGP qui a le plus probablement généré les données. La démarche
statistique, classique ou bayésienne, est une démarche d'inversion
--- il s'agit de remonter des observations aux paramètres du DGP.
Un modèle (la vraisemblance) définit la densité d'un ensemble de
variables conditionnellement à des paramètres inconnus.
L'observation de l'échantillon donne en retour de l'information
sur les paramètres. La notation définie par l'équation
(\ref{equ:bayes:lik}) résume le principe de l'inférence~; la
vraisemblance est la densité de l'échantillon $\sample$ sachant
les paramètres $\theta$ mais nous écrivons habituellement la
vraisemblance comme une fonction des paramètres, \textit{ie.}
formellement nous échangeons les rôles de $\sample$ et
$\theta$.}\newline

\par{Nous disposons des densités $p_0(\theta_{\mathcal M}|\mathcal M)$, qui caractérise
l'information postulée \textit{a priori}, et $p(\sample
|\theta_{\mathcal M},\mathcal M)$, qui caractérise l'information amenée par les données. On croise
ces deux sources d'informations orthogonales, en utilisant le théorème de
Bayes, pour obtenir la densité de $\theta_M$ sachant les données
$\sample$, \textit{ie} la densité postérieure~:
\begin{equation}\label{equ:bayes:posterior}
p_1\left(\theta_{\mathcal M}|\sample,\mathcal M\right) =
\frac{p_0\left(\theta_{\mathcal M} | \mathcal M\right)
p(\sample|\theta_{\mathcal M},\mathcal M)}{p(\sample | \mathcal M)}
\end{equation}
avec
\begin{equation}\label{equ:bayes:marginaldensity}
p\left(\sample | \mathcal M\right) = \int_{\Theta_{\mathcal M}}
p_0\left(\theta_{\mathcal M} | \mathcal M\right)
p(\sample|\theta_{\mathcal M},\mathcal M)\mathrm d \theta_{\mathcal
M}
\end{equation}
la densité marginale. Ainsi, la densité postérieure est
proportionnelle à la densité \textit{a priori} multipliée par la
vraisemblance~:
\[
\begin{split}
p_1\left(\theta_{\mathcal M}|\sample,\mathcal M\right) & \propto
p_0\left(\theta_{\mathcal M} | \mathcal M\right)
p\left(\sample|\theta_{\mathcal M},\mathcal M\right)\\
& \equiv \mathcal K \left(\theta_{\mathcal M} | \sample,{\mathcal M}
\right)
\end{split}
\]
Puisque le dénominateur dans \equaref{equ:bayes:posterior}, la
densité marginale, ne dépend pas de $\theta_{\mathcal M}$,
l'inférence sur les paramètres, par exemple l'évaluation de
l'espérance postérieure, peut être mise en oeuvre à l'aide du seul
noyau postérieur, $\mathcal K \left(\theta_{\mathcal M} |
\sample,{\mathcal M} \right)$. On représente nos croyances \textit{a
posteriori} en exhibant les propriétés de la distribution \textit{a
posteriori}. Nous pouvons représenter graphiquement la densité
postérieure marginale de chaque paramètre $\theta$, construire des
intervalles contenant $\alpha \%$ de la distribution postérieure, ou
encore calculer des moments \textit{a posteriori}. Par exemple, la
comparaison des variances \textit{a priori} et \textit{a posteriori}
peut nous renseigner sur l'information apportée par les données,
relativement à celle contenue dans nos croyances \textit{a priori}.
Les variances de chaque paramètre sont définies à partir des
éléments diagonaux des matrice suivantes~:
\[
\mathbb V_0[\theta] = \int_{\Theta} \theta\theta'p_0(\theta)\mathrm d\theta -
\left(\int_{\Theta} \theta p_0(\theta)\mathrm d\theta\right)\left(\int_{\Theta} \theta p_0(\theta)\mathrm d\theta\right)'
\]
et
\[
\mathbb V_1[\theta] = \int_{\Theta} \theta\theta'p_1(\theta|\sample)\mathrm d\theta -
\left(\int_{\Theta} \theta p_1(\theta|\sample)\mathrm d\theta\right)\left(\int_{\Theta} \theta p_1(\theta|\sample)\mathrm d\theta\right)'
\]
Si la variance postérieure d'un paramètre est plus faible que sa
variance \textit{a priori} cela signifie que les données apportent
une information supplémentaire sur ce paramètre, relativement à
l'information \textit{a priori}. Dans certains cas, il est possible
d'obtenir analytiquement la densité postérieure et ses
moments\footnote{Une condition nécessaire est que nous disposions
d'une expression analytique de la vraisemblance.}~; nous verrons un
exemple dans la section suivante. Plus généralement, il est
nécessaire de recourir à des algorithmes numériques, pour
caractériser la distribution postérieure, \textit{ie} pour évaluer
les intégrales nécessaires au calcul des moments.}\newline

\par{Pour communiquer nos croyances \textit{a posteriori} on désire souvent recourir à un média plus synthétique en
résumant, à l'image de l'approche classique, la distribution postérieure par un point. On parle alors d'estimation
ponctuelle. Réduire la distribution postérieure à un point s'apparente à un choix en univers incertain. Il est donc
naturel de construire une estimation ponctuelle en minimisant l'espérance postérieure d'une fonction de perte~:
\begin{equation}\label{equ:bayes:pointestimate}
\hat{\theta} = \arg\underset{a}{\min}\int_{\Theta} p_1(\theta|\sample,\mathcal M)L(a,\theta)\mathrm d\theta
\end{equation}
où $L(a,\theta)$ est une fonction associant une perte au choix $a$ si la vraie valeur du paramètre est $\theta$. Si,
par exemple, la fonction de perte est quadratique\footnote{Nous supposons un instant qu'il n'y a qu'un paramètre
dans le modèle.}~:
\[
L(a,\theta) = (a-\theta)^2
\]
alors on montre \cite[page 24]{ZellnerBook} que l'estimation ponctuelle doit être l'espérance postérieure de
$\theta$. D'autres fonctions de perte aboutiront à d'autres estimations ponctuelles. La médiane postérieure
peut être rationnalisée en considérant la fonction de perte $L(a,\theta) = |a-\theta|$~; plus généralement on
peut exhiber une fonction de perte pour justifier l'utilisation de chaque quantile de la distribution comme une
estimation ponctuelle.}\newline

\par{Tant que l'inférence porte sur les paramètres d'un modèle, nous
pouvons mettre de côté la constante d'intégration, $p(\mathcal
Y_T^{\ast}|\mathcal M)$. Cependant, la densité marginale
con\-tient une information pertinente si nous désirons comparer
différents modèles. En effet, celle-ci nous renseigne sur la
densité de l'échantillon conditionnellement au modèle.
L'interprétation de sa définition
(\ref{equ:bayes:marginaldensity}) est directe~: la densité
marginale est une moyenne des vraisemblances, obtenues pour
différentes valeurs des paramètres, pondérées par nos croyances
\textit{a priori} sur les valeurs de ces paramètres. Comme cette
quantité ne dépend pas des paramètres, puis\-que nous avons
intégré le noyau postérieur par rapport aux paramètres pour
l'obtenir, elle autorise facilement la comparaison de modèles non
emboîtés\footnote{La comparaison est moins simple dans le cas de
l'inférence classique.}. Par exemple, si nous disposons de deux
modèles alternatifs, $\mathcal A$ et $\mathcal B$, pour expliquer
les données et si nous obtenons $p(\mathcal Y_T^{\ast}|\mathcal
A)>p(\mathcal Y_T^{\ast}|\mathcal B)$, alors cela signifie que les
données à notre disposition ont plus probablement été générées par
le modèle $\mathcal A$ que par le modèle $\mathcal B$. Cette
approche ne fournit pas un test formel, pour prendre une décision
il faut spécifier un critère de perte\footnote{Se reporter à
Zellner \cite[chapitre 10]{ZellnerBook}, en particulier la
première section pages 292 à 298. Le choix d'un modèle parmi une
collection de modèles s'apparente à l'estimation d'un paramètre
dont la distribution est discrète. Il y a donc une analogie entre
le choix d'un modèle et l'estimation ponctuelle de $\theta$.}.
Ici, nous supposons implicitement que nous n'avons pas de << 
préférence >> \textit{a priori} pour un des deux modèles. En
toute généralité nous pourrions définir une densité de probabilité
(discrète) \textit{a priori} pour les modèles $\mathcal I =
\mathcal A,\mathcal B$. Par exemple, nous pourrions supposer que
$p_0(\mathcal A)>p_0(\mathcal B)$ \textit{a priori}. Par le
théorème de Bayes, nous obtenons la probabilité \textit{a
posteriori} du modèle $\mathcal I$~:
\[
p_1(\mathcal I|\sample) = \frac{p_0(\mathcal I)p(\sample|\mathcal
  I)}{\sum_{\mathcal I} p_0(\mathcal I)p(\sample|\mathcal I)}
\]
L'interprétation est directe, il s'agit d'une version discrète des
équations (\ref{equ:bayes:posterior}) et
(\ref{equ:bayes:predictivedensity}). Si nous devons sélectionner
un modèle, notre préférence ira au modèle qui maximise la densité
postérieure. L'élicitation d'une densité de probabilité \textit{a
priori} sur la collection de modèles n'est pas une chose aisée~;
on ne peut associer une probabilité à un modèle de la même façon
que l'on pose une densité \textit{a priori} sur le paramètre de
Calvo d'une courbe de Phillips. C'est pourquoi une densité de
probabilité uniforme est souvent utilisée dans la littérature
concernée par l'estimation des modèles DSGE. Pour une revue des
enjeux de la comparaison de modèles, on peut lire \citet[chapitre
7]{RobertBook} ou \citet{Sims2003}. Enfin, notons que la comparaison de différents
modèles, par l'intermédiaire de leurs densités marginales, ne doit
pas nécessairement aboutir à un choix. Dans certaines situations,
il peut être souhaitable de combiner plusieurs modèles, en les
pondérant par leurs densités marginales respectives\footnote{Pour
une description des méthodes \emph{Bayesian Model Averaging}, voir
\cite{BookKoop03}, chapitre 11.}.}\newline

\par{Après l'estimation, le modèle peut être utilisé pour construire des prévisions et des fonctions
de réponse. À l'image de l'estimation des paramètres, le paradigme bayésien ne fournit pas en premier lieu
des prévisions ponctuelles mais des densités prédictives. Supposons que nous souhaitions établir des
prédictions sur $y_{T+1}^*$, un vecteur d'endogènes non encore observées, le but est de construire la
densité (prédictive) de ce vecteur. Cette densité peut être obtenue en intégrant par rapport à $\theta$
la densité postérieure jointe de $y_{T+1}^*$ et $\theta$~:
\[
p\left(y_{T+1}^*|\sample\right) = \int_{\Theta}p\left(y_{T+1}^*,\theta|\sample\right)\mathrm d \theta
\]
où la densité sous l'intégrale est définie par~:
\[
p\left(y_{T+1}^*,\theta|\sample\right) = p\left(y_{T+1}^*|\theta,\sample\right)p_1\left(\theta|\sample\right)
\]
par le théorème de Bayes. La densité jointe de $y_{T+1}^*$ et $\theta$ est le produit de la densité
postérieure de $\theta$ et de la densité de $y_{T+1}^*$ conditionnelle à $\theta$. Cette dernière
densité est directement obtenue à partir de la définition du
modèle.\newline Donnons un exemple dans le cas scalaire. Si le
modèle est un processus autorégressif d'ordre un~: $y^*_t = \theta  y^*_{t-1} + \varepsilon_t$ avec $t=1,\dots,T$,
 $\varepsilon_t \underset{iid}{\sim}\normal{0}{\sigma^2}$ et $\sigma^2$, la variance de l'innovation, connue.
La distribution de $y_{T+1}^*$ conditionnellement à $\theta$ et
$\sample$\footnote{Dans le cas du modèle AR(1) l'information
apportée par l'échantillon est résumée par la dernière observation
$y^*_{T}$.} est gaussienne~: $y_{T+1}|
\sample,\theta\sim\normal{\theta y^*_T}{\sigma^2}$. La densité
prédictive s'écrit finalement~:
\begin{equation}\label{equ:bayes:predictivedensity}
p\left(y_{T+1}^*|\sample\right) = \int_{\Theta} p\left(y_{T+1}^*|\theta,\sample\right)
p_1\left(\theta|\sample\right)\mathrm d \theta
\end{equation}
et s'interpréte comme une moyenne des densités conditionnelles de $y_{T+1}^*$ sachant $\theta$, pondérées par la
densité postérieure de $\theta$.\newline À partir de cette densité prédictive, on peut construire une prédiction
ponctuelle des variables en se donnant une fonction de perte, représenter un intervalle contenant $\alpha\%$ de
la distribution de $y_{T+1}^*$, ou encore représenter la densité prédictive (marginale) de chaque variable
composant le vecteur des observables. En confrontant la densité prédictive aux réalisations effectives des
variables, on peut alors évaluer dans quelle mesure notre modèle tend à sur-estimer ou sous-estimer, par exemple, le
taux de croissance à un trimestre du PIB par tête. Cette comparaison peut fournir un critère d'évaluation du modèle.
Si on se rend compte que les réalisations effectives d'une variable se situent systématiquement dans les queues de la
densité prédictive, alors on peut conclure que le modèle est mal spécifié vis-à-vis de cette variable.}

\subsection{Le choix des croyances \textit{a priori}}\label{sec:1:2}

\par{On comprend déjà que le choix des croyances \textit{a priori} est essentiel, dans la mesure où
il détermine partiellement les résultats (surtout pour un échantillon de taille réduite comme
nous le verrons par la suite). La subjectivité de l'économètre ne peut
intervenir que dans la première étape d'élicitation de l'\textit{a priori}, les étapes (l'évaluation
de la vraisemblance,...) suivantes sont automatiques et nécessairement
objectives. La question du choix des croyances \textit{a priori} est donc crucial,
d'autant plus qu'il s'agit du point qui engendre le plus de méfiance parmi les économistes, de
formation -- très généralement -- classique. Il est donc important de bien comprendre le rôle de la
densité \textit{a priori} dans les résultats, par exemple en menant des exercices de sensibilité aux
croyances \textit{a priori}. Ces expériences, en donnant une idée du rôle des priors, dévoilent implicitement la forme
de la vraisemblance. L'expérience la plus extrème\footnote{Une expérience moins extrême serait de considérer des
densités \textit{a priori} plus générales. Supposons que notre \textit{a priori} sur un
paramètre $\mu$ soit caractérisé par une loi normale centrée en $\mu_0$ et de variance $\sigma_0^2$. Nous pourrions
évaluer la sensibilité des résultats à ce choix en reprenant l'estimation avec une densité \textit{a priori} de student~:
\[
  p_0(\mu) \propto \left(\nu s+(x-\mu_0)^2\right)^{-\frac{\nu+1}{2}}
\]
L'espérance \textit{a priori} serait alors $\mu_0$ mais la variance \textit{a priori} serait $\frac{\nu}{\nu-2}s$
(pour $\nu$ strictement supérieur à 2). En faisant varier le nombre de degré de liberté $\nu$ on s'écarte ou se
rapproche du prior gaussien.} est de considérer un \textit{a priori} non informatif, c'est-à-dire
le cas où nous n'avons aucune croyance \textit{a priori} sur la
paramétrisation du modèle. De façon assez surprenante,
les statisticiens bayésiens ne parviennent pas à s'accorder sur une chose aussi essentielle
que la caractérisation du non savoir.}\newline

\par{Dans la section \ref{sec:1:1} nous avons examiné le cas d'un \textit{a priori} informatif sur le paramètre
de Calvo définissant le degré de rigidité de l'inflation. Dans ce cas notre connaissance \textit{a priori}
provient de l'observation de données microéconomiques, différentes de celles utilisées pour l'estimation
du modèle. Lorsque l'information \textit{a priori} est basée sur des données, celles-ci doivent être différentes
des données utilisées pour identifier le modèle. Dans le cas contraire la démarcation entre vraisemblance et densité
\textit{a priori} devient plus ambiguë, ce qui paraît inacceptable pour de nombreux statisticiens. Notons néanmoins que
de non moins nombreux statisticiens utilisent l'échantillon pour définir les croyances \textit{a priori}. Par exemple
quand il s'agit de spécifier la densité \textit{a priori} de façon à optimiser les capacités prédictives d'un
modèle (voir plus loin la section \ref{sec:1:4:pratique}). Les croyances \textit{a priori} peuvent aussi être basées
sur des considérations purement théoriques (voir la section \ref{sec:bvardsge}). Dans la littérature concernée par
l'estimation des modèles DSGE (et aussi des VAR), les croyances \textit{a priori}, indépendamment de l'origine de ces
croyances, sont généralement représentées par des densité paramétrées (distribution gaussienne, gamma,...). Dans
certains cas, on parle alors d'\textit{a priori} conjugués, elles sont choisies de façon que la densité
\textit{a posteriori} soit de la même famille paramètrique (voir l'exemple du modèle VAR plus loin). La motivation est
essentiellement technique, l'utilisation de formulations conjuguées ou plus généralement paramètrées des
\textit{a priori} est un héritage du passé. Aujourd'hui, la technologie ne nous empêche pas d'adopter une formulation
non paramétrique plus générale. Par exemple nous pourrions caractériser nos croyances \textit{a priori} sur chaque
paramètre en spécifiant les quantiles de chaque distribution. Il est vrai que nos croyances sont rarement aussi
précises.}\newline

\par{Dans certaines situations nos connaissances \textit{a priori} sont faibles ou nous désirons faire comme si c'était
le cas (par exemple pour évaluer la sensibilité des résultats aux \textit{a priori}). Malheureusement la caractérisation
de l'ignorance est toujours sujet à débat. Un exemple frappant est donné par \citet{SimsUhlig91} puis
\citet{PhillipsJAE91a}, \citet{PhillipsJAE91b} et \citet{SimsJAE91}, qui débatent de la caractérisation de
l'ignorance dans un modèle autorégressif d'ordre un et des
conséquences sur la détection de racines unitaires.\newline Une première approche est de considérer un prior plat. Pour
un paramètre $\mu$ qui peut prendre des valeurs entre $-\infty$ et $\infty$, \citet{JeffreyBook61} propose d'adopter une
distribution uniforme entre $-\infty$ et $\infty$~:
\[
p_0(\mu) \propto \mathrm 1
\]
Évidemment cette densité est impropre dans le sens où $\int p_0(\mu)\mathrm d\mu$ est indéfini.
Mais c'est précisément cette propriété qui, pour \citeauthor{JeffreyBook61}, rend ce prior non informatif.
En effet, pour tout $a<b<c<d$ on ne peut pas dire que $\mu \in [a,b]$ soit \textit{a priori} plus probable que
$\mu \in [c,d]$, puisque les probabilités de ces événements sont nulles. Pour un paramètre $\sigma$, par exemple
un écart-type, qui peut prendre des valeurs entre 0 et $\infty$, \citeauthor{JeffreyBook61} propose d'adopter
une distribution uniforme pour le logarithme de $\sigma$ entre $-\infty$ et $\infty$~:
\[
\begin{split}
p_0(\log \sigma) &\propto 1\\
\Leftrightarrow p_0(\sigma) &\propto \frac{1}{\sigma}
\end{split}
\]
Comme dans le cas précédent l'intégrale de cette densité est
impropre. En particulier, on ne peut définir $\int_0^c
p_0(\sigma)\mathrm d \sigma$ et $\int_c^{\infty} p_0(\sigma)\mathrm
d \sigma$, nous ne pouvons dire s'il est plus probable que $\sigma$
soit supérieur ou inférieur à $c$\footnote{Cette propriété est
indispensable pour représenter l'ignorance. Dans la littérature
DSGE, la distribution inverse-gamma avec un moment d'ordre deux
infini est souvent utilisée pour représenter le peu d'information
dont nous disposons sur la variance des chocs structurels (voir par
exemple \citet{SmetsWouters2004}). Cette distribution est
informative dans le sens où, même si le moment d'ordre deux n'est
pas défini, il est possible de comparer les probabilités qu'une
variance soit supérieure ou inférieure à $c>0$.}. On note en passant
que cette densité a l'heureuse propriété d'être invariante à une
transformation puissance\footnote{Si le logarithme de $\sigma$ est
uniformément distribué sur $]-\infty,\infty[$ alors le logarithme de
$\sigma^a$ (avec $a>0$) est aussi uniformément distribué sur
$]-\infty,\infty[$ car $\log \sigma^a = a\log\sigma$.}, si le prior
est non informatif sur l'écart type, il en va de même pour la
variance ($\sigma^2$).\newline Plus tard, Jeffrey généralisa ce
résultat d'invariance et proposa un prior non informatif (le plus
souvent impropre) plus général basé sur la matrice d'information de
Fisher~:
\[
p_0(\theta) \propto |I(\theta)|^{\frac{1}{2}}
\]
avec
\[
I(\theta) = \mathbb E \left[\left(\frac{\partial p(\sample|\theta)}{\partial \theta}\right)
\left(\frac{\partial p(\sample|\theta)}{\partial \theta}\right)'\right]
\]
La matrice d'information de Fisher quantifie l'information amenée par le modèle et les données sur
le paramètre $\theta$. En favorisant les valeurs de $\theta$ pour lesquelles l'information de Fisher
est plus grande, on diminue l'influence de la loi \textit{a priori} puisque l'information véhiculée
par celle-ci est peu différente de l'information provenant de la vraisemblance. La définition de la densité \textit{a
priori} est donc liée à la courbure de la vraisemblance. Cette densité \textit{a priori} est invariante à toute
reparamétrisation (continue) du modèle (voir \citet[annexe du chapitre 2]{ZellnerBook} pour une description plus
détaillée des propriétés d'invariance).}\newline

\par{L'utilisation d'un prior plat ou d'un prior dérivé de la matrice
d'information de Fisher pour caractériser l'absence d'information affecte généralement l'inférence. Par exemple, dans
un modèle AR(1), voir \citet{PhillipsJAE91a}, un prior basé sur l'information de Fisher n'est pas équivalent à
un prior uniforme (plat). En effet, dans un modèle dynamique, la quantité d'information véhiculée par les données
(\textit{ie} la vraisemblance) dépend de la valeur du paramètre autorégressif ($\rho$). Si le paramètre est proche de
l'unité, voire égal ou supérieur à un, les données sont plus informatives. Ainsi, pour \citeauthor{PhillipsJAE91a},
l'utilisation d'un prior plat, à l'instar de \citet{ZellnerBook} ou \citet{SimsUhlig91}, biaise la distribution
postérieure de $\rho$ en faveur de la stationnarité. En donnant autant de poids aux valeurs explosives de $\rho$
qu'aux valeurs stationnaires, le prior plat ne prend pas en compte le fait que des données générées par un modèle
à racine unitaire ou explosif sont plus informatives. Il existe d'autres approches pour caractériser l'ignorance, on
peut lire le chapitre 3 de \citet{RobertBook} et plus spécialement la section 5.}\newline

\par{Le choix d'une densité \textit{a priori} et ses conséquences sur l'inférence sont l'objet de toutes les
critiques de la part des statisticiens ou économistes classiques. Il ne faudrait pourtant pas oublier que
le paradigme classique n'est pas plus exempt de choix aux conséquences non négligeables sur l'inférence. Par exemple,
le choix d'une métrique (minimiser la somme des carrés des résidus ou la somme des valeurs absolues des résidus), le
choix des variables instrumentales, modèles auxiliaires ou des conditions de moments, sont rarement discutés même s'ils
déterminent les résultats. Dans une certaine mesure, nous n'avons même plus conscience des choix effectués. L'approche
bayésienne est de ce point de vu bien plus transparente.}\newline


\subsection{Comportement asymptotique et approximations} \label{sec:1:3}

\par{Même si l'approche bayésienne ne repose pas sur des arguments asymptotiques, comme généralement
l'approche classique, il est utile de s'interroger sur le comportement asymptotique de ses estimateurs.
Le résultat rassurant est que si les conditions de normalité asymptotique de l'estimateur du maximum de
vraisemblance sont réunies\footnote{Se reporter à \citet[chapitre 7]{GourierouxMonfortStatBook1}.}, alors
la distribution postérieure tend vers une gaussienne multivariée. Asymptotiquement, la distribution
postérieure est centrée sur l'estimateur du maximum de vraisemblance. Ce résultat, avancé par Laplace,
est intuitif puisque lorsque la taille de l'échantillon tend vers l'infini, le poids de l'information \textit{a priori}
relativement à l'information contenue dans l'échantillon devient marginal.}\newline

\par{Plus formellement, si on note $\theta^*$ l'unique mode de la distribution postérieure obtenu en maximisant
le noyau postérieur $\mathcal{K}(\theta) \equiv \mathcal K \left(\theta_{\mathcal A} | \sample,{\mathcal A}
\right)$, et s'il est possible d'écrire une approximation de Taylor à l'ordre deux du noyau postérieur autour
de $\theta^*$, alors nous avons~:
\[
\begin{split}
\log \mathcal{K}(\theta) = \log \mathcal{K}(\theta^*)
&+ (\theta-\theta^*)'\left.\frac{\partial \log \mathcal{K}(\theta) }{\partial \theta}\right|_{\theta = \theta^*}\\
&+ \frac{1}{2}(\theta-\theta^*)'\left.\frac{\partial^2 \log \mathcal{K}(\theta) }
{\partial \theta\partial \theta'}\right|_{\theta = \theta^*}(\theta-\theta^*) + \mathcal O (||\theta-\theta^*||^3)
\end{split}
\]
Puisque les dérivées premières sont, par définition, nulles en $\theta^*$, nous avons de façon équivalente~:
\[
\log \mathcal{K}(\theta) = \log \mathcal{K}(\theta^*)
- \frac{1}{2}(\theta-\theta^*)'[\mathcal H(\theta^*)]^{-1}(\theta-\theta^*) + \mathcal O (||\theta-\theta^*||^3)
\]
où $\mathcal H(\theta^*)$ est l'opposé de l'inverse de la matrice hessienne évaluée au mode. Ainsi, en ne considérant
que le terme quadratique, le noyau postérieur peut être approximé par~:
\[
\mathcal{K}(\theta)\ \dot{=}\ \mathcal{K}(\theta^*)e^{-\frac{1}{2}(\theta-\theta^*)'[\mathcal H(\theta^*)]^{-1}(\theta-\theta^*)}
\]
on reconnaît, à une constante d'intégration près\footnote{Notons $c$ cette constante d'intégration, c'est-à-dire la
constante telle que $\int c^{-1}\mathcal K (\theta)\mathrm d \theta = 1$. Cette constante (voir les équations
(\ref{equ:bayes:posterior}) et (\ref{equ:bayes:marginaldensity})) est une approximation de la densité marginale,
$p(\sample)$. Par définition de la densité d'une loi normale, on a~:
\[
c = \mathcal{K}(\theta^*)(2\pi)^{\frac{q}{2}}
|\mathcal H(\theta^*)|^{\frac{1}{2}}
\]
On dit que $c$ est l'approximation de Laplace de la densité marginale.
L'erreur d'approximation est d'ordre $\mathcal O (T^{-1})$.} la densité d'une loi normale multivariée. En complétant pour la
constante d'intégration, nous obtenons finalement une approximation de la densité postérieure $p_1(\theta)\equiv p_1\left(\theta_{\mathcal A}|\sample,\mathcal A\right)$~:
\begin{equation}
p_1\left(\theta\right)\ \dot{=}\ (2\pi)^{-\frac{q}{2}}|\mathcal H(\theta^*)|^{-\frac{1}{2}}
e^{-\frac{1}{2}(\theta-\theta^*)'[\mathcal H(\theta^*)]^{-1}(\theta-\theta^*)}
\end{equation}
Généralement, la matrice hessienne est d'ordre $\mathcal O(T)$,
lorsque la taille de l'échantillon augmente la distribution
postérieure se concentre autour du mode. \`{A} partir de cette
approximation asymptotique on peut alors très facilement calculer,
par exemple, des moments postérieurs ou approximer les densités
prédictives. Par exemple, l'espérance postérieure de $\varphi
(\theta)$ est définie par~:
\[
\mathbb E\left[\varphi (\theta)\right] = \frac{\int_{\Theta}\varphi(\theta)p(\sample|\theta)p_0(\theta)\mathrm d\theta}
{\int_{\Theta}p(\sample|\theta)p_0(\theta)\mathrm d\theta}
\]
\citet{TierneyKadane1986} montrent que si l'on approxime à l'ordre
deux le numérateur autour du mode de
$\varphi(\theta)p(\sample|\theta)p_0(\theta)$ et le dénominateur
autour du mode de $p(\sample|\theta)p_0(\theta)$, alors l'erreur
d'approximation de l'espérance est d'ordre $\mathcal O (T^{-2})$.
Les erreurs d'approximation du numérateur et du dénominateur, qui
sont d'ordre $\mathcal O(T^{-1})$, se compensent favorablement.
L'approche de \citeauthor{TierneyKadane1986} ne va pas sans poser certains problèmes.
Si on cherche à calculer $\mathbb E\left[\varphi (\theta)\right]$ pour différentes
fonctions $\varphi$, alors il est nécessaire de recourir à une nouvelle maximisation pour
chaque paramètre et chaque fonction $\varphi$. Par exemple, si on
désire calculer les espérances et écart-types \emph{a posteriori}
pour chacun des $k$ paramètres, il faut recourir à $2k+1$ maximisations, auxquelles il faut rajouter le calcul des
matrices hessiennes. Il est alors évident que si $k$ est élevé, une
telle approximation peut devenir coûteuse en temps de calculs.
\cite{TierneyKassKadane89} propose différentes méthodes pour
pallier cette difficulté\footnote{Pour plus de détails, voir
\cite{CarlinLouis00}, \cite{Poirier95}, et
\cite{TierneyKadane1986}.}. Notons néanmoins qu'une approche basée sur des simulations
(voir la section \ref{sec:dsge}) devient aussi plus coûteuse lorsque le nombre de paramètres
augmente.}


\subsection{Un modèle linéaire : le modèle VAR}\label{sec:1:4}

\par{Dans cette section, nous considérons un exemple où les résultats
peuvent être obtenus analytiquement. Le modèle VAR gaussien se
prête, comme tout modèle linéaire gaussien, à cet exercice et a l'avantage d'être un outil
couramment utilisé en macroéconomie (voir par exemple la contribution de
Fabrice Collard et Patrick Fève dans ce numéro).}\newline

\par{Nous considérons un modèle VAR($p$) pour caractériser le vecteur $1
\times m$ de variables endogènes $y^*_t$ observées~:
\[
y^*_t = \sum_{i=1}^p y^*_{t-i} \mathbf{A}_i + \varepsilon_t
\]
où $\{\mathbf{A}_i\}$ est une suite de matrice $m\times m$ et
$\varepsilon_t$ est un bruit blanc gaussien, de dimension $1\times
m$ d'espérance nulle et de variance
$\mathbb{V}\left[\varepsilon_t\right]=\Sigma$. Nous pourrions
compléter le modèle avec des variables exogènes, une constante par
exemple, mais nous allons à l'essentiel en omettant cette
possibilité.}\newline

\par{On note $\sample \equiv \{y^*_t\}_{t=-p+1}^T$ les données à notre
disposition et on note $z_t$ la concaténation horizontale des
vecteurs lignes $y^*_{t-1}$, $y^*_{t-2}$, ..., $y^*_{t-p}$. En concaténant
verticalement les vecteurs lignes $y^*_t$, $z_t$ et $\varepsilon_t$,
pour $t=1,\dots,T$, on obtient la représentation matricielle
suivante du modèle VAR($p$)~:
\[
Y = Z \mathcal{A} + E
\]
où $Y$ et $E$ sont des matrices $T\times m$, $Z$ est une matrice
$T\times(mp)$ et $\mathcal A = \left(\mathbf{A}_1',\mathbf{A}_2',\dots,\mathbf{A}_p'\right)'$
%\[
%\mathcal{A} = \left(%
%\begin{array}{c}
%  \mathbf{A}_1 \\ \mathbf{A}_2 \\ \vdots \\ \mathbf{A}_p \\
%\end{array}%
%\right)
%\]
la matrice $k\times m$ (avec $k = mp$) regroupant les coefficients
auto-régressifs. La vraisemblance associée à ce modèle
linéaire gaussien est donnée par~:
\[
\begin{split}
\mathcal L(\mathcal{A},\Sigma;\sample) =
&(2\pi)^{-\frac{mT}{2}}\left|\Sigma \right|^{-\frac{T}{2}}\\
&\times
e^{-\frac{1}{2}\trace\{(Y-Z\mathcal{A})\Sigma^{-1}(Y-Z\mathcal{A})'\}}
\end{split}
\]
L'estimateur du maximum de vraisemblance (MCO) est défini par~:
\[
\widehat{\mathcal{A}} = (Z'Z)^{-1}Z'Y
\]
et
\[
\widehat{\Sigma} =
T^{-1}(Y-Z\widehat{\mathcal{A}})'(Y-Z\widehat{\mathcal{A}})
\]
Nous verrons plus loin qu'il est profitable de réécrire la
vraisemblance en faisant apparaître l'estimateur des MCO~:
\[
\begin{split}
\mathcal L(\mathcal{A},\Sigma;\sample) &=
(2\pi)^{-\frac{mT}{2}}\\\times &\left|\Sigma \right|^{-\frac{k}{2}}
e^{-\frac{1}{2}\trace
\left\{\Sigma^{-1}(\mathcal{A}-\widehat{\mathcal{A}})'Z'Z(\mathcal{A}-\widehat{\mathcal{A}})\right\}}\\
\times &\left|\Sigma \right|^{-\frac{T-k}{2}} e^{-\frac{1}{2}\trace
\left\{\Sigma^{-1}(Y-Z\widehat{\mathcal{A}})'(Y-Z\widehat{\mathcal{A}})\right\}}.\\
\end{split}
\]
Aux constantes d'intégration près on reconnaît ici les fonctions de
densité de probabilité d'une gaussienne matricielle et d'une inverse
Wishart (voir l'annexe \ref{annex:A}). La vraisemblance se réécrit
donc sous la forme suivante~:
\[
\begin{split}
\mathcal L(\mathcal{A},\Sigma;\sample) =
&(2\pi)^{-\frac{mT}{2}}\times
(2\pi)^{\frac{km}{2}}\left|Z'Z\right|^{-\frac{m}{2}}\\
&\times f_{MN_{k,m}}(\mathcal{A};\widehat{\mathcal{A}},(Z'Z)^{-1},\Sigma)\\
&\times \frac{2^{\frac{\nu m}{2}}\pi^{\frac{m(m-1)}{4}}\prod_{i=1}^m
\Gamma
\left(\frac{\nu+1-i}{2}\right)}{|\widehat{S}|^{\frac{\nu}{2}}}\\&\times
f_{i\mathcal{W}_{m}}(\Sigma;\widehat{S},\nu)
\end{split}
\]
avec $\nu = T-k-m-1$ les degrés de liberté et $\widehat{S} =
T\widehat{\Sigma}$. Cette écriture nous apprend que la vraisemblance
du VAR($p$) est proportionnelle au produit de la densité d'une
normale matricielle et d'une loi inverse Wishart~:
\begin{equation}\label{equ:var:likprop}
\begin{split}
\mathcal L(\mathcal{A},\Sigma;\sample) &\propto
f_{MN_{k,m}}(\mathcal{A};\widehat{\mathcal{A}},(Z'Z)^{-1},\Sigma)\\
\times &
f_{i\mathcal{W}_{m}}(\Sigma;\widehat{S},\nu)
\end{split}
\end{equation}
Cette propriété va nous aider à poser une forme de la densité
\textit{a priori} telle que nous puissions obtenir une expression analytique
de la densité postérieure.}

\subsubsection{A priori non informatif}\label{sec:1:4:noninf}
\par{ Dans cette section nous supposons
que nos croyances sont non informatives en adoptant un \textit{a priori} plat à la Jeffrey~:
\begin{equation}\label{equ:var:jeffrey}
    p_0\left(\mathcal A,\Sigma\right) = \left|\Sigma\right|^{-\frac{m+1}{2}}
\end{equation}
On note que dans le cas scalaire, $m=1$, on retrouve le prior suggéré par Jeffrey ($1/\sigma^2$) décrit plus haut.
La densité \textit{a posteriori} satisfait donc~:
\[
p\left(\mathcal{A},\Sigma|\sample\right) \propto \left|\Sigma\right|^{-\frac{m+1}{2}}
\times\mathcal L(\mathcal{A},\Sigma;\sample)
\]
La densité jointe postérieure est donc proportionnelle au produit
d'une loi normale multivariée et d'une  loi inverse Wishart~:
\begin{equation}\label{equ:var:jeffrey:post1}
\begin{split}
 p(\mathcal{A},\Sigma;\sample) &\propto
f_{MN_{k,m}}(\mathcal{A};\widehat{\mathcal{A}},(Z'Z)^{-1},\Sigma)\\
\times &
f_{i\mathcal{W}_{m}}(\Sigma;\widehat{S},\tilde{\nu})
\end{split}
\end{equation}
avec $\tilde{\nu} = T-k$. Ainsi, la densité postérieure s'écrit sous
la forme suivante~:
\begin{equation}\label{equ:var:jeffrey:post2}
\begin{split}
 \mathcal{A}|\Sigma,\sample &\sim MN_{k,m}\left(\widehat{\mathcal{A}},\Sigma,(Z'Z)^{-1}\right)\\
 \Sigma|\sample &\sim i\mathcal{W}_m\left(\widehat{S},\tilde{\nu}\right)
\end{split}
\end{equation}
Il n'est pas surprenant de constater que la distribution postérieure de $\mathcal A$ (conditionnelle
à la matrice de variance covariance) est centrée sur l'estimateur du maximum de vraisemblance, puisque notre
\textit{a priori} est non informatif. Nous pourrions montrer, en intégrant par rapport à $\Sigma$, que la
distribution postérieure (marginale) de $\mathcal A$ est une version matricielle de la loi de Student
(voir \cite[chapitre 8]{ZellnerBook}). L'\textit{a priori} de Jeffrey n'affecte que le nombre de degré de
liberté de la distribution postérieure de $\mathcal A$. On obtient la densité marginale postérieure de $\sample$
en intégrant le noyau postérieur successivement par rapport à $\Sigma$ et $\mathcal{A}$~:
\begin{equation}\label{equ:var:jeffrey:marginaldensity}
\begin{split}
p\left(\sample\right) =
&(2\pi)^{-\frac{mT}{2}}\times
(2\pi)^{\frac{km}{2}}\left|Z'Z\right|^{-\frac{m}{2}}|\widehat{S}|^{-\frac{\tilde{\nu}}{2}}\\
&\times 2^{\frac{\tilde{\nu} m}{2}}\pi^{\frac{m(m-1)}{4}}\prod_{i=1}^m
\Gamma
\left(\frac{\tilde{\nu}+1-i}{2}\right).\\
\end{split}
\end{equation}
Cette quantité nous renseigne sur la qualité d'ajustement du modèle
VAR($p$). On note que la densité marginale de $\sample$ est une
fonction décroissante de la taille des erreurs ($|\hat{S}|$). Dans
cet exemple, nous pouvons caractériser la distribution postérieure
analytiquement. Notons néanmoins que même si nous connaissons
l'expression analytique de la distribution de $\mathcal A$ et
$\Sigma$, la construction des densités prédictives nécessite une
approche par simulations\footnote{Se reporter, par exemple, à
\cite{KadiyalaKarlsson1997} qui comparent différentes spécifications
des croyances \textit{a priori} et étudient les conséquences sur les
prévisions.}, puisque les prévisions sont des fonctions non
linéaires des matrices auto-régressives (dont nous connaissons la
distribution postérieure). L'intérêt pratique de l'approche
bayésienne peut paraître faible dans ce cas, dans la mesure où
l'espérance postérieure n'est pas différente de l'estimateur du
maximum de vraisemblance.}

\subsubsection{Un exemple d'a priori informatif}\label{sec:1:4:inf}

\par{Nous considérons maintenant un prior plus informatif qui va écarter l'espérance
de la distribution \textit{a posteriori} de l'estimateur du maximum de vraisemblance~;
dans un modèle linéaire gaussien, l'espérance \textit{a posteriori} est une combinaison convexe de l'estimateur
du maximum de vraisemblance et de l'espérance \textit{a priori}. Afin d'aller à l'essentiel\footnote{Nous
pourrions choisir un prior conjugué, c'est-à-dire une densité \textit{a priori} qui confrontée aux données
via la vraisemblance induit une densité postérieure de la même forme. Les propriétés des densités gaussienne
et Wishart, ainsi que l'équation
\equaref{equ:var:likprop}, suggèrent la densité jointe \textit{a priori} conjuguée suivante~:
\[
\left\{%
\begin{array}{ll}
    \mathcal{A}|\Sigma &\sim MN_{k,m}\left(\mathcal{A}_0,\Sigma,M_0^{-1}\right)\\
    \Sigma &\sim i\mathcal{W}_m\left(S_0,\nu_0\right)\\
\end{array}
\right.
\]
où $A_0$ est une matrice réelle de même dimension que $\mathcal A$, $\Sigma$ et $M_0$ sont des matrices symétriques
définies postives respectivement de dimensions $m \times m$ et $p\times p$, $S_0$ est une matrice symétrique définie
positive. On montre alors facilement que la densité postérieure est encore Normale-Wishart. Ce résultat est direct si
on couple le prior non informatif de la section \ref{sec:1:4:noninf} et un pré-échantillon pour former le prior normal
Wishart, voir \cite{TiaoZellner1964a} et la section \ref{sec:bvardsge:theil}.}, nous adoptons une
densité \textit{a priori} dégénérée
pour la matrice de variance-covariance des erreurs, en supposant que la matrice $\Sigma$ est connue
(on posera $\Sigma = \widehat{\Sigma}$). Enfin nous spécifions le prior sur $\mathcal A$ de la façon suivante~:
\begin{equation}\label{equ:var:gaussian:prior}
p_0(\VEC\ \mathcal A) \sim \normal{a_0}{\Omega_0}
\end{equation}
où $\Omega_0$ est une matrice symétrique définie positive de dimension $mp \times mp$. En multipliant la
vraisemblance par (\ref{equ:var:gaussian:prior}), on établit facilement que le noyau postérieur est~:
\begin{subequations}\label{equ:var:gaussian:posterior}
\begin{equation}\label{equ:var:gaussian:posterior:kernel}
\begin{split}
\mathcal K (\mathcal A|\sample) &= \exp\left\{ -\frac{1}{2}(a-a_1)'\Omega_1^{-1}(a-a_1)\right\}\\
&\times \exp\left\{ -\frac{1}{2}\left[a_0'\Omega_0^{-1}a_0 + \widehat{a}' \left(\Sigma^{-1}\otimes Z'Z\right)  \widehat{a}
-a_1'\Omega_1^{-1}a_1 \right]\right\}\\
&\times(2\pi)^{-\frac{km}{2}}|\Omega_0|^{-\frac{1}{2}}(2\pi)^{-\frac{mT}{2}}|\Sigma|^{-\frac{T}{2}}e^{-\frac{1}{2}\trace \Sigma^{-1}\widehat{S}}
\end{split}
%\begin{split}
%\mathcal K\left(\mathcal A\right|\sample ) = (2\pi)^{-\frac{T+k}{2}m}\left|\Omega_0\right|^{-\frac{1}{2}}
%\left|\Sigma\right|^{-\frac{T}{2}}&\times e^{-\frac{1}{2}\trace\ \left\{\Sigma^{-1}\widehat{S}\right\}}\\
%&\times e^{-\frac{1}{2}(\VEC\ \mathcal A - a_1)'\Omega_1^{-1}(\VEC\ \mathcal A - a_1)}
%\end{split}
\end{equation}
\begin{equation}\label{equ:var:gaussian:posterior:variance}
\Omega_1 = \left(\Omega_0^{-1} + \Sigma^{-1} \otimes Z'Z\right)^{-1}
\end{equation}
\begin{equation}\label{equ:var:gaussian:posterior:expectation}
a_1 = \Omega_1 \left[\Omega_0^{-1}a_0 + \left(\Sigma^{-1} \otimes Z'Z\right) \VEC \widehat{\mathcal A}\right]
\end{equation}
\end{subequations}
La distribution postérieure de $\mathcal A$ est donc gaussienne
$\normal{a_1}{\Omega_1}$, son interprétation est immédiate.
L'inverse de la variance postérieure ($\Omega_1^{-1}$, que l'on peut
interpréter comme une quantification de l'information \textit{a
posteriori}) est égale à la somme de l'inverse de la variance
\textit{a priori} ($\Omega_0^{-1}$, l'information \textit{a priori})
et de l'inverse de la variance de l'estimateur du maximum de
vraisemblance de $\mathcal A$ ($\Sigma^{-1} \otimes Z'Z$,
l'information apportée par les données). \textit{Ceteris paribus},
quand l'information \textit{a priori} est importante, la matrice de
variance-covariance $\Omega_0$ est petite, la variance \textit{a
posteriori} est faible. L'espérance postérieure est une combinaison
linéaire convexe de l'espérance \textit{a priori}, $a_0$, et de
l'estimateur du maximum de vraisemblance, $\VEC\ \widehat{\mathcal
A}$. Les pondérations respectives sont définies par le contenu
informatif des croyances \textit{a priori} et de l'échantillon.
Lorsque l'information \textit{a priori} tend vers l'infini,
\textit{ie} $\Omega_0 \rightarrow 0$, l'espérance postérieure tend
vers l'espérance \textit{a priori}. Lorsque l'information amenée par
les données tend vers l'infini, \textit{ie} $\Sigma^{-1}\otimes Z'Z
\rightarrow 0$, l'espérance \textit{a posteriori} tend vers
l'estimateur du maximum de vraisemblance. On peut donc interpréter
le paradigme bayésien comme un pont entre la calibration et
l'estimation par le maximum de vraisemblance. En notant que $Z'Z$
est généralement, si le modèle est stationnaire\footnote{La présence
d'une racine unitaire ne ferait qu'accroître l'ordre de divergence,
ce qui ne change pas qualitativement la conclusion.}, d'ordre
$\mathcal O(T)$, on retrouve un résultat conforme aux considérations
asymptotiques introduites dans la section \ref{sec:1:3}~:
l'espérance postérieure tend vers l'estimateur du maximum de
vraisemblance lorsque $T$ tend vers l'infini.}\newline

\par{\`A nouveau, en intégrant le noyau postérieur ~(\ref{equ:var:gaussian:posterior}) par rapport aux paramètres
auto-régressifs, on obtient une expression analytique de la densité marginale postérieure~:
\begin{equation}\label{equ:var:gaussian:posterior:marginaldensity}
\begin{split}
p(\sample) &= (2\pi)^{\frac{km}{2}}|\Omega_1|^{\frac{1}{2}}\\
&\times \exp\left\{ -\frac{1}{2}\left[a_0'\Omega_0^{-1}a_0 + \widehat{a}' \left(\Sigma^{-1}\otimes Z'Z\right)  \widehat{a}
-a_1'\Omega_1^{-1}a_1 \right]\right\}\\
&\times(2\pi)^{-\frac{km}{2}}|\Omega_0|^{-\frac{1}{2}}(2\pi)^{-\frac{mT}{2}}|\Sigma|^{-\frac{T}{2}}e^{-\frac{1}{2}\trace \Sigma^{-1}\widehat{S}}
\end{split}
\end{equation}
qui mesure la qualité d'ajustement du modèle et permet de comparer le VAR à d'autres modèles estimés à l'aide du
même échantillon.}

\subsubsection{La pratique}\label{sec:1:4:pratique}

\par{L'intérêt pratique de l'approche bayésienne des modèles VAR s'explique par l'équation
(\ref{equ:var:gaussian:posterior:variance}). Celle-ci établit que la variance postérieure de $\mathcal A$ est
inférieure à la variance de l'estimateur du maximum de vraisemblance, $\widehat{\mathcal A}$, dès lors que l'on
apporte de l'information \textit{a priori}. L'estimation des modèles VAR sur des données macroéconomiques pose
souvent des problèmes de précision. En effet, un modèle avec cinq variables et quatre retards demande l'estimation
de vingt paramètres alors que les échantillons sont habituellement de l'ordre de la centaine d'observations.
En incorporant de l'information à l'aide d'une densité \textit{a priori}  tout se passe comme si nous augmentions
le nombre de degrés de liberté. Ce gain en variance sur les paramètres du modèle, permettra d'obtenir des
prévisions ou des fonctions de réponses plus précises.}\newline

\par{On peut faire l'analogie avec l'incorporation de contraintes sur les paramètres
d'un modèle estimé dans le paradigme classique. Par exemple, si nous pensons que $\mathcal A$ doit satisfaire
les contraintes linéaires définies par $R\times\VEC\mathcal A = b$ (où $R$ est une matrice $
r\times mp$, $b$ est un vecteur $r\times 1$ et $r$ le nombre de restrictions linéaires), l'incorporation de
ces contraintes lors de l'estimation, \textit{ie} l'utilisation de moindres carrés contraints, permet de réduire
la variance des estimateurs et aussi l'erreur quadratique moyenne (dans la mesure où la contrainte n'est pas
en contradiction avec le processus générateur des données)
\footnote{Voir \citet[pages 52-54]{JudgeGriffithsHillLutkepohlLeeBook} la
section 3.2.1 intitulée << Exact Nonsample Information >>.}. Le paradigme bayésien est plus souple, dans le
sens où il ne pose pas des contraintes déterministes. Dans certains cas\footnote{Se reporter à
\citet[pages 670-673]{TheilBook}}, nous verrons une application dans la
section \ref{sec:bvardsge:theil}, la définition de croyances \textit{a priori} revient à poser une contrainte
probabiliste de la forme $R\times\VEC\mathcal A - b = \epsilon$, où $\epsilon$ est une variable aléatoire gaussienne.
Plus la variance de $\epsilon$ est importante, moins la contrainte sur $\VEC \mathcal A$ est forte
(plus l'information \textit{a priori} est floue).}\newline

\par{La formalisation de l'information \textit{a priori} ne se limite pas au choix de la forme d'une distribution.
Dans le cas du modèle BVAR de la section précédente, nous devons
aussi choisir les paramètres $a_0$ et $\Omega_0$. Dans cette
perspective, un prior qui s'est montré particulièrement efficace
quand on cherche à modéliser des séries macroéconomiques est le
prior de Minnesota\footnote{Voir, par exemple, \citet{Litterman1986}
ou \citet{KadiyalaKarlsson1997}. Pour d'autres priors on peut se
reporter à \citeauthor{KadiyalaKarlsson1997}.}. Celui-ci correspond
à la croyance \textit{a priori} que les séries observées sont des
marches aléatoires indépendantes. L'espérance \textit{a priori} de
$\VEC \mathcal A$ est alors telle que $\mathbb E \left[\mathbf
A_1\right] = I_m$ et $\mathbb E \left[\mathbf A_i\right] = 0_m$
pour $i=2,\dots,p$. La variance \textit{a priori} de $\VEC \mathcal
A$ est supposée diagonale. En notant $\omega_{i,j,k}$
($i,j=1,\dots,m$, $k=1,\dots,p$) la variance associée au paramètre
correspondant à la variable $j$ dans l'équation $i$ au retard $k$,
la variance \textit{a priori} est définie par~:
\[
\begin{split}
  \omega_{i,i,k} &= \frac{\pi_1}{k^{\pi_3}}\quad i=1,\dots,m\text{ et }k=1,\dots,p \\
  \omega_{i,j,k} &= \frac{\pi_2}{k^{\pi_3}}\frac{\sigma_i}{\sigma_j}\quad i=1,\dots,m\quad j\neq i\text{ et }k=1,\dots,p
\end{split}
\]
où les hyperparamètres $\pi_h$ pour $h=1,2,3$ sont positifs,
$\{\sigma_i^2\}$ est l'estimateur de la variance des résidus dans
l'estimation d'un AR(p) pour la variable $i$. Le ratio des
écarts-types permet de prendre en compte les différences
d'échelles entre les différentes variables composant le vecteur
des observables. La variance \textit{a priori} décroît lorsque le
retard $k$ augmente, ce qui traduit l'idée que plus le retard est
important plus nous croyons que la matrice $\mathbf A_k$ est
nulle. L'hyperparamètre $\pi_3$ indique à quelle vitesse la
variance \textit{a priori} tend vers zéro. Des valeurs fréquemment
utilisées pour $\pi_1$ et $\pi_2$ sont respectivement 0,05 et
0,005. Cela revient à dire, dans la mesure où $\sigma_i$ et
$\sigma_j$ sont proches, que nos croyances \textit{a priori} sont
plus fortes sur la nullité des termes hors des diagnonales de
$\mathbf A_k$ ($k=1,\dots,p$), c'est-à-dire sur l'absence de
causalité\footnote{Nous devons aussi poser un \textit{a priori}
sur la matrice de variance-covariance de l'innovation du VAR,
$\Sigma$. \cite{Litterman1986} considère que celle-ci est
diagonale et donnée (variance \textit{a priori} nulle). Cela ne
correspond pas à notre hypothèse de la section \ref{sec:1:4:inf},
où nous avons supposé que cette matrice était pleine (égale à
l'estimateur du maximum de vraisemblance). Nous pourrions, sans
changer qualitativement les résultats, adopter une matrice
diagonale dans la section \ref{sec:1:4:inf}. Nous choisissons de
poursuivre avec une matrice pleine, comme \cite{Phillips1996},
mais notre prior ne correspondra plus à des marches aléatoires
indépendantes. En fait la motivation principale de
\citeauthor{Litterman1986} était de justifier une estimation
équation par équation, car à l'époque l'estimation d'un système
était trop coûteuse numériquement, et ne reposait pas sur la
croyance que les séries macroéconomiques sont réellement
indépendantes.}. Notons que l'\textit{a priori} de Minnesota
suppose l'absence de relations de cointégration entre les
variables, il y a, \textit{a priori}, autant de racines unitaires
que de variables. Néanmoins rien n'empêche l'apparition de
relations de long terme dans la distribution postérieure.}\newline

\par{Il nous reste à choisir les valeurs des hyperparamètres du prior de Minnesota. Cette étape est importante
car l'expérience montre que l'inférence postérieure, en particulier l'évaluation de la densité marginale qui nous
permet d'évaluer le modèle, est très sensible à ce choix. Si, comme souvent dans la littérature (voir par exemple
\cite{SmetsWouters2004} ou \cite{VillaverdeRubio}), l'estimation d'un BVAR ne sert qu'à titre de comparaison
afin d'évaluer la qualité d'ajustement d'un modèle DSGE, le choix des hyperparamètres devient crucial. Ce point n'est
malheureusement jamais abordé dans la littérature. Le con\-tenu économique d'un modèle VAR étant faible, il paraît
difficile de recourir à la théorie pour spécifier la densité \textit{a priori}. Un critère objectif à notre
disposition est de choisir les hyperparamètres ($\pi_1,\pi_2$ et $\pi_3$) qui maximisent les performances en
prévisions du modèle BVAR. En spécifiant ainsi les priors de notre BVAR, nous savons au moins que nous ne comparons
par notre DSGE avec un BVAR aux performances prédictives médiocres. Dans cet esprit \citet{Phillips1996} propose
le critère PIC (\textit{Posterior Information Criterion}) que l'on peut minimiser par rapport aux
hyperparamètres. Ce critère peut être vu comme une généralisation, au cas non stationnaire, du bien connu
critère BIC\footnote{Voir \cite{Kim1998}.}. Dans le cas du modèle considéré dans la section \ref{sec:1:4:inf}
on choisit les hyperparamètres de la façon suivante~:
\begin{equation}
(\pi^*_1,\pi^*_2,\pi^*_3 ) = \arg\underset{\pi_1,\pi_2,\pi_3}{\min} \quad \log |\tilde{\Sigma}|
+\frac{1}{T}\log\frac{|\Omega_0^{-1}+\tilde{\Sigma}^{-1}\otimes Z'Z|}
{|\Omega_0^{-1}+\tilde{\Sigma}^{-1}_{T_0}\otimes Z_{T_0}'Z_{T_0}^{ }|}
\end{equation}
où $\tilde{\Sigma}$ est la matrice de variance covariance des
innovations au mode postérieur, les matrices indicées par $T_0$ sont
obtenues à partir du sous échantillon $1,\dots,T_0$ (où $T_0$ est
supérieur au nombre de paramètres estimés). Ici nous avons considéré
le nombre de retards comme une donnée, mais nous pourrions aussi
optimiser par rapport à $p$ le critère PIC (voir
\cite{Phillips1996}). À notre connaissance, l'utilisation
d'\textit{a priori} objectif pour les modèles BVAR, tel que
l'optimisation du critère PIC proposé par \citet{Phillips1996},
demeure inappliqué dans la littérature. On peut donc légitimement
douter de la pertinence des comparaisons entre BVAR et DSGE
effectuées jusqu'à présent\footnote{Quand
\citeauthor*{SmetsWouters2004} établissent que leur modèle DSGE
surpasse \textit{un} modèle BVAR en terme de densité marginale de
l'échantillon, on ne peut véritablement conclure à la bonne la
qualité du DSGE puisque nous n'avons aucune idée des performances du
modèle BVAR.}.}



\section{Modèles DSGE}\label{sec:dsge}

\par{Dans cette section nous présentons de façon générale les
modèles DSGE, puis soulignons les problèmes que peut poser
l'estimation. En particulier, nous expliquons pourquoi, à la
différence des BVAR, il n'est pas possible d'obtenir une expression
analytique de la distribution postérieure. Nous terminons en
présentant les méthodes de Monte Carlo utilisées pour caractériser
les croyances postérieures.}

\subsection{Résolution et vraisemblance}\label{sec:dsge:1}

\par{Nous limitons notre attention aux modèles \textsc{dsge} que nous
pouvons écrire sous la forme suivante~:
\begin{equation}\label{equ:dsge:def}
\mathbb E_t \left[\mathcal
F_{\theta}(y_{t+1},y_{t},y_{t-1},\varepsilon_t)\right] = 0
\end{equation}
avec $\varepsilon_t \sim \iid{0}{\Sigma}$, une variable aléatoire
dans $ \mathbb R^r $, les innovations structurelles, $ y_t \in
\Lambda\subseteq\mathbb R^n $ un vecteur regroupant les variables endogènes,
$\mathcal F: \Lambda^{3} \times \mathbb R^r \rightarrow \Lambda$ une
fonction réelle dans $\mathcal C^2$ paramétrée par un vecteur réel
$\theta \in \Theta \subseteq \mathbb R^q$ regroupant l'ensemble des
paramètres structurels du modèle. La fonction $\mathcal F$ est
simplement l'ensemble des équations qui définissent un modèle~; on a
autant d'équations que de variables endogènes. Le vecteur des
variables endogènes, $y_t$, inclut des variables d'état (endogènes
ou exogènes), des variables de choix et des variables statiques
\footnote{La première catégorie correspond aux variables prédéterminées, les suivantes
aux variables non prédéterminées.}. On supposera qu'il est possible d'exhiber une unique
solution stable et invariante au modèle décrit par l'équation \equaref{equ:dsge:def}~:
\begin{equation}\label{equ:dsge:sol}
y_t = \mathcal H_{\theta} \left(y_{t-1},\varepsilon_t\right)
\end{equation}
qui exprime les variables endogènes en fonction du passé et des chocs structurels
contemporains. La fonction paramétrée $\mathcal H_{\theta}$ regroupe les \textit{policy rules} et
les équations de transition (voir la contribution de Michel Juillard et Tarek Ocaktan dans ce numéro). La
solution \equaref{equ:dsge:sol}, en décrivant une récurrence stochastique non
linéaire, définit la distribution jointe d'un ensemble de variables.}\newline

\par{Pour estimer les paramètres $\theta$ du modèle, ou d'un sous ensemble des paramètres,
nous devons évaluer la vraisemblance associée au modèle
\equaref{equ:dsge:def} ou à sa forme réduite \equaref{equ:dsge:sol}. Même si la fonction
$\mathcal H_{\theta}$ est linéaire en $y_{t-1}$ et $\varepsilon_t$, cette évaluation ne peut être
 directe comme dans le cas du modèle VAR examiné plus haut. En effet, l'équation \equaref{equ:dsge:sol}
décrit la distribution jointe d'un ensemble de variables qui ne sont pas toutes observées.
Afin d'amener le modèle aux données on peut l'écrire sous une forme état-mesure~:
\begin{subequations}\label{equ:dsge:sol:ss}
\begin{equation}\label{equ:dsge:sol:measure}
        y^*_t = Z  y_t + \eta_t
\end{equation}
\begin{equation}\label{equ:dsge:sol:state}
        y_t = \mathcal H_{\theta} \left(y_{t-1},\varepsilon_t\right)
\end{equation}
\end{subequations}
où $y^*_t$ est un vecteur $m \times 1$, avec $r \leq m < n$, regroupant les variables observées et
$Z$ est une matrice de sélection $m \times n$. On peut éventuellement augmenter l'équation
de mesure d'un bruit blanc multivarié, $\eta_t$, représentant l'inadéquation des variables
théoriques avec les variables observées, ou plus simplement une erreur de mesure. On note $\mathcal
Y_T^* = \left\{y_t^*\right\}_{t=1}^T$ l'échantillon à notre disposition et
$\psi\in\Psi\subseteq \mathbb R^{q+\frac{n(n+1)}{2}+\frac{r(r+1)}{2}}$ le vecteur des paramètres
du modèle état-mesure ($\theta$, $\Sigma$ et éventuellement la matrice de variance-covariance de
$\eta_t$). La vraisemblance est la densité de l'échantillon, conditionnellement aux paramètres $\psi$ et
au modèle défini par \equaref{equ:dsge:sol:ss}~:
\begin{equation}\label{equ:dsge:likelihood}
    \mathcal L(\psi ; \sample) = p\left(\sample | \psi\right) = \prod_{t=1}^T
    p\left(y_t^*|\samplet{t-1},\psi \right)
\end{equation}
L'évaluation de la densité de $y_t^*$ conditionnellement à $\samplet{t-1}$ n'est généralement pas
directe, %Nous pourrions imaginer un cas où $y_t^*$ regroupe l'ensemble des variables
%d'état du modèle \equaref{equ:dsge:def}, dans ce cas il ne serait pas nécessaire d'avoir des idées
%sur la distribution des variables latentes.
dans la mesure où $y_t^*$ dépend de variables endogènes inobservables. Nous pouvons néanmoins
utiliser la relation suivante~:
\begin{equation}\label{equ:prediction}
    p\left(y_t^*|\samplet{t-1},\psi \right) = \int_{\Lambda} p\left(y_t^* | y_t , \psi \right)p(y_t
    | \samplet{t-1},\psi) \mathrm d y_t
\end{equation}
La densité de $y_t^*$ conditionnellement à $\samplet{t-1}$ est obtenue comme la moyenne de la densité
de $y^*_t$ sachant $y_t$, pondérée par la densité de $y_t$ sachant $\samplet{t-1}$. La première
densité sous l'intégrale est specifiée par l'équation de mesure (\ref{equ:dsge:sol:measure}).
L'évaluation de la densité de la prévision des variables latentes, conditionnellement l'information disponible
en $t-1$, est moins directe, et on doit utiliser un filtre de Kalman. Il s'agit d'une procédure
récursive. À chaque date (entre 1 et T) on forme une prévision des variables latentes ($y_t$
sachant $\samplet{t-1}$), en utilisant l'équation d'état ($y_{t}$ sachant $y_{t-1}$,
équation \ref{equ:dsge:sol:state}) et une estimation initiale des variables latentes ($y_{t-1}$
sachant $\samplet{t-1}$), puis on corrige cette prévision quand une nouvelle observation ($y^*_t$)
augmente l'ensemble d'information. On peut interpréter cette démarche comme une estimation récursive
bayésienne du vecteur des variables latentes. En initialisant les variables latentes avec la densité
associée à la distribution ergodique des variables latentes définie par l'équation d'état
(\ref{equ:dsge:sol:state}),
$p\left(y_0|\samplet{0},\psi\right)=p\left(y_0|\psi\right)=p\left(y_{\infty}|\psi\right)$, la
récursion s'écrit de la façon suivante~:
\begin{subequations}\label{equ:kalman}
\begin{equation}\label{equ:kalman:prediction}
        p\left(y_t|\samplet{t-1},\psi\right) = \int_{\Lambda}
        p\left(y_t|y_{t-1},\psi\right)p\left(y_{t-1}|\samplet{t-1},\psi\right)\mathrm d y_{t-1}
\end{equation}
\begin{equation}\label{equ:kalman:update}
        p\left(y_t|\samplet{t},\psi\right) =
        \frac{p\left(y^*_t|y_t,\psi\right)p\left(y_t|\samplet{t-1},\psi\right)}
        {\int_{\Lambda} p\left(y^*_t|y_t,\psi\right)p\left(y_t|\samplet{t-1},\psi\right)\mathrm d
y_t}
\end{equation}
\end{subequations}
L'interprétation de l'équation de prédiction (\ref{equ:kalman:prediction}) est directe~: la densité
de la prédiction des variables latentes en $t$ est la moyenne de la densité de
$y_t$ sachant $y_{t-1}$, définie par l'équation d'état (\ref{equ:dsge:sol:state}), pondérée par la
densité de $y_{t-1}$ sachant $\samplet{t-1}$. Cette dernière densité est définie par l'équation de
mise à jour (\ref{equ:kalman:update}) ou la condition initiale. L'équation de mise à jour
(\ref{equ:kalman:update}) est, à l'instar de l'équation (\ref{equ:bayes:posterior}), une application
directe du théorème de Bayes. Le premier terme au numérateur, $p\left(y_t|\samplet{t-1}\right)$, est
la densité \textit{a priori} du vecteur des variables latentes. Le second terme,
$p\left(y^*_t|y_t\right)$, la densité de l'observation sachant l'état obtenu via l'équation de
mesure (\ref{equ:dsge:sol:measure}), est la vraisemblance. Le dénominateur est la densité marginale
de la nouvelle observation.}\newline

\par{Puisque nous pouvons, au moins théoriquement, évaluer la vraisemblance associée au modèle
DSGE, nous devrions être capable d'estimer ses paramètres.
Malheureusement, les équations (\ref{equ:prediction}) et
(\ref{equ:kalman}) nécessitent l'évaluation
d'intégrales\footnote{Plusieurs approches sont envisageables~:
quadrature, Monte Carlo, quasi Monte Carlo... voir \cite{JuddBook}.}
dans l'espace des variables d'état. Quand le nombre de variables
latentes agmente il devient très coûteux d'évaluer ces intégrales
(on parle de \textit{curse of dimensionality}). La dérivation de la
forme réduite du modèle (\ref{equ:dsge:sol}) nécessite également
l'évaluation d'intégrales. En pratique, même pour des modèles de
dimensions modestes, l'évaluation de la vraisemblance est difficile.
Nous devons donc approximer celle-ci. Dans le cas où le modèle
état-mesure (\ref{equ:dsge:sol:ss}) est linéaire et gaussien,
l'évaluation des intégrales devient très simple car les variables
latentes et observées sont normalement distribuées à chaque date.
Ainsi la dynamique de la distribution des variables latentes est
complétement caractérisée par la dynamique de l'espérance et de la
variance des variables latentes. On peut trouver une présentation du
filtre de Kalman dans ce cas simple dans \citet[chapitre
13]{GourierouxMonfortStatBook1}, \citet[chapitre 3]{HarveyBook} ou
encore dans la contribution de Fabrice Collard et Patrick Fève à ce
même numéro. Ceci explique pourquoi les modèles DSGE estimés sont
généralement (log-) linéarisés autour de l'état
stationnaire.}\newline

\par{Malgré l'approximation (log-) linéaire de la forme réduite du modèle ($\mathcal
H_{\theta}$) l'évaluation de la vraisemblance est numérique. Nous ne
disposons pas d'une expression analytique, comme dans le cas du
modèle VAR, et ne pouvons donc écrire formellement la densité
postérieure ou les moments postérieurs. Deux possibilités s'offrent
à nous.\newline La première est de considérer une approximation
asymptotique de la densité postérieure. Il est alors possible
d'approximer, voir la section \ref{sec:1:3} et
\citet{TierneyKadane1986}, tout moment \textit{a posteriori} ou la
densité marginale. Nous avons vu que l'erreur d'approximation des
moments est d'ordre $\mathcal O(T^{-2})$ et que l'erreur
d'approximation de la densité marginale est d'ordre $\mathcal
O(T^{-1})$. L'expérience suggère, pour les dimensions d'échantillon,
$T$, généralement considérées dans la littérature, que
l'approximation de Laplace de la densité marginale est
satisfaisante\footnote{Par exemple, \citet[tableau
2]{SmetsWouters2004} estiment la densité marginale de leur modèle
DSGE à l'aide de l'approximation de Laplace et d'une méthode <<
exacte >> (dans le sens où elle ne repose pas sur une
approximation asymptotique, voir plus bas) basées sur des
simulations. Avec l'approximation de Laplace ils obtiennent (en
logarithme) -269,59 à comparer au -269,20 obtenu avec un exercice de
Monte-Carlo. Ces deux évaluations sont très proches, on retrouve
généralement cette proximité dès lors que l'échantillon est assez
grand.}.\newline La deuxième possibilité est d'évaluer les moments
en recourant à des simulations par Monte-Carlo. L'intuition de cette
approche repose sur la loi des grands nombres. Supposons, par
exemple, que nous souhaitions évaluer l'espérance d'une variable
aléatoire $X$ de distribution $\mathcal  G$. Si l'on génère une
suite de variables aléatoires $X_1, X_2, \dots, X_n$  indépendantes
et distribuées selon $\mathcal{G}$, alors une approximation de
l'espérance est donnée par la moyenne empirique de ces
variables\footnote{En fait on peut montrer, dans certains cas, qu'il
s'agit d'un estimateur du maximum de vraisemblance.}~; la loi des
grands nombres assure que l'erreur d'approximation tend vers zéro
presque sûrement lorsque le nombre de tirages $n$ tend vers
l'infini. Si on admet de plus que le moment du second ordre existe,
la vitesse de convergence est $\mathcal{O}(n^{1/2})$ par application
du théorème \textit{central limit}. En pratique nous pouvons être
intéressés par les moments postérieurs de $\theta$. En notant que~:
\begin{equation}\label{equ:dsge:postmom}
\mathbb E \left[ \varphi(\theta) \right] = \int_{\Theta} \varphi (\theta) p_1(\theta|\sample)\mathrm d\theta
\end{equation}
il semble alors naturel d'utiliser la moyenne empirique de
$\left(\varphi(\theta^{(1)}),\varphi(\theta^{(2)}),\dots,\varphi(\theta^{(n)})\right)$,
où les $\theta^{(i)}$ sont des tirages indépendants dans la
distribution postérieure, pour évaluer l'espérance de $\varphi
(\theta)$. L'erreur d'approximation tend vers zéro lorsque le
nombre de simulations ($n$) tend vers l'infini. il convient de
noter que $p_1$ est généralement d'une forme inconnue et que l'on
ne peut donc pas définir un générateur pseudo-aléatoire
reproduisant la distribution \emph{a posteriori}. Comme nous le
verrons par la suite, la méthode de Monte Carlo dite de fonction
d'importance permet de remédier à cette difficulté sous certaines
conditions. Ce principe de Monte Carlo se généralise au cas où les
variables simulées ne sont pas indépendantes. Il est ainsi
possible dans certains cas de construire, moyennant certaines
conditions, une chaîne de Markov $\theta^{(1)},
\theta^{(2)},\dots, \theta^{(n)}$ de loi stationnaire (ergodique)
$p_1$ telle que la moyenne empirique des $\varphi(\theta^{(j)})$
($\varphi$ est la statistique d'intérêt) converge presque sûrement
vers la quantité d'intérêt comme dans le cas des tirages
indépendants. Tout l'art de l'exercice est alors de déterminer une
chaîne de Markov (et plus précisément son noyau de transition)
telle que sa loi ergodique corresponde à la loi a posteriori
désirée et d'évaluer le temps d'arrêt des simulations, \textit{ie}
de diagnostiquer la convergence de la chaîne de
Markov\footnote{Dans ce qui suit, nous omettons les méthodes
d'échantillonnage de Gibbs. Cette méthode consiste à générer
chaque paramètre conditionnellement à tous les autres paramètres.
Il est donc nécessaire de pouvoir écrire toutes les distributions
conditionnelles. C'est pourquoi cette méthode n'est généralement
pas priviligiée pour l'estimation des modèles DSGE. Cependant, il
est à noter que les algorithmes de Metropolis-Hasting et
d'échantillonnage de Gibbs peuvent être combinés, on parle alors
d'algorithme \emph{Metroplis-Within-Gibbs}.}.}


\subsection{L'échantillonnage bayésien par fonction d'importance}

\par{Idéalement, on souhaiterait générer les paramètres d'intérêt suivant la distribution \emph{a
posteriori}. Cependant, à l'exception de rares modèles, ceci n'est pas possible. On peut néanmoins exploiter le
fait qu'il n'est pas nécessaire de générer une suite de tirage $\{\theta^{(s)}\}$, suivant la distribution
\emph{a posteriori}, pour obtenir une évaluation correcte des moments postérieurs. En effet, sous certaines conditions
de régularité, on peut utiliser une densité de probabilité $h$, définie sur le même espace et appelée fonction
d'importance, suffisamment proche de $p_1$ (dans un sens à préciser) et échantillonner à partir de cette dernière. Il
est alors possible de montrer par la loi des grands nombres que l'intégrale (\ref{equ:dsge:postmom}) définissant le
moment postérieur est approchée par~:
\[
  \mathbb E\left[\varphi(\theta) \right] \approx \frac{\sum_{s=1}^n
  \varphi(\theta^{(s)})w(\theta^{(s)})}{\sum_{s=1}^n w(\theta^{(s)})}.
\]
avec
\begin{equation}\label{equ:dsge:weight}
w(\theta^{(s)}) = \frac{p(\mathcal{Y}^*_T \mid
\theta^{(s)})p_0(\theta^{(s)})}{h(\theta^{(s)})}
\end{equation}
Le choix de la fonction d'importance est crucial~: elle doit être suffisamment proche de la loi \emph{a posteriori},
ce qui n'est pas toujours simple en pratique. En effet, si $h$ est une mauvaise approximation de $p_1$, alors les
poids sont généralement faibles pour la plupart des valeurs échantillonnées de $\theta$, la somme est alors dominée
par quelques termes dont les poids sont très élevés. Il en résulte une estimation peu fiable, voir
\cite{CasellaRobertBook} pour plus de détails. L'algorithme se résume comme suit :
\begin{algorithm}

\item[(1)] Maximiser le noyau postérieur par rapport à $\theta$. On obtient le mode de la densité postérieure,
$\theta^m$, et le hessien au mode qui caractérise la courbure de la densité postérieure au mode et dont l'inverse de
l'opposé, noté $\Sigma(\theta^m)$, approxime la variance postérieure.

\item[(2)] Générer $\theta^{(s)}$, suivant une fonction d'importance,$h$, dont les moments du premier
et second ordre dépendent de $\theta^{m}$ et $ \Sigma(\theta^{m})$

\item[(3)] Déterminer les poids $w(\theta^{(s)})$ selon (\ref{equ:dsge:weight}).

\item[(4)] Reprendre (2-3) pour $s=1,\dots,n$.

\item[(5)] Calculer~:
\[
\frac{\sum_{s=1}^n \varphi(\theta^{(s)})w(\theta^{(s)})}{\sum_{s=1}^n w(\theta^{(s)})}
\]

\end{algorithm}
La première étape n'est pas spécifique à l'algorithme par fonction d'importance~: il s'agit
de calculer les moments associés à la fonction d'importance, $h$.  Cette <<  calibration >>
des moments de la fonction d'importance est généralement faite à partir de la maximisation du logarithme du noyau
\emph{a posteriori}. Étant données les propriétés asymptotiques de la distribution postérieure, ces choix sur les
moments d'ordre un et deux associés à $h$ sont d'autant plus satisfaisants que la taille de l'échantillon est
importante.}\newline

\par{Bien que très populaire en statistique, cette méthode est peu utilisée dans le cadre de l'estimation de modèles
DSGE. A titre d'exemples, \cite{DeJongIngramWhiteman00} estiment
avec cette méthode un modèle de croissance stochastique
linéarisée. An et Schorfheide (2007) comparent l'algorithme
d'importance avec celui de Metropolis (à pas aléatoire) dans une
version simplifiée du modèle de \citet{SmetsWouters2004}. Pour ce
faire, ils retiennent comme fonction d'importance une distribution
de Student multivariée.}

\subsection{Les méthodes de Monte-Carlo à chaînes de Markov}\label{sec:dsge:simul:MCMC}

\par{Cette seconde classe d'algorithmes permet de générer des variables aléatoires suivant approximativement
la loi \emph{a posteriori}, lorsque cette dernière n'est pas disponible. Elle évite donc l'appel à une fonction
d'importance, $h$, souvent difficile à déterminer pour les modèles DSGE\footnote{Nous verrons par la suite qu'un
algorithme à chaînes de Markov indépendantes peut s'interpéter comme un algorithme par fonction
d'importance. Par ailleurs, les méthodes de Monte-Carlo à chaînes de Markov nécessitent de déterminer le noyau de
transition de la chaîne de Markov, dont on sait seulement qu'il vérifie certaines propriétés d'ergodicité, etc.
Le choix de la fonction définissant les changements d'état peut être assimilée, \emph{toutes choses égales par ailleurs},
à celui de la fonction d'importance.}. On cherche ainsi à définir une chaîne de Markov dont la distribution
ergodique est approximativement le noyau \emph{a posteriori}. Si cette chaîne existe, la méthode d'échantillonnage est
grossièrement définie comme suit. Dans un premier temps, on initialise (arbitrairement) la chaîne de Markov. Dans
un second temps, on génère les $\theta^{(s)}$ à partir de cette chaîne. A l'issue d'un certain nombre de tirages
(disons $n_0$), on dispose de réalisations de variables aléatoires $\{\theta^{(s)}, s =n_0, \cdots, n\}$ approximativement
distribuées comme la distribution \emph{a posteriori}.}

\subsubsection{Chaînes de Markov}\label{sec:dsge:simul:MCMC:CM}
\par{Une chaîne de Markov est une suite de variables aléatoires continues à valeurs dans $\Theta$,
$\left(\theta^{(0)}, \dots,\theta^{(n)}\right)$, générée par un
processus de Markov. Une suite de variables aléatoires est générée
par un processus de Markov (d'ordre 1) si la distribution de
$\theta^{(s)}$ ne dépend que de $\theta^{(s-1)}$. Une chaîne de
Markov est caractérisée par un noyau de transition qui spécifie la
probabilité de passer de $\eta \in \Theta$ à $S \subseteq \Theta$.
Nous noterons $P(\eta,S)$ le noyau de transition, il vérifie
$P(\eta,\Theta) = 1$ pour tout $\eta$ dans $\Theta$. Si la chaîne de
Markov définie par le noyau $P$ converge vers une distribution
invariante $\pi$, alors le noyau doit satisfaire l'identité
suivante~:
\[
 \pi(S) = \int_{\Theta} P(\eta,S)\pi\left(\mathrm d\eta\right)
\]
pour tout sous ensemble mesurable $S$ de $\Theta$. Plus généralement, avant d'atteindre la distribution ergodique $\pi$,
si nous notons $P^{(s)}(\eta,S)$ la probabilité que $\theta^{(s)}$ soit dans $S$ sachant que $\theta^{(s-1)} = \eta$, nous
avons~:
\[
P^{(s)}(\eta,S) = \int_{\Theta} P(\nu,S)P^{(s-1)}\left(\eta,\mathrm d\nu\right)
\]
la distribution de $\theta$ s'ajuste d'itération en itération puis rejoint la distribution ergodique,
$\lim_{s\rightarrow\infty}P^{(s)}(\eta,S) = \pi(S)$. L'idée est alors de choisir le noyau de transition qui nous aménera
vers la distribution invariante désirée.}\newline

\par{Définissons $p(\eta,\nu)$ et $\tilde{\pi}$ les densités associées au noyau $P$ et à la distribution
$\pi$\footnote{\label{note:mc:dens}Le noyau $P(\eta,S)$ spécifie la probabilité d'aller de $\eta$ à $S$. Dans un
cas favorable,
$\theta$ est en $S$ à l'itération suivante, nous pouvons envisager deux possibilités~: (\textit{i}) $\theta$
se déplace effectivement et rejoint la région $S$ à l'itération suivante, (\textit{ii}) $\theta$ ne se déplace pas
mais $\eta$ est déjà dans $S$. La densité associée au noyau est donc une densité continue - discrète,
\citeauthor{Tierney1994} adopte la définition suivante~:
\[
P(\eta,\mathrm d \nu) = p(\eta,\nu)\mathrm d\nu + (1-r(\eta))\delta_{\eta}(\mathrm d \nu)
\]
où $p(\eta,\nu)\equiv p(\nu|\eta)$ est la densité associée à la
transition de $\eta$ à $\nu$,  $r(\eta) = \int p(\eta,\nu)d\nu < 1$,
$1-r(\eta)$ s'interprète comme la probabilité de ne pas quitter la
position $\theta = \eta$, $\delta_{\eta}(S)$ est une fonction
indicatrice égale à un si et seulement si $\eta \in S$.}.
\citet{Tierney1994} montre que si la densité $p(\eta,\mu)$ vérifie
la condition de réversibilité\footnote{Il ne s'agit pas à proprement
parler de la condition de réversibilité, mais d'une implication de
la propriété de réversibilité.}~:
\[
\tilde{\pi}(\eta)p(\eta,\nu) = \tilde{\pi}(\nu)p(\nu,\eta)
\]
alors $\pi$ est la distribution invariante associée au noyau $P$\footnote{(Suite de la note \ref{note:mc:dens})
Techniquement, il suffit de substituer la définition du noyau dans $\int_{\Theta} P(\eta,S)
\pi(\eta)\mathrm d \eta$ qui, si la chaîne est réversible, se réduit alors à $\pi(S)$.}. De façon équivalente~:
\[
\frac{\tilde{\pi}(\eta)}{\tilde{\pi}(\nu)} = \frac{p(\nu,\eta)}{p(\eta,\nu)} >1
\]
Cette condition nous dit simplement que si la densité de $\theta = \eta$,
$\tilde{\pi}(\eta)$, domine la densité associée à $\theta = \nu$, $\tilde{\pi}(\nu)$, alors il doit être plus
<<  facile >> de passer de $\nu$ à $\eta$ que de $\eta$ à $\nu$.}\newline

\par{Cette propriété nous aidera à construire une chaîne de Markov dont la distribution invariante est la distribution
postérieure des paramètres $\theta$ dans le modèle DSGE. On comprend bien que le noyau de cette chaîne est difficile
à définir. Supposons que l'on puisse choisir un noyau de transition $Q(\eta,S)$~; alors il est presque sûr que la
condition de réversibilité ne sera pas vérifiée, c'est-à-dire que nous aurons
$p_1\left(\eta|\sample\right)q(\eta,\nu) \neq p_1\left(\nu|\sample\right)q(\nu,\eta)$. L'algorithme de
Metropolis-Hastings est une approche générale qui permet de <<  corriger >> ce noyau, de façon à respecter
la condition de réversibilité.}


\subsubsection{L'algorithme de Métropolis-Hasting}\label{sec:dsge:simul:MCMC:MH}

\par{Supposons que l'on puisse définir une densité instrumentale, qui permette d'approcher le noyau
de transition de la chaîne de Markov dont la densité ergodique est la loi \emph{a posteriori} de notre modèle.
Cette densité est définie par $q(\eta,\nu) \equiv q(\nu \mid \eta)$.}

\begin{algorithm}[Metropolis-Hastings]

\item[(1)] Se donner une condition initiale $\theta^{(0)}$ telle que $\mathcal K\left(\theta^{(0)}|\sample\right)>0$ et
poser $s=1$.

\item[(2)] Générer un candidat (une proposition) $\theta^*$ à partir d'une densité $q(\theta^{(s-1)}, \theta^*)$.

\item[(3)] Générer $u$ dans une loi uniforme entre $[0,1]$

\item[(4)] Appliquer la règle suivante~:
\[
\theta^{(s)} =
\begin{cases}
\theta^*   &\text{ si } \alpha\left(\theta^{(s-1)},\theta^*\right) > u\\
\theta^{(s-1)}&\text{ sinon.}
\end{cases}
\]
où
\[
\alpha(\theta^{(s-1)},\theta^*) = \min \left\{1,
\frac{\mathcal K\left(\theta^* \mid \sample\right)}{\mathcal K\left(\theta^{(s-1)} \mid \sample\right)}
\frac{q\left(\theta^{(s-1)} \mid \theta^*\right)}{q\left(\theta^* \mid \theta^{(s-1)}\right)} \right\}
\]

\item[(5)] Reprendre (2-4) pour $s=2,\dots,n$.

\end{algorithm}
\par{\noindent Notons qu'il suffit de pouvoir évaluer le noyau postérieur pour mettre en oeuvre cet algorithme~; la
connaissance de la densité postérieure à une constante près est suffisante. L'algorithme de Metropolis-Hasting
requiert le choix d'une fonction instrumentale $q$ à partir de laquelle on génère des transitions dans l'espace des
paramètres. La densité conditionnelle $q$ permet de générer un vecteur candidat $\theta^{*}$. Puisqu'elle n'est
pas nécessairement la densité conditionnelle associée au noyau de transition dont la distribution ergodique est la
distribution \emph{a posteriori} recherchée, la condition de réversibilité n'est pas vérifiée\footnote{Le noyau
de transition du MH, $Q(\eta,S)$, est défini de la même façon que $P$ dans la section \ref{sec:dsge:simul:MCMC:CM}
et la note \ref{note:mc:dens}.}. L'algorithme de MH corrige cette erreur\footnote{Par exemple, nous pourrions avoir~:
\[
\mathcal K \left(\eta \mid \sample \right)
q(\eta,\nu) > \mathcal K\left(\nu \mid
\sample \right) q(\nu , \eta)
\]
Dans ce cas, l'échantillonage à partir de $q$ ne propose pas assez souvent des transitions de $\theta=\nu$ à
$\theta=\eta$ ou trop souvent des mouvements de $\theta=\eta$ à $\theta=\nu$. L'algorithme de MH corrige cette
erreur en n'acceptant pas systématiquement les propositions de $q$. En introduisant une probabilité d'acceptation
de la transition proposée, $\alpha$, élevée (faible) quand il s'agit de rejoindre une région où la densité
\emph{a posteriori} est élevée (faible), on peut rétablir l'équilibre et finalement vérifier la condition de
réversibilité. Dans notre exemple, la probabilité d'acceptation de la transition de $\nu$ à $\eta$ devrait être la
plus grande possible puisque $q$ ne propose pas assez souvent cette transition, nous poserons donc
$\alpha(\nu,\eta) = 1$. À l'inverse la densité conditionnelle $q$ propose trop de transitions de $\eta$ vers $\nu$,
la probabilité d'acceptation associée à cette proposition, $\alpha(\eta,\nu)$, doit donc être inférieure à 1. Pour
équilibrer les deux transitions, elle doit être telle que~:
\[
  \mathcal K \left(\eta | \sample \right)
  q(\eta,\nu)\alpha(\eta,\nu) = \mathcal K\left(\nu |
  \sample \right)q(\nu,\eta) \alpha(\nu , \eta)
\]
soit, puisque $\alpha(\nu,\eta) = 1$, de façon équivalente~:
\[
  \alpha(\eta,\nu) = \frac{\mathcal K\left(\nu |
  \sample \right)}{\mathcal K \left(\eta | \sample \right)} \frac{q(\nu,\eta)}{q(\eta,\nu)}
\]
On ne rejette donc pas systématiquement la transition proposée par $q$. En considérant l'exemple en renversant
l'inégalité, on comprend la règle donnée dans l'étape 4 de l'algorithme 2.} en n'acceptant pas systématiquement les
propositions de $q$. En introduisant une probabilité d'acceptation de la transition proposée, $\alpha$, on peut
finalement vérifier la condition de réversibilité. Pour cela, la probabilité d'acceptation doit être telle que~:
\[
  \mathcal K \left(\eta | \sample \right)
  q(\eta,\nu)\alpha(\eta,\nu) = \mathcal K\left(\nu |
  \sample \right)q(\nu,\eta) \alpha(\nu , \eta)
\]
soit
\[
\alpha(\eta,\nu) = \min \left\{1,\frac{\mathcal K\left(\nu |
  \sample \right)}{\mathcal K \left(\eta | \sample \right)} \frac{q(\nu,\eta)}{q(\eta,\nu)}\right\}
\]
Il nous reste à déterminer (\textit{i}) comment nous devons
initialiser la chaîne et (\textit{ii}) la longueur de la chaîne.
Nous reviendrons par la suite, lors de la présentation de
l'algorithme de Metropolis à pas aléatoires, sur le premier point.
Nous aborderons la question du nombre de simulations nécessaires,
c'est-à-dire de la longueur de la chaîne, dans la section 3.3.4.
Nous reviendrons plus loin sur cette question, pour l'instant nous
supposerons que pour tout $s>n_0$ les $\theta^{(s)}$ sont tirés dans
la distribution ciblée. Afin de s'assurer que les résultats sont
indépendants des conditions initiales, on ne considère pas les
simulations indicées par $s=0,\dots,n_0$. Ainsi, pour évaluer
$\mathbb E [\varphi(\theta)]$ nous calculons~:
$(n-n_0)^{-1}\sum_{s=n_0+1}^n\varphi(\theta^{(s)})$ qui converge
vers le moment postérieur recherché lorsque le nombre de
simulations, $n$, tend vers l'infini.}

\subsubsection{Deux variantes de l'algorithme MH}

\par{\noindent\textbf{L'algorithme de MH à pas aléatoires.} Comme nous l'avons expliqué plus haut, l'utilisation
de l'algorithme de Metropolis-Hastings repose sur le fait qu'il est aisé d'échantillonner à partir de la densité
instrumentale $q$. L'inconvénient est que cette dernière n'est pas toujours facile à déterminer.  Dans cette
perspective, l'algorithme de Metropolis à pas aléatoires est utile lorsqu'il est difficile d'obtenir une
bonne approximation de la densité \emph{a posteriori}. Une proposition à l'itération $s$ est définie par~:
\[
\theta^{*} = \theta^{(s-1)}+z
\]
où $z$ est le pas aléatoire. Le choix de la densité de $z$ détermine la forme précise de la densité
instrumentrale, $q$. Un choix standard est la distribution gaussienne multivariée~: $z\sim \normal{0}{\Sigma}$. Ainsi
la densité de $\theta^{*}$ conditionnel à $\theta^{(s-1)}$ est gaussienne~:
\[
q\left(\theta^{(s-1)},\theta^*\right) \equiv q\left(\theta^* \left| \theta^{(s-1)}\right.\right) \sim
\normal{\theta^{(s-1)}}{ \Sigma}
\]
Par symétrie de la loi normale, la densité instrumentale vérifie $q(\eta,\nu) = q(\nu,\eta)$, ainsi la probabilité
d'acceptation ne dépend que du noyau postérieur~:
\[
\alpha(\eta,\nu) = \min \left\{1,\frac{\mathcal K(\nu|\sample)}{\mathcal K (\eta|\sample)} \right\}
\]
Autrement dit, si $\mathcal K(\theta^* |\sample) \geq \mathcal K(\theta^{(s-1)} |\sample)$, la chaîne de Markov se
déplace en $\theta^*$. Si ce n'est pas le cas, la chaîne se déplace avec une probabilité égale au rapport des densités
\emph{a posteriori}. On accepte avec une probabilité unitaire la proposition dans une phase ascendante
(c'est-à-dire lorsque la probabilité \emph{a posteriori} croît) et avec une probabilité non nulle la proposition dans
une phase descendante (si nous décidions de rejeter systématiquement ces propositions défavorables la chaîne ne Markov
ne visiterait pas complétement l'espace des paramètres et \textit{a fortiori} elle ne pourrait converger vers la
distribution postérieure recherchée).}\newline


\par{La probabilité d'acceptation peut être ajustée à partir de la matrice de variance-covariance $\Sigma$. En effet, si
celle-ci est <<  grande >> alors il y a de fortes chances pour que la transition proposée nous amène dans les queues
de la distribution \textit{a posteriori} c'est à dire dans une région où la densité est faible et où donc la probabilité
d'acceptation est proche de zéro (si on vient d'une zone plus dense). Si la matrice $\Sigma$ est <<  petite >> alors
les transitions proposées ne seront que des petits pas, dans ce cas il n'y aura pas de grandes variations de la densité
postérieure et donc la probabilité d'acceptation sera proche de un. Nous écrivons $\Sigma$ sous la forme $c \Omega$ où
$\Omega$ est une estimation de la matrice de variance-covariance \emph{a posteriori} et $c$ est un paramètre
d'échelle qui permet de jouer sur le caractère plus ou moins diffus de la matrice de variance-covariance et donc
d'ajuster le taux d'acceptation. On peut expérimenter différentes valeurs de $c$ afin de déterminer une  probabilité
d'acceptation raisonnable\footnote{Il n'existe pas une règle universelle. Un taux d'acceptation de l'ordre de 0.25-0.40
est généralement considéré comme approprié. Dans le même temps, il est important de noter que ce n'est pas tant le taux
d'acceptation qui est crucial mais plutôt la garantie que la chaîne de Markov a effectivement convergé. Le taux
d'acceptation peut néanmoins influencer le temps qu'il faudra à la chaîne de Markov pour rejoindre sa distribution
invariante.}.}\newline

\par{Cette version de l'algorithme de MH est généralement utilisée dans la littérature concernée par l'estimation
bayésienne des modèles DSGE.}

\begin{algorithm}

\item[(1)] Maximiser le noyau postérieur par rapport à $\theta$. On obtient le mode de la densité postérieure,
$\theta^m$, et le hessien au mode qui caractérise la courbure de la densité postérieure au mode et dont l'inverse de
l'opposé, noté $\omega(\theta^m)$, approxime la variance postérieure. On pose $\Sigma = c \Omega(\theta^m)$ avec $c>0$,
$s = 1$ et $\theta^{(0)} = \theta^m$.

\item[(2)] Générer $\theta^{*}$ à partir d'une gaussienne d'espérance $\theta^{(s-1)}$ et de variance $\Sigma$

\item[(3)] Générer $u$ dans une loi uniforme entre $[0,1]$

\item[(4)] Appliquer la règle suivante~:
\[
\theta^{(s)} =
\begin{cases}
\theta^*   &\text{ si } \alpha\left(\theta^{(s-1)},\theta^*\right) > u\\
\theta^{(s-1)}&\text{ sinon.}
\end{cases}
\]
où
\[
\alpha(\theta^{(s-1)},\theta^*) = \min \left\{1,
\frac{\mathcal K\left(\theta^* | \sample\right)}{\mathcal K\left(\theta^{(s-1)} | \sample\right)}
\right\}
\]

\item[(5)] Reprendre (2-4) pour $s=2,\dots,n$.

\end{algorithm}


\par{\noindent\textbf{L'algorithme à chaînes indépendantes.} Si la proposition est indépendante de l'état
courant, l'algorithme est dit à chaines indépendantes (\citet{Tierney1994}).  La probabilité d'acceptation se
simplifie alors comme suit~:
\[
\alpha(\theta^{(s-1)},\theta^{*}) = \min \left\{1,
\frac{\mathcal K (\theta^* | \sample}{\mathcal K (\theta^{(s-1)} | \sample)}
\frac{q(\theta^{(s-1)})}{q(\theta^*)} \right\}.
\]
Cet algorithme est particulièrement adapté au cas où
il existe une approximation naturelle de la densité \emph{a
posteriori}. En effet, l'algorithme à chaînes indépendantes est
alors similaire à l'algorithme par fonction d'importance. Pour
s'en convaincre, il suffit de remarquer que l'on peut définir des
poids analogues à ceux vus précédemment comme suit~:
\[
w\left(\theta\right) = \frac{p_1(\theta^* \mid \mathcal{Y}_T,
\mathcal{A})}{q(\theta)}.
\]
La probabilité d'acceptation est alors donnée par
$\alpha(\theta^{(s-1)},\theta^{(s)}) = min \left\{1,
w(\theta^*)/w(\theta^{(s-1)})\right\}$. Autrement dit, il s'agit du
ratio des poids d'échantillonnage par importance entre le vecteur
candidat et le vecteur précédent.}

\subsubsection{Les diagnostics de convergence}\label{sec:mh:conv}

\par{Un certain nombre de résultats théoriques sont disponibles pour l'étude de la convergence des
chaînes de Markov. Il est cependant extrêmement compliqué d'énoncer des règles pratiques. Ainsi, il n'existe
aucune réponse simple à la question~: Quel est le nombre optimal de simulations pour garantir la convergence de la chaîne
de Markov vers la distribution ergodique~? Nous discutons brièvement quelques diagnostics de convergence. Le lecteur
intéressé pourra consulter, par exemple, \citet[chapitre 12]{CasellaRobertBook}.}\newline

\par{La littérature bayésienne distingue généralement trois types de convergence : (\textit{i}) la convergence vers la
distribution stationnaire, (\textit{ii}) la convergence des moments empiriques (ou approchés) vers les moments
théoriques et (\textit{iii}) la convergence vers un échantillonage i.i.d.. Nous nous intéressons ici aux deux
dernières formes de convergence\footnote{Pour plus de détails sur la convergence vers
la distribution stationnaire et l'hypothèse de stationnarité, voir \cite{GelfandSmith90}, \cite{Roberts92}, et
\cite{Liu92} }. Quelle que soit la convergence étudiée, les résultats reposent soit sur des méthodes d'évaluation
graphiques dont il est difficile d'en déduire des règles générales, soit des tests statistiques formels.
Avant de présenter certaines de ces méthodes, il convient de noter qu'il est important de distinguer les méthodes
qui font appel à $M$ chaînes de Markov parallèles et indépendantes et celles basées sur une seule chaîne
(\emph{on-line Markov chain}). L'utilisation de plusieurs chaînes est coûteuse en temps mais elle réduit la
dépendance aux conditions initiales et accroît la possibilité de parcourir efficacement l'espace des paramètres,
$\Theta$. En particulier, si une chaîne de Markov est mélangeante au sens faible --- elle reste coincée dans
des régions (de mesure dominée) de l'espace des paramètres --- une explication possible est la présence d'une
distribution \emph{a posteriori} multimodale (notamment lorsque les densités \emph{a priori} sont en conflit
avec la vraisemblance du modèle). Dans ce cas, la mise en oeuvre de chaînes de Markov en parallèle et indépendantes,
très dispersées, peut permettre de résoudre cette difficulté. L'argument des chaînes multiples est aussi avancé
pour s'assurer de la convergence. Si on se donne plusieurs vecteurs initiaux de paramètres, suffisamment dispersés,
et que l'on obtient les mêmes résultats, la convergence serait assurée. L'argument est valide si et seulement
si on a bien  vérifié que chaque chaîne a convergée. Nous revoilà donc au point de départ... Il existe une
littérature abondante discutant des avantages et inconvénients respectifs de chaque
méthode\footnote{Pour plus de détails, voir \cite{RafteryLewis1996}, \cite{CowlesCarlin96},  et
\cite{BrooksRoberts98}.}.}\newline


\par{La convergence des moments empiriques ou approchés vers les moments théoriques peut
s'appréhender à partir de méthodes graphiques. \cite{YuMykland94}
se basent sur les sommes cumulatives des moments d'intérêt dans le
cas d'une seule chaîne de Markov. Au contraire,
\cite{GelmanRubin92} proposent un test formel qui repose sur des
méthodes de chaînes de Markov en parallèle. La convergence est
diagnostiquée si les différences entre $J$ des $M$ chaines de
Markov restent dans un intervalle raisonnable.
\citeauthor{GelmanRubin92} formalisent cette idée en
recourant à des statistiques de type ANOVA. Pour chaque
statistique d'intérêt $\phi$, ils déterminent la variance intra et
inter-chaînes. L'intuition du test est alors la suivante. Si
l'effet des valeurs initiales de chaque chaîne de Markov a été
supprimé, les chaînes en parallèle doivent être relativement
proches. En d'autres termes, la variance inter-chaîne ne devrait
pas être trop grande par rapport à la variance
intra-chaîne\footnote{Plus formellement, notons
$\hat{\phi}_{n_1}^{(i)}$ l'estimateur de
$\mathbb{E}\left[\phi(\theta)\right]$ obtenu à partir du vecteur
initial $\theta^{i}$ lorsque les $n_1 \equiv n-n_0$ dernières
valeurs de la chaîne sont prises en compte. La variance intra
d'une chaîne, obtenue à partir du vecteur initial $\theta^{i}$,
est définie par $s_i= \frac{1}{n_1-1}\sum_{s=n_0+1}^{n}
\left[\phi\left(\theta^{(s,i)}\right)-\hat{\phi}_{S_1}^{(i)}
\right]^2$. La moyenne des variance-intra est alors donnée par $W
= \frac{1}{m}\sum_{i=1}^m s_i^2$ où $m$ est le nombre de chaînes
en parallèle ou de vecteurs initiaux. De la même manière, on peut
montrer que la variance-inter est estimée par $B =
\frac{n_1}{m-1}\sum_{i=1}^m \left(\hat{\phi}_{n_1}^{(i)} -
\hat{\phi}\right)^2$ où $\hat{\phi}$ est donnée par $\hat{g}=
\frac{1}{m}\sum_{i=1}^m \hat{\phi}_{n_1}^{(i)}$. Un estimateur de
la variance a posteriori de $\phi$ est alors définie comme
$\frac{n_1-1}{n_1}W+\frac{1}{n_1}B$.}. La statistique de test est
alors définie à partir d'un estimateur de la variance \textit{a
posteriori} de $\phi$. Plus précisément, ce dernier est une
moyenne pondérée de la variance intra-chaîne et de la variance
inter-chaîne. Le critère de convergence est ainsi le rapport de
cet estimateur à la variance intra-chaîne.  En utilisant une
approximation de ce critère, les auteurs montrent que si sa valeur
excède 1.2, on peut en conclure qu'il n'y a pas convergence. Dans
une autre optique, \cite{Geweke92} propose de comparer la moyenne
de deux sous-échantillons disjoints, $\mathcal S_1$ et $\mathcal
S_3$, d'une chaîne de Markov (après avoir éliminé les $n_0$
premières valeurs). On choisit $\mathcal S_1$ (resp. $\mathcal
S_3$) au début (resp. à la fin) de la chaîne de Markov. Si la
chaîne de Markov a atteint la distribution stationnaire, la
moyenne des deux sous-échantillons doit être égale. Une version
modifiée de la statistique $z$ est alors élaborée par
Geweke\footnote{Plus formellement, supposons que l'on dispose
d'une chaîne de Markov $\left(\theta^{s}\right)_{s=1,\cdots,n}$ et
que l'on subdivise cette chaîne en sous-ensembles, $\mathcal{S}_0
= \left(\theta^{s}, s=1, \cdots, n_0 \right)$, $\mathcal{S}_1 =
\left(\theta^{s}, s=n_0+1, \cdots, n_0 + n_a \right)$,
$\mathcal{S}_2 = \left(\theta^{s}, s=n_0+n_a + 1, \cdots, n_0 +
n_a + n_b\right)$ et $\mathcal{S}_3 = \left(\theta^{s},
s=n_0+n_a+n_b+1, \cdots, n_0 + n_a + n_b + n_c\right)$. On choisit
généralement, $n_a = 0.1n_1$, $n_b = 0.5n_1$ et $n_c = 0.4n_1$, où
$n_1 \equiv n-n_0$. Le test de Geweke revient à déterminer la
variance a posteriori de $\phi$, $\hat{s}_{1}$ et $\hat{s}_{3}$,
pour les sous-ensembles $\mathcal{S}_1$ et $\mathcal{S}_3$ et à
évaluer $\hat{\phi}_{\mathcal{S}_1}$ et
$\hat{\phi}_{\mathcal{S}_3}$. La statistique de test est alors
définie par
$\frac{\hat{\phi}_{\mathcal{S}_1}-\hat{\phi}_{\mathcal{S}_3}}{\frac{\hat{s}_{1}}{\sqrt{n_a}}+\frac{\hat{s}_{3}}{\sqrt{n_c}}}$
.}. Une valeur de la statistique de test supérieure à 2 indique
qu'un nombre plus élevé d'itérations est sans doute nécessaire. Le
test de \cite{RafteryLewis92a} (voir aussi \cite{RafteryLewis92b})
est plus informatif. Il se base sur les quantiles de la
statistique d'intérêt. L'idée est de construire une chaîne de
Markov à deux états  à partir d'un quantile (par exemple, 2,5\% et
97,5\%) qui permette d'estimer les probabilités de transition et
ainsi d'estimer le nombre de simulations nécessaires pour
approcher la stationnarité.}\newline

\par{Finalement, les méthodes de Monte-Carlo présentées dans les sections précédentes
ne sont valides que si les éléments de la chaîne de Markov sont
i.i.d.  Or, l'intuition suggère que les valeurs adjacentes d'une
chaîne de Markov devraient être corrélées positivement. De manière
plus générale, le fait que des autocorrélations d'ordre élevé
puissent subsister est problématique si la taille de la chaîne de
Markov n'est pas suffisamment grande. Dans cette perspective, il
est nécessaire de vérifier cette propriété ou tout du moins
d'éviter une corrélation trop élevée de la chaîne de Markov à
partir de laquelle on déduit les quantités ou statistiques
d'intérêt. Plusieurs procédures ont été suggérées dans la
littérature. Nous mentionnons ici deux stratégies. La première
repose sur un facteur correctif à appliquer à la dimension de la
chaîne de Markov en présence d'un degré observé $k$
d'autocorrélations\footnote{L'intuition repose sur un théorème
fondamental de l'analyse des séries temporelles qui nous indique
que si les $\theta^{(s)}$ sont échantillonnées à partir d'un
processus stationnaire et corrélé, les  réalisations des tirages
(qui sont donc corrélés) fournissent encore une information non
biaisée de la distribution si la taille de l'échantillon est
suffisamment grande.}. Une autre stratégie est de conserver
seulement chaque $k^{\text{ième}}$ élément de la chaîne de Markov
(après avoir supprimé les $n_0$ premiers éléments de la chaine de
Markov). Cette technique est connue sous le nom de
sous-échantillonnage, voir \cite{Schmeiser89},
\cite{RafteryLewis92a} ou \cite{RafteryLewis92b}.}

\subsection{Estimation de la densité marginale}

\par{Nous disposons d'une suite de vecteurs de paramètres $\{\theta^{(s)}\}_{s=n_0+1,\dots,n}$ où chaque $\theta^{(s)}$
est extrait de la distribution postérieure. À partir de cette suite nous pouvons estimer les moments postérieurs,
les densités prédictives, et finalement la densité marginale de l'échantillon, $p(\sample)$. Cette
densité marginale, comme nous l'avons vu en section \ref{sec:1:1}, pemet de quantifier la capacité du modèle à
expliquer l'échantillon à notre disposition et éventuellement de comparer différents modèles. Par exemple,
\citet{RabanalRubioramirez2005} évaluent différentes spécifications des rigidités nominales sur les salaires et
les prix dans le cadre d'un modèle DSGE, en comparant des densités marginales. Il existe de nombreuses méthodes pour
estimer $p(\sample)$. Dans cette section nous présentons la méthode
généralement utilisée pour les modèles DSGE.}\newline

\par{L'estimateur par la moyenne harmonique est motivé par
la propriété suivante de l'espérance postérieure~:
\[
    \mathbb E\left[\frac{f(\theta)}{p_0(\theta)p(\sample|\theta)}\right]
    = \int_{\Theta}
    \frac{f(\theta)p_1(\theta|\sample)}{p_0(\theta)p(\sample|\theta)}\mathrm d\theta
 \]
 où $f$ est une fonction de densité quelconque et $\mathbb E$ est l'espérance postérieure. Le membre de droite de
l'égalité, en utilisant la définition de la densité postérieure,
s'écrit alternativement~:
\[
  \int_{\Theta}
  \frac{f(\theta)}{p_0(\theta)p(\sample|\theta)}\frac{p_0(\theta)p(\sample|\theta)}
  {\int_{\Theta}p_0(\theta)p(\sample|\theta)\mathrm d\theta}\mathrm d\theta
  \]
On obtient donc~:
\[
  \mathbb E \left[\frac{f(\theta)}{p_0(\theta)p(\sample|\theta)}\right]
  =
  \frac{\int_{\Theta}f(\theta)\mathrm d\theta}{\int_{\Theta}p_0(\theta)
  p(\sample|\theta)\mathrm d\theta}
\]
Puisque l'intégrale de $f$ somme à un, nous obtenons finalement~:
\[
  \mathbb E \left[\frac{f(\theta)}{p_0(\theta)p(\sample|\theta)}\right]
  =
  \frac{1}{\int_{\Theta}p_0(\theta)
  p(\sample|\theta)\mathrm d\theta}
\]
Ainsi, un estimateur de la densité marginale (l'intégrale du noyau postérieur qui apparaît au dénominateur du
second membre), est l'inverse de l'espérance postérieure de $f(\theta)/\mathcal K(\theta|\sample)$. Ceci suggère
l'estimateur suivant de la densité marginale~:
\[
\hat{p}\left(\sample\right) = \frac{1}{n-n_0}\sum_{s=n_0+1}^n \frac{f\left(\theta^{(s)}\right)}
{\mathcal K\left(\theta^{(s)}|\sample\right)}
\]
Ce résultat est valable pour toute densité $f$. \citet{Geweke92} propose d'utiliser une gaussienne tronquée centrée
sur l'espérance postérieure. L'idée est d'accorder moins de poids, voire d'éliminer, les simulations éloignées
du centre de la distribution postérieure. Cela permet de diminuer la variance de l'estimateur de la densité
marginale.}\newline


\section{Un DSGE pour le B du VAR}\label{sec:bvardsge}

\par{Dans cette section, nous illustrons en quoi les modèles VAR et DSGE sont des outils complémentaires que l'on
ne doit pas nécessairement chercher à opposer. Nous avons vu dans la
section \ref{sec:1:4} que la spécification des croyances \textit{a
priori} sur la paramètrisation d'un VAR ne va pas de soi. En effet,
dans la mesure où le contenu économique d'un modèle VAR est ténu,
l'interprétation des paramètres du VAR est délicate, ce qui rend
l'élicitation des priors ardue. \cite{IngramWhiteman1994} proposent
d'utiliser un modèle DSGE afin de construire le prior d'un modèle
VAR. Ils montrent qu'en utilisant les restrictions définies par un
modèle RBC pour définir le prior d'un modèle VAR, on peut produire
avec ce dernier des prévisions comparables, en termes de précision,
à celles que nous obtiendrions avec un prior Minnesota. Ce résultat
est remarquable, car même si le modèle RBC canonique est mal
spécifié dans de nombreuses directions, il impose des restrictions
utiles pour améliorer les prévisions du VAR.}\newline

\par{Plus récemment, \cite{DelNegroSchorfheide2004} ont repris cette idée sous une forme plus simple à mettre en
oeuvre. Leur approche permet d'estimer simultanément les paramètres structurels du DSGE et les paramètres du modèle
VAR. Nous présentons leur approche dans cette section.}\newline

\subsection{Les régressions mixtes}\label{sec:bvardsge:theil}

\par{Dans la section \ref{sec:1:4:pratique} nous avons noté, au moins dans le cadre d'un modèle linéaire gaussien,
une analogie entre les priors du paradigme bayésien et les contraintes linéaires sur les paramètres de l'approche
classique. \citeauthor{DelNegroSchorfheide2004} utilisent cette analogie (voir \citet{TheilGoldberger1961},
\citet{TiaoZellner1964a} et \citet{TheilBook}) pour mettre en oeuvre le prior DSGE. Dans le modèle VAR, on peut
définir un prior sur $\mathcal A$ en utilisant des observations artificielles, cohérentes avec nos croyances, et
un prior diffus à la Jeffrey. Par exemple, si ces observations artificielles sont générées par un modèle DSGE, alors
l'estimation sur la base de l'échantillon augmenté sera attirée vers la projection du DSGE dans l'espace des VAR.}\newline

\par{Plus formellement supposons que nous disposions des observations artificielles
$(\check Y(\eta),\check Z(\eta))$, où $\eta$ est un vecteur de paramètres qui définit le processus générateur
des données artificielles (\textit{ie} les croyances \textit{a priori}). Comme l'échantillon artificiel est
indépendant de $\sample$, la vraisemblance de l'échantillon augmenté s'écrit de la façon suivante~:
\begin{equation}\label{equ:bvardsge:lik}
  p\left(\check Y(\eta),\sample|\mathcal A, \Sigma\right)
  = p\left(\check Y(\eta)|\mathcal A, \Sigma\right)\times
  p\left(\sample|\mathcal A, \Sigma\right)
\end{equation}
Le premier terme du membre de droite, si l'échantillon artificiel est de dimension $[\lambda T]$ où
$\lambda\in\mathbb R$, s'écrit~:
\begin{equation}\label{equ:bvardsge:lik:simulateddata}
   p\left(\check Y(\eta)|\mathcal A,\Sigma\right) \propto |\Sigma|^{-\frac{[\lambda T]}{2}}
   e^{-\frac{1}{2}\trace \left[\Sigma^{-1}
       \left(\check Y' \check Y - \mathcal A' \check Z' \check Y - \check Y' \check X \mathcal A +
         \mathcal A' \check Z'\check Z \mathcal A \right)\right]}
\end{equation}
et, à la lumière de l'avant dernière équation, s'interprète comme un prior pour $\mathcal A$ et $\Sigma$. La
croyance \textit{a priori} est d'autant plus informative que l'échantillon artificiel est de grande taille. Quand
$\lambda$ tend vers l'infini, le poids de la vraisemblance (le second terme sur le membre de droite de
(\ref{equ:bvardsge:lik})) devient négligeable par rapport au prior (le premier terme sur le membre de droite de
(\ref{equ:bvardsge:lik})). En complétant le prior, défini avec les données artificielles, par un prior diffus (ou plat)
à la Jeffrey~:
\[
  p_0(\mathcal A,\Sigma) \propto |\Sigma|^{-\frac{m+1}{2}}
\]
le prior est au final de type normal-Wishart, le prior conjugué dans un modèle linéaire
gaussien. En particulier, $\mathcal A$ est \textit{a priori} normalement distribué~:
\[
\VEC \mathcal A | \Sigma \sim \normal{\VEC \tilde{\mathcal A}(\eta)}{\Sigma\otimes (\check Z'\check Z)^{-1}}
\]
où $\tilde{\mathcal A}(\eta) = (\check Z'\check Z)^{-1}\check Z'\check Y$, est l'estimateur des MCO (MV) des paramètres
autorégressifs pour l'échantillon artifiel.}\newline

\par{On voit immédiatement, en considérant la vraisemblance de l'échantillon augmenté (\ref{equ:bvardsge:lik}),
le prior diffus à la Jeffrey et les résultats de la section \ref{sec:1:4:noninf}, que la distribution
\textit{a posteriori} est de type normale-Wishart~:
\begin{equation}\label{equ:var:mixed:post}
\begin{split}
  \mathcal{A}|\Sigma,\sample,\check Y (\eta) &\sim MN_{k,m}\left(\widehat{\mathcal{A}}(\eta),\Sigma,(Z'Z+
\check Z'\check Z)^{-1}\right)\\
  \Sigma|\sample,\check Y (\eta) &\sim i\mathcal{W}_m\left(\widehat{S}(\eta),\tilde{\nu}\right)
\end{split}
\end{equation}
où $\widehat{\mathcal{A}}(\eta)$ et $\widehat{S}(\eta)$ sont respectivement les estimateurs du maximum de vraisemblance
de $\mathcal A$ et $(T+[\lambda T])\Sigma$, pour l'échantillon augmenté des données artificielles, $\tilde \nu =
[(1+\lambda)T]-k$. En intégrant la densité jointe postérieure par rapport à $\Sigma$, on montre que la distribution
postérieure marginale de $\mathcal A$ est une distribution de student matricielle, centrée en
$\widehat{\mathcal{A}}(\eta)$. Lorsque $\lambda$ augmente, $\widehat{\mathcal{A}}(\eta)$ se rapproche de
$\tilde{\mathcal{A}}(\eta)$, en effet, on établit facilement que~:
\[
\widehat{\mathcal{A}}(\eta) = (\check Z'\check Z + Z'Z)^{-1}\left(\check Z'\check Z \tilde{\mathcal A}(\eta)+ Z'Z \widehat{\mathcal A}\right)
\]
Ainsi, lorsque le poids du prior augmente, la distribution postérieure du VAR se rapproche de la projection
dans l'espace des VAR du modèle générateur des données cohérent avec nos croyances \textit{a priori}.}\newline


\subsection{Le modèle BVAR-DSGE}

\par{\cite{DelNegroSchorfheide2004} proposent, à la suite de \cite{IngramWhiteman1994}, d'utiliser un modèle
DSGE pour spécifier le prior d'un modèle VAR. Contrairement à ces derniers, \citeauthor{DelNegroSchorfheide2004}
utilisent les régressions mixtes décrites dans la section \ref{sec:bvardsge:theil}, même si en pratique ils ne
simulent pas des données. Afin d'éviter que les résultats puissent varier, à cause des simulations, ils préfèrent
remplacer les moments empiriques dans (\ref{equ:bvardsge:lik:simulateddata}) par des moments théoriques calculés
à partir d'une approximation de Taylor d'ordre un de la forme réduite (\ref{equ:dsge:sol}) du modèle DSGE. Par
exemple ils remplacent $\check Y' \check Y$ par la matrice de variance-covariance des endogènes observées,
c'est-à-dire une sous matrice de $\Gamma_{yy}(\theta) = \mathbb E \left[(y_t-\mathbb E y_t)'(y_t-\mathbb E y_t)\right]$,
multipliée par la taille de l'échantillon artificiel, $[\lambda T]$. Pour tout vecteur de paramètres structurels,
$\theta \in \Theta$, la définition du prior du VAR est pratiquement immédiate, il suffit d'écrire la forme réduite
du modèle DSGE et de calculer ses moments asymptotiques. \citeauthor{DelNegroSchorfheide2004} ne se contentent pas
d'estimer les paramètres du VAR, ils estiment simultanément les paramètres du modèle DSGE. Ils spécifient donc un
prior joint sur les paramètres du modèle VAR et les paramètres structurels du modèle DSGE:
\[
p_0\left(\mathcal A, \Sigma, \theta | \lambda\right) = p_0(\mathcal A, \Sigma | \theta,\lambda)
\times p_0\left(\theta\right)
\]
Le prior est conditionnel au paramètre $\lambda$ qui spécifie la taille de l'échantillon artificiel relativement
à l'échantillon d'origine, c'est-à-dire la quantité relative d'information structurelle \textit{a priori}. On peut
alors utiliser l'algorithme de Metropolis-Hastings pour obtenir la distribution postérieure de $\theta$ (et
indirectement de $\mathcal A$ et $\Sigma$) en utilisant la densité postérieure du modèle BVAR spécifiée par
(\ref{equ:var:mixed:post}). Ici, la vraisemblance du modèle DSGE n'a pas à être calculée, ce qui simplifie
considérablement l'estimation puisque le filtre de Kalman n'est plus nécessaire. Les paramètres du modèle DSGE
sont identifiés grâce à la vraisemblance, plus exactement la densité postérieure, de son approximation VAR. Le modèle
VAR joue ici en quelque sorte le même rôle qu'un modèle auxiliaire en inférence indirecte (voir
\citet{GourierouxMonfortSimulation}).}\newline

\par{L'estimation de $\theta$ (et donc de $\mathcal A$ et $\Sigma$) est conditionnelle aux choix de $p$, le nombre de
retards dans le VAR, et $\lambda$, la quantité relative d'information structurelle \textit{a priori} dans le VAR. Il
convient de choisir un nombre de retard assez grand pour que le modèle VAR puisse être une approximation acceptable
du modèle DSGE. En effet la forme réduite (\ref{equ:dsge:sol}) approximée du modèle DSGE n'appartient pas à la
famille des modèles VAR, il faudrait un nombre de retard infini pour approximer au mieux le modèle
DSGE\footnote{Par exemple, \cite{Campbell1994}
établit, en écrivant analytiquement la forme réduite du modèle RBC linéarisé, que le produit par tête est un processus
ARMA(2,1). Ce modèle prédit donc que le produit par tête est un AR($\infty$) que l'on pourrait approximer avec
un AR($p$) pour un nombre de retards, $p$, assez grand.}. \citeauthor{DelNegroSchorfheide2004} estiment un VAR décrivant
l'inflation, le taux d'intérêt et le taux de croissance du produit. Ils affirment qu'un VAR(4) permet une approximation
satisfaisante de leur modèle DSGE. Le choix de $\lambda$ est plus délicat, en variant ce paramètre de zéro à l'infini,
on passe d'un prior diffus (l'espérance postérieure de $\mathcal A$ est alors l'estimateur du MV) à un prior très
informatif (l'espérance postérieur de $\mathcal A$ tend vers $\Gamma_{zz}(\theta)^{-1}\Gamma_{yz}(\theta)$ les contraintes
DSGE sur les paramètres du modèle VAR). \citeauthor{DelNegroSchorfheide2004} proposent d'estimer plusieurs modèles
pour une grille de valeurs de $\lambda$. Ils choisissent alors le modèle, c'est-à-dire la valeur de $\lambda$, qui
maximise la densité marginale. Ils sélectionnent le modèle dont la qualité d'ajustement est la meilleure.
\cite{DelNegroSchorfheideSmetsWouters2004} utilisent le BVAR-DSGE pour estimer le modèle de \cite{SmetsWouters2004}, ils
obtiennent $\lambda=0,75$. Ils montrent ainsi que les restrictions apportées par le modèle de Smets et Wouters
sont utiles pour améliorer les performances du modèle VAR. Cette procédure est relativement compliquée à mettre en
oeuvre. Pour chaque valeur de $\lambda$ il faut s'assurer de la convergence de l'algorithme de Metropolis-Hastings,
afin d'estimer la densité marginale\footnote{On pourrait se contenter de l'estimation du mode postérieur et d'une
approximation de Laplace, mais cette possibilité n'est pas évoquée par \citeauthor{DelNegroSchorfheide2004}.}. Plus
haut nous avons noté l'analogie entre le choix d'un modèle dans une collection de modèles et l'estimation d'un
paramètre dont les valeurs seraient discrètes. Une approche plus directe est d'associer une distribution
\textit{a priori} à $\lambda$ puis d'estimer ce paramètre (avec les paramètres structurelles $\theta$). Il faut alors
définir un prior joint sur $\mathcal A$, $\Sigma$, $\theta$ et $\lambda$~:
\[
p_0(\mathcal,\Sigma,\theta,\lambda) = p_0(\mathcal A, \Sigma | \theta,\lambda)
\times p_0\left(\theta\right) \times p_0(\lambda)
\]
\cite{AdjemianDarracqPariesTwoCountry} estiment une version deux pays du modèle de \citeauthor{SmetsWouters2004}, avec
le modèle BVAR-DSGE, en posant un prior uniforme (entre 0 et 10) pour le paramètre $\lambda$. Ils obtiennent une
distribution postérieure de $\lambda$ centrée sur $2,5$. Il n'est pas surprenant d'obtenir dans ce cas une valeur de
$\lambda$ largement supérieure. La version deux pays du modèle de Smets et Wouters est estimée avec un VAR comprenant
dix-huit variables observables, alors que \citet{DelNegroSchorfheideSmetsWouters2004}, pour la version un pays, ne
considèrent que sept variables. Avec dix-huit variables, les restrictions structurelles deviennent plus nécessaires, même
si le modèle n'est pas mieux spécifié.}

\subsection{Usages et avantages du BVAR-DSGE}

\par{\citet{DelNegroSchorfheide2004} et surtout \citet{DelNegroSchorfheideSmetsWouters2004} présentent le modèle
BVAR-DSGE comme un outil d'évaluation de la qualité d'ajustement d'un modèle DSGE. Pour ces derniers la valeur
de $\lambda$, le poids du prior structurel, donne une idée de l'intérêt empirique du modèle. Si les restrictions
structurelles définies par le modèle DSGE sont pertinentes, alors la procédure sélectionne une valeur élevée de
$\lambda$. Si le modèle apporte des informations totalement incohérentes avec les données alors la procèdure
sélectionne une valeur proche de zéro. Malheureusement cette mesure n'a pas d'échelle et nous ne savons pas à
partir de quelle valeur de $\lambda$ on peut dire que le modèle apporte des informations pertinentes. Un autre
problème est que ce paramètre ne mesure pas la qualité d'ajustement du modèle DSGE, il nous donne la quantité
optimale, au sens du \textit{fit} du modèle BVAR, d'information DSGE qu'il faut incorporer dans le prior du VAR.
\citet{DelNegroSchorfheideSmetsWouters2004} utilisent le BVAR-DSGE afin de dévoiler les éventuelles mauvaises,
spécifications d'un modèle DSGE. Or le niveau optimal de $\lambda$ ou la densité marginale, $p(\sample)$, du modèle
DSGE ne sauraient donner une idée précise des directions dans lesquels le modèle est insatisfaisant puisque ces deux
indicateurs donnent des informations trop agrégées. Les quatre auteurs recherchent les directions dans lesquels le
modèle DSGE est mal spécifié en comparant les fonctions de réponses (IRF) du modèle BVAR-DSGE avec celles du modèle
DSGE. Ils identifient les chocs structurels dans le modèle BVAR-DSGE en se fondant sur le modèle DSGE
(\ref{equ:dsge:def}). À partir de la forme réduite (\ref{equ:dsge:sol}) il est possible de calculer l'impact instantanné
de chaque choc structurel sur les variables observables~:
\[
Z\frac{\partial \mathcal H_{\theta}}{\partial \varepsilon}
\]
où $Z$ est une matrice de sélection définie dans l'équation de
mesure (\ref{equ:dsge:sol:measure}).
\citeauthor{DelNegroSchorfheide2004} utilisent cette information
pour indentifier les innovations structurelles dans le modèle BVAR
(se reporter à \cite{DelNegroSchorfheide2004} pour les détails).
Même si le BVAR-DSGE est construit sur la base (au moins
partiellement) d'une information structurelle provenant du DSGE, ce
modèle est moins contraint que le modèle DSGE. Ainsi, l'observation
d'une différence significative entre les IRFs du BVAR-DSGE et celles
du modèle DSGE conduit \citet{DelNegroSchorfheideSmetsWouters2004} à
identifer les directions dans lesquelles le modèle DSGE est mal
spécifié. Par exemple, les quatre auteurs observent que les réponses
du produit, de la consommation et des heures face à un choc de
préférence (sur la désutilité du travail) sont plus persistentes
dans le BVAR-DSGE que dans le DSGE. Ils concluent alors que le
modèle DSGE manque de mécanismes de propagation des chocs sur
l'offre de travail. Cet exercice de comparaison entre BVAR-DSGE et
DSGE peut être mis en oeuvre en utilisant des statistiques autres
que des fonctions de réponse~: décompositions de variance des
variables observées, moments théoriques des variables observées,...
La limite de l'exercice est que les conditions d'identification des
chocs dans le BVAR sont directement dérivées du modèle DSGE. Si nous
n'observons pas de grandes différences entre les IRFs du BVAR-DSGE
et celles du DSGE c'est peut être parce que nous utilisons les mêmes
conditions d'identifications. Ce problème ne se pose pas si on
compare des statistiques qui ne reposent pas sur des conditions
d'identification, par exemple si on compare des moments (variances,
fonction d'autocorrélation,...).}\newline

\par{L'avantage du modèle BVAR-DSGE est plus évident en termes de prévisions. Tout modèle est, par nature, mal spécifié
dans une multitude de directions. Malgré cette limite intrinsèque les modèles apportent souvent des informations
utiles et pertinentes. L'expérience de \cite{IngramWhiteman1994} est, à cet égard, des plus éclairantes. Ils montrent
que même le plus stylisé des modèles DSGE (le modèle de cycle réel canonique) est suffisament riche pour aider un
BVAR à produire des prévisions plus précises et moins biaisées. Un modèle que personne ne voudrait
utiliser pour produire des prévisions peut aider un modèle a-théorique (plus souple) à produire de meilleures
prévisions. Cette idée pourrait être développée dans de nombreuses directions. Nous pourrions par exemple utiliser
plusieurs modèles DSGE pour construire le prior d'un modèle VAR (ou de tout autre modèle a-théorique, par exemple
un modèle à facteurs communs) et optimiser les parts de chaque modèle dans le prior du VAR.}



\end{spacing}

\bibliography{ecoprev3}
%% À décocher pour éco et prev
%\theendnotes


\appendix

\section{Densités pour le modèle BVAR}\label{annex:A}

\subsection{Distribution normale matricielle}
 
\begin{definition}
La matrice $p \times q$ aléatoire $\mathbf{X}$ est distribuée
conformément à une loi normale matricielle
\[
\mathbf{X} \sim MN_{p,q} (\mathbf{M},\mathbf{Q},\mathbf{P})
\]
où $\mathbf{M}$ est une matrice $p \times q$, $\mathbf{Q}$ et
$\mathbf{P}$ sont respectivement des matrices $q\times q $ et $p\times p$ 
symétriques et définies positives, si et seulement si
$vec (\mathbf{X})$ est distribué comme une v.a. normale multivariée
\[
vec (\mathbf{X}) \sim \mathcal{N}_{pq}
(vec(\mathbf{M}),\mathbf{Q}\otimes\mathbf{P})
\]
\end{definition}
\noindent Ainsi, la fonction de densité associée à $\mathbf{X}$ est
donnée par~:
\[
\begin{split}
f_{MN_{p,q}}(X;\mathbf{M},\mathbf{P},\mathbf{Q}) =&
(2\pi)^{-\frac{pq}{2}}|\mathbf{Q}|^{-\frac{p}{2}}|\mathbf{P}|^{-\frac{q}{2}}\\&e^{-\frac{1}{2}\trace
\left\{\mathbf{Q}^{-1}(X-\mathbf{M})'\mathbf{P}^{-1}(X-\mathbf{M})\right\}}
\end{split}
\]

\subsection{Distributions de Wishart}
La loi de Wishart est une version multivariée de la loi du $\chi^2$. Soit $\{X_i\}_{i=1}^{\nu}$ une suite de variables aléatoires gaussiennes indépendantes et identiquement distribuées $\mathcal N(0,Q)$, avec $Q$ une matrice symétrique définie positive $q\times q$. Par définition $Y = \sum_{i=1}^{\nu} X_iX_i'$ est distribué selon une loi de Wishart. Les définitions suivantes caractérisent cette loi et la densité de l'inverse d'une v.a. de Wishart.
\begin{definition}
La matrice aléatoire, de dimension $q \times q$, symétrique et semi
définie positive $\mathbf{Y}$ est distribuée selon une loi de
Wishart, $\mathbf{Y} \sim \mathcal{W}_q (\mathbf{Q},\nu)$, si et
seulement si sa densité est donnée par
\[
f(Y;\mathbf{Q},\nu) =
\frac{|\mathbf{Q}|^{-\frac{\nu}{2}}|Y|^{\frac{\nu-q-1}{2}}}{2^{\frac{\nu
q}{2}}\pi^{\frac{q(q-1)}{4}}\prod_{i=1}^q \Gamma
\left(\frac{\nu+1-i}{2}\right)}e^{-\frac{1}{2}\trace
\left\{Y\mathbf{Q}^{-1}\right\}}
\]
pour $\mathbf{Q}$ une matrice symétrique semie définie positive, et
$\nu \leq q$ le degré de liberté.
\end{definition}

\begin{definition}
Une matrice aléatoire, de dimension $q \times q$, $\mathbf{X}$ est
distribuée selon une loi inverse Wishart,
\[
\mathbf{X} \sim i\mathcal{W}_q (\mathbf{Q},\nu)
\]
si et seulement si $\mathbf{X}^{-1} \sim \mathcal{W}_q
(\mathbf{Q}^{-1},\nu)$.
\end{definition}
\noindent Ainsi la fonction de densité associée à $\mathbf{X}$ est
définie par~:
\[
f_{i\mathcal{W}_{q}}(X;\mathbf{Q},\nu) =
\frac{|\mathbf{Q}|^{\frac{\nu}{2}}|X|^{-\frac{\nu+q+1}{2}}}{2^{\frac{\nu
q}{2}}\pi^{\frac{q(q-1)}{4}}\prod_{i=1}^q \Gamma
\left(\frac{\nu+1-i}{2}\right)}e^{-\frac{1}{2}\trace
\left\{X^{-1}\mathbf{Q}\right\}}
\]

\end{document}
