Bayesian source separation with mixture of Gaussians prior for sources and Gaussian prior for mixture coefficients
Hichem Snoussi
} \def\eabs{
beginend.tex
% TEX 7(ascii) bits
% BEGINEND.tex                          LaTeX document
% Author: J.I.          Date : 07/96    Derniere modif le 10/07/96
% Raccourcis des \begin et \end, d'apr\`es les Grands Anciens
% Modification provisoire le 12/02/97 sur \beq#1\eeq --->
% ---> Compatibilit\'e avec Amsmath.sty

\def\babs{\begin{abstract}}             \def\eabs{\end{abstract}}
\def\barr{\begin{array}}                \def\earr{\end{array}}
\def\bcc{\begin{center}}                \def\ecc{\end{center}}
\def\cl#1{\centerline{#1}}
\def\bdes{\begin{description}}          \def\edes{\end{description}}
\def\bdoc{\begin{document}}             \def\edoc{\end{document}}
\def\ben{\begin{enumerate}}             \def\een{\end{enumerate}}
\def\beqn{\begin{eqnarray}}             \def\eeqn{\end{eqnarray}}
\def\beqnl#1{\beqn\label{#1}}           \def\eeqnl#1{\label{#1}\eeqn}
\def\beqnx{\begin{eqnarray*}}           \def\eeqnx{\end{eqnarray*}}
\def\bseqn{\begin{subeqnarray}}         \def\eseqn{\end{subeqnarray}}
%
\def\beq#1\eeq{\begin{equation}#1\end{equation}}
%\def\beq{\begin{equation}}             \def\eeq{\end{equation}}
%\def\beql#1{\beq\label{#1}}            \def\eeql#1{\label{#1}\eeq}
%
\def\bal#1\eal{\begin{align}#1\end{align}}
\def\balx#1\ealx{\begin{align*}#1\end{align*}}
%
\def\beqx{$$}                           \def\eeqx{$$}
\def\bfig{\protect\begin{figure}}       \def\efig{\protect\end{figure}}
\def\bfigx{\protect\begin{figure*}}     \def\efigx{\protect\end{figure*}}
\def\bfigt{\protect\begin{figurette}}   \def\efigt{\protect\end{figurette}}
\def\bfl{\begin{flushleft}}             \def\efl{\end{flushleft}}
\def\bfr{\begin{flushright}}            \def\efr{\end{flushright}}
\def\bit{\begin{itemize}}               \def\eit{\end{itemize}}
\def\bmi{\begin{minipage}}              \def\emi{\end{minipage}}
\def\bfmi{\begin{fminipage}}            \def\efmi{\end{fminipage}}
\def\bpic{\begin{picture}}              \def\epic{\end{picture}}
\def\bqu{\begin{quote}}                 \def\equ{\end{quote}}
\def\bqun{\begin{quotation}}            \def\equn{\end{quotation}}
\def\bsl{\begin{slide}}                 \def\esl{\end{slide}}
\def\btabb{\begin{tabbing}}             \def\etabb{\end{tabbing}}
\def\btabl{\begin{table}}               \def\etabl{\end{table}}
\def\btablx{\begin{table*}}             \def\etablx{\end{table*}}
\def\btab{\begin{tabular}} %\def\etab{\end{tabular}} Conflit avec eta gras... Elle est pas grasse, ma soeur !
\def\btabu{\begin{tabular}}             \def\etabu{\end{tabular}}
\def\btabx{\begin{tabular*}}            \def\etabx{\end{tabular*}}
\def\bbib{\begin{thebibliography}{999}} \def\ebib{\end{thebibliography}}
\def\bver{\begin{verbatim}}             \def\ever{\end{verbatim}}
\def\bca{\begin{cases}}                          \def\eca{\end{cases}}

abrmath.tex
% TEX 7(ascii) bits
%
% ABRMATH.tex           LaTeX document
% Author: J.I.          Date  : Juillet 1996
% Raccourcis d'expressions tres usitees en math
% Derniere modif importante le 12-08-96


%-- Accolades, parenth�ses, etc --------------------------------------
        
\def\pth#1{\left(#1\right)}                \def\stdpth#1{(#1)}
\def\acc#1{\left\{#1\right\}}              \def\stdacc#1{\{#1\}}
\def\cro#1{\left[#1\right]}                \def\stdcro#1{[#1]}
\def\bars#1{\left|#1\right|}               \def\stdbars#1{|#1|}
\def\norm#1{\left\|#1\right\|}             \def\stdnorm#1{\|#1\|}
\def\scal#1{\left\langle#1\right\rangle}   \def\stdscal#1{\langle#1\rangle}
 
\def\bigpth#1{\bigl(#1\bigr)}              \def\biggpth#1{\biggl(#1\biggr)}
\def\bigacc#1{\bigl\{#1\bigr\}}            \def\biggacc#1{\biggl\{#1\biggr\}}
\def\bigcro#1{\bigl[#1\bigr]}              \def\biggcro#1{\biggl[#1\biggr]}
\def\bigbars#1{\bigl|#1\bigr|}             \def\biggbars#1{\biggl|#1\biggr|}
\def\bignorm#1{\bigl\|#1\bigr\|}           \def\biggnorm#1{\biggl\|#1\biggr\|}
\def\bigscal#1{\bigl\langle#1\bigr\rangle} \def\biggscal#1{\biggl\langle#1\biggr\rangle}

\def\Bigpth#1{\Bigl(#1\Bigr)}              \def\Biggpth#1{\Biggl(#1\Biggr)}
\def\Bigacc#1{\Bigl\{#1\Bigr\}}            \def\Biggacc#1{\Biggl\{#1\Biggr\}}
\def\Bigcro#1{\Bigl[#1\Bigr]}              \def\Biggcro#1{\Biggl[#1\Biggr]}
\def\Bigbars#1{\Bigl|#1\Bigr|}             \def\Biggbars#1{\Biggl|#1\Biggr|}
\def\Bignorm#1{\Bigl\|#1\Bigr\|}           \def\Biggnorm#1{\Biggl\|#1\Biggr\|}
\def\Bigscal#1{\Bigl\langle#1\Bigr\rangle} \def\Biggscal#1{\Biggl\langle#1\Biggr\rangle}


% Quelques fonctions classiques + Arguments entre [] ----------------
%
% Laisser les \mathrm entre {}, sinon �a d�conne dans le style "ieeetran". 
%
\def\diag{{\mathrm{diag}}}		\def\Diag#1{{\mathrm{diag}}\bigcro{#1}}	\def\Diagold#1{{\mathrm{diag}}\cro{#1}}
\def\tr{{\mathrm{tr}}\,}		\def\Tr#1{{\mathrm{tr}}\bigcro{#1}}			\def\Trold#1{{\mathrm{tr}}\cro{#1}}
\def\rg{{\mathrm{rg}}\,}		\def\Rg#1{{\mathrm{rg}}\bigcro{#1}}			\def\Rgold#1{{\mathrm{rg}}\cro{#1}}
\def\esp{{\mathrm{E}}\,}		\def\Esp#1{{\mathrm{E}}\bigcro{#1}}			\def\Espold#1{{\mathrm{E}}\cro{#1}} 
\def\var{{\mathrm{var}}\,}		\def\Var#1{{\mathrm{var}}\bigcro{#1}}		\def\Varold#1{{\mathrm{var}}\cro{#1}}
\def\cov{{\mathrm{Cov}}\,}		\def\Cov#1{{\mathrm{Cov}}\bigcro{#1}}		\def\Covold#1{{\mathrm{Cov}}\cro{#1}}
\def\Cos#1{\cos\bigcro{#1}}	\def\Cosold#1{\cos\cro{#1}}
\def\Sin#1{\sin\cro{#1}}		\def\Sinold#1{\sin\cro{#1}}
\def\Exp#1{\exp\bigcro{#1}}	\def\Expold#1{\exp\cro{#1}}
\def\Log#1{\log\bigcro{#1}}	\def\Logold#1{\log\cro{#1}}
\def\Ln#1{\ln\cro{#1}}			\def\Det#1{\det\bigcro{#1}}					\def\Detold#1{\det\cro{#1}}

\def\cond{\textrm{Cond}\,}		\def\Cond#1{\cond\bigpth{#1}}

\def\sinc{{\mathrm{sinc}}\,}	\def\Sinc#1{{\mathrm{sinc}}\bigcro{#1}}	\def\Sincold#1{{\mathrm{sinc}}\cro{#1}}
\def\rang{{\mathrm{rang}}\,}	\def\Rang#1{\rang\bigcro{#1}}					\def\Rangold#1{\rang\cro{#1}}
\def\ker{\textrm{Ker}\,}		\def\Ker#1{\ker\bigcro{#1}}					\def\Kerold#1{\ker\cro{#1}}
\def\img{\textrm{Im}\,}			\def\Img#1{\img\bigcro{#1}}					\def\Imgold#1{\img\cro{#1}}
\def\vect{{\mathrm{Vect}}\,}	\def\Vect#1{\vect\bigcro{#1}}					\def\Vectold#1{\vect\cro{#1}}
\def\sgn{{\mathrm{sgn}}}		\def\Sgn#1{\sgn\bigcro{#1}}					\def\Sgnold#1{\sgn\cro{#1}}

% Signe proba
\def\Pr{\mathop{\textrm{Pr}}}

% Variante pour \Re:
\def\reel{{\mathrm{Re}}} 		\def\Reel#1{{\mathrm{Re}}\cro{#1}}
\def\card{{\mathrm{Card}}}		\def\Card#1{{\mathrm{Card}}\cro{#1}}

\def\sh{{\mathrm{sh}}}                % sin hyperbolique en FR
\def\ch{{\mathrm{ch}}}                % cos hyperbolique en FR
\def\th{{\mathrm{th}}}                % th hyperbolique en FR
\def\coth{{\mathrm{coth}}}            % cot hyperbolique en FR
\def\coeffbin#1#2{\pth{\setlength{\arraycolsep}{.1em}\barr{c}#1\\#2\earr}}
% Remarque : \exist \binom{}{}

% Op�rateurs sur champs
\def\Div{{\mathrm{div}}}											% divergence
\def\Rotv{\overrightarrow{\mathop{{\mathrm{rot}}}}}		% rotationnel avec fleche
\def\Gradv{\overrightarrow{\mathop{{\mathrm{grad}}}}}		% gradient avec fleche
%-- Textes (if, si,...) droit en math ---------------------------

\def\IF{\text{if\:}}             \def\SI{\text{si\:}}
\def\If{\text{If\:}}             \def\Si{\text{Si\:}}
\def\AND{\text{and\:}}           \def\ET{\text{et\:}}
\def\OR{\text{or\:}}             \def\OU{\text{ou\:}}
\def\THEN{\text{then\:}}         \def\ALORS{\text{alors\:}}
                                 \def\DOU{\text{d'o\`u\:}}
\def\WHERE{\text{where\:}}       \def\Ou{\text{o\`u\:}}
\def\WHEN{\text{when\:}}         \def\QUAND{\text{quand\:}}
\def\FOR{\text{for\:}}           \def\POUR{\text{pour\:}}
\def\FORALL{\text{for all\:}}    \def\POURTOUT{\text{pour tout\:}}
\def\ST{\text{s.t.\:}}           \def\SC{\text{s.c.\:}}
\def\SUBJTO{\text{subject to\:}} \def\SOUSC{\text{sous contraintes\:}}
\def\OTHERWISE{\text{otherwise}} \def\SINON{\text{sinon}}
\def\WITH{\text{with\:}}         \def\AVEC{\text{avec\:}}
\def\IN{\text{in\:}}             \def\DANS{\text{dans\:}}

%-- Triture param�tres tableaux... -----------------------------

\def\arrayp{\renewcommand{\arraystretch}{.7}\setlength{\arraycolsep}{2pt}}
\def\tabp{\renewcommand{\arraystretch}{.7}\setlength{\tabcolsep}{2pt}}

%-- Environnement encadr\'e ------------------------------------
% Ces d\'efinitons introduisent un nouvel environnement appele "fminipage". Il 
% permet d'encadrer une partie de texte.

\newsavebox{\fminibox}
\newlength{\fminilength}
\newenvironment{fminipage}[1][\linewidth]
  {\setlength{\fminilength}{#1}%-2\fboxsep-2\fboxrule}%
   \begin{lrbox}{\fminibox}\begin{minipage}{\fminilength}}
  {\end{minipage}\end{lrbox}\noindent\fbox{\usebox{\fminibox}}}

%-- Divers : "^{-1}", "dag", "," ... %---------------------------

\def\M{^{-1}} \def\T{^\tD} \def\+{^\dagger}
\def\I{\,|\,}           % "sachant" bien espac\'e pour les formules
\def\J{\mathop{\,;\,}}  % "point virgule" bien espac\'e pour les formules.
\def\w{,\thinspace}
\def\W{,\thickspace}
\def\ldotsv{,\,\ldots,\,}
\def\V{,}               % Virgule pour les nombres decimaux
%\def\V{.}              % Point pour les anglo-saxons
\def\e#1{.10^{#1}}      % Notation scientifique a la francaise
%\def\e#1{\rm{e}^{#1}}  % Notation scientifique anglo-saxonne

\def\egdef{\stackrel{\Delta}{=}}
\def\nequiv{\not\kern-.05em\equiv}
\def\egal{\kern-.5em=\kern-.5em}        % Moins d'espace autour de "="
\def\propt{\kern-.2em\propto\kern-.2em} % Idem
\def\dans{\in\!}                        % Trop d'espace apres "appartient a"
\def\pourtt{\forall\,}                  % Pas assez d'espace
\def\wh#1{\widehat{#1}}                 % Sombrero !
\def\wt#1{\widetilde{#1}} 

%\def\arg{\mathop{\mathrm{arg}}}
\def\argmax{\mathop{\mathrm{arg\,max}}} % Mieux que \def\argmax{\arg\max}
\def\argmin{\mathop{\mathrm{arg\,min}}} % car l'indice est reparti
\def\EL{\mathrm{EL}}
\def\Vect{\mathop{\text{Vect}}}

\def\froc#1#2{{#1/#2}}                  % Frac en toc
\def\fric#1#2{\frac1{#2}#1}
\def\fracds#1#2{\frac{\displaystyle#1}{\displaystyle#2}}
\def\diff#1#2{{\frac{d#1}{d#2}}}
\def\dd{\,d}                            % doit etre en italique en anglais
%\def\dd{\,\dD}                         % doit etre droit en francais
\def\derpar#1#2{{\frac{\partial #1}{\partial #2}}}
\def\derpor#1#2{{\froc{\partial #1}{\partial #2}}}
\def\parsec#1#2#3{{\frac{\partial^2 #1}{\partial #2\,\partial #3}}}
\def\parsecd#1#2{{\frac{\partial^2 #1}{{\partial #2}^2}}}
\def\porsec#1#2#3{{\froc{\partial^2 #1}{\partial #2\,\partial #3}}}
\def\porsecd#1#2{{\froc{\partial^2 #1}{{\partial #2}^2}}}
\def\intdouble{\int\kern-0.3em\int}
\def\inttriple{\int\kern-0.3em\int\kern-0.3em\int}
\def\prods{\mathop{\text{\footnotesize$\displaystyle\prod$}}}

% Pour mettre un \No{} au-dessus d'une lettre (\TM, 02/97) (package{amsmath.sty})
\def\rond#1{\overset{\kern-0.33em~_\circ}{#1}}
\def\rondit[#1]#2{\overset{\kern#1~_\circ}{#2}}

% Pour mettre un texte dans un rond : (utilise pstricks ; \JI, 08/97)
\def\incirc#1{\pscirclebox[framesep=1pt]{\scriptsize#1}}
\def\incircp#1{\pscirclebox[framesep=1pt]{\tiny#1}}
\def\outcirc#1{\raisebox{-3pt}{\huge$\times$}\hspace*{-.45cm}\text{\incirc{#1}}}




abrege.tex
% TEX 7(ascii) bits
%
% ABREGE.tex
% LaTeX document
% Author: J.I. 
% Date  : Juillet 1996
%
% Raccourcis d'expressions tr\`es usit\'ees
%
% Derniere modif le 10-07-96
%

% Raccourci pour \XS
\def\XS{\xspace}

%%%%%%%%%%%%%%%%%%% Pour definir nouvelle commande de changement de taille
\RequirePackage{ifthen}
\RequirePackage{calc}

\newcommand{\taille}[1][\scad]{%
\ifthenelse{#1 = -5}{}{}%
\ifthenelse{#1 = -4}{\tiny}{}%
\ifthenelse{#1 = -3}{\scriptsize}{}%
\ifthenelse{#1 = -2}{\footnotesize}{}%
\ifthenelse{#1 = -1}{\small}{}%
\ifthenelse{#1 = 0}{\normalsize}{}%
\ifthenelse{#1 = 1}{\large}{}%
\ifthenelse{#1 = 2}{\Large}{}%
\ifthenelse{#1 = 3}{\LARGE}{}%
\ifthenelse{#1 = 4}{\huge}{}%
\ifthenelse{#1 = 5}{\Huge}{}}

\newboolean{@serif}
\setboolean{@serif}{true}
\def\scad{-5} % scadefaultsize

\newcounter{taille}
\newcommand{\sca}[2][\scad]{\setcounter{taille}{#1}%
  \ifthenelse{\boolean{@serif}}
  {{\taille[\thetaille]\textsc{#2}}}
  {\setcounter{taille}{\value{taille}-1}{\uppercase{\taille[\thetaille]#2}}}}

% Mise en commentaire rapide !
\def\rem#1{}                    % rem{bla bla bla}

% References automatiquement entre parenthese
\def\reff#1{\XS(\ref{#1})}
\def\refs#1{section~\ref{#1}}
\def\refS#1{Section~\ref{#1}}

\newcommand{\clearemptydoublepage}{\newpage{\pagestyle{empty}\cleardoublepage}}  
% Num�ro (pour biblio)
\def\numero{n�}
\def\Numero{N�}

% Majuscules
\def\UP#1{\uppercase{#1}}       % majuscules
\def\up#1{\uppercase{#1}}       % majuscules

% Liste des abreges
\def\adhoc{\textit{ad hoc}\XS}
\def\apost{\textit{a posteriori}\XS}
\def\Apost{\textit{A posteriori}\XS}
\def\map{\textsc{map}\XS}
\def\post{\textit{posterior}\XS}
\def\Post{\textit{Posterior}\XS}
\def\afortio{\textit{a fortiori}\XS}
\def\Afortio{\textit{A fortiori}\XS}
\def\aprio{\textit{a priori}\XS}
\def\Aprio{\textit{A priori}\XS}
\def\Cad{C'est-\`a-dire\XS}
\def\cad{c'est-\`a-dire\XS}
\def\cf{\textit{cf.}\XS}
\def\cv{\textit{curriculum vit\ae}\XS}
\def\Cv{\textit{Curriculum vit\ae}\XS}
\def\CV{\textit{Curriculum Vit\ae}\XS}
\def\eg{\textit{e.g.,}\XS}
\def\etal{\textit{et al.}\XS}
\def\etc{\textit{etc\dots\XS}}
\def\etcoll{{\textit et coll.}\XS}
\def\idem{\textit{idem}\XS}
\def\ie{\textit{i.e.,}\XS}
\def\IFF{\textit{if and only if}\XS}
\def\iid{\textit{i.i.d.}\XS}
\def\infine{\textit{in fine}\XS}
\def\oeuvre{{\oe}uvre\XS}
\def\oeuvres{{\oe}uvres\XS}
\def\recto{\textit{recto}\XS}
\def\Toeplitz{T{\oe}plitz\XS}
\def\vav{vis-\`a-vis\XS}
\def\via{\textit{via}\XS}
\def\vic{\textit{vice-versa}\XS}
\def\verso{\textit{verso}\XS}
%\def\vitae{\textit{vit\ae}\XS}
\def\vs{\textit{vs}\XS}
\def\versus{\textit{versus}\XS}
\def\ssi{\textit{si et seulement si}\XS}
\def\strsen{\textit{stricto sensu}\XS}
\def\wrt{w.r.t.\XS}


\newcommand{\art}[1][\scad]{\sca[#1]{art}\XS}
\newcommand{\cnd}[1][\scad]{\sca[#1]{cnd}\XS}
\newcommand{\ndt}[1][\scad]{\sca[#1]{ndt}\XS}
\newcommand{\gnc}[1][\scad]{\sca[#1]{gnc}\XS}
\newcommand{\mfa}[1][\scad]{\sca[#1]{mfa}\XS}

\newcommand{\icm}[1][\scad]{\sca[#1]{icm}\XS}
\newcommand{\icd}[1][\scad]{\sca[#1]{icd}\XS}
\newcommand{\rsb}[1][\scad]{\sca[#1]{rsb}\XS}
\newcommand{\snr}[1][\scad]{\sca[#1]{snr}\XS}

\newcommand{\dess}[1][\scad]{\sca[#1]{dess}\XS}
\newcommand{\dea}[1][\scad]{\sca[#1]{dea}\XS}
\newcommand{\deug}[1][\scad]{\sca[#1]{deug}\XS}
\newcommand{\deust}[1][\scad]{\sca[#1]{deust}\XS}
\newcommand{\iut}[1][\scad]{\sca[#1]{iut}\XS}
\newcommand{\iup}[1][\scad]{\sca[#1]{iup}\XS}

\newcommand{\eea}[1][\scad]{\sca[#1]{eea}\XS}
\def\EEA{\'Electronique \'Electrotechnique Automatique\XS}
\newcommand{\ats}[1][\scad]{\sca[#1]{ats}\XS}
\def\ATS{Automatique et Traitement du Signal\XS}

\newcommand{\te}[1][\scad]{\sca[#1]{te}\XS}
\def\TE{travaux d'�tudes\XS}
\newcommand{\td}[1][\scad]{\sca[#1]{td}\XS}
\def\TD{travaux dirig�s\XS}
\newcommand{\tp}[1][\scad]{\sca[#1]{tp}\XS}
\def\TP{travaux pratiques\XS}
%\newcommand{\tdao}[1][\scad]{\sca[#1]{tdao}\XS}
\def\TDAO{travaux dirig�s assist�s sur ordinateur\XS}
\newcommand{\tdao}[1][\scad]{\sca[#1]{tdao}\XS}
\def\ts{traitement du signal\XS}
\def\TS{Traitement du Signal\XS}

\newcommand{\gbm}[1][\scad]{\sca[#1]{gbm}\XS}
\newcommand{\irmf}[1][\scad]{\sca[#1]{irm-f}\XS}
\newcommand{\irm}[1][\scad]{\sca[#1]{irm}\XS}
\newcommand{\rmn}[1][\scad]{\sca[#1]{rmn}\XS}
\newcommand{\eeg}[1][\scad]{\sca[#1]{eeg}\XS}
\newcommand{\meg}[1][\scad]{\sca[#1]{meg}\XS}

\newcommand{\pc}[1][\scad]{\sca[#1]{pc}\XS}
% Operating Systemes 
\newcommand{\linux}[1][\scad]{\sca[#1]{Linux}\XS}
\newcommand{\unix}[1][\scad]{\sca[#1]{Unix}\XS}
\newcommand{\msdos}[1][\scad]{\sca[#1]{msdos}\XS}
\newcommand{\sunos}[1][\scad]{\sca[#1]{Sunos}\XS}
\newcommand{\hupux}[1][\scad]{\sca[#1]{hp-ux}\XS}
\newcommand{\vms}[1][\scad]{\sca[#1]{vms}\XS}
\newcommand{\solaris}[1][\scad]{\sca[#1]{Solaris}\XS}
\newcommand{\winnt}[1][\scad]{\sca[#1]{Windows-nt}\XS}
\newcommand{\vmware}[1][\scad]{\sca[#1]{vmware}\XS}
\newcommand{\win}[1][\scad]{\sca[#1]{Windows}\XS}

% Calcul scientifique
\newcommand{\labview}[1][\scad]{\sca[#1]{Labview}\XS}
\newcommand{\matlab}[1][\scad]{\sca[#1]{Matlab}\XS}
\newcommand{\maple}[1][\scad]{\sca[#1]{Maple}\XS}
\newcommand{\mathe}[1][\scad]{\sca[#1]{Mathematica}\XS}
\newcommand{\idl}[1][\scad]{\sca[#1]{Idl}\XS}

% langages
\newcommand{\cmex}[1][\scad]{\sca[#1]{Cmex}\XS}
\newcommand{\pascal}[1][\scad]{\sca[#1]{Pascal}\XS}
\newcommand{\java}[1][\scad]{\sca[#1]{Java}\XS}
\newcommand{\lisp}[1][\scad]{\sca[#1]{Lisp}\XS}
\newcommand{\forts}[1][\scad]{\sca[#1]{Fortran}{\footnotesize 77}\XS}
\newcommand{\fortq}[1][\scad]{\sca[#1]{Fortran}{\footnotesize 90}\XS}
\newcommand{\ada}[1][\scad]{\sca[#1]{Ada}\XS}
\newcommand{\eiffel}[1][\scad]{\sca[#1]{Eiffel}\XS}
\newcommand{\tcsh}[1][\scad]{\sca[#1]{tcsh}\XS}
\newcommand{\bash}[1][\scad]{\sca[#1]{bash}\XS}
\newcommand{\lgc}[1][\scad]{\sca[#1]{C}\XS}
\newcommand{\cpp}[1][\scad]{\sca[#1]{C}{\footnotesize ++}\XS}

\def\crec{charg� de recherches\XS}
\def\Crec{Charg� de recherches\XS}
\def\dr{directeur de recherches\XS}
\def\Dr{Directeur de recherches\XS}


% congres
\newcommand{\gretsi}[1][\scad]{\sca[#1]{gretsi}\XS}
\newcommand{\maxent}[1][\scad]{\sca[#1]{maxent}\XS}
\newcommand{\icip}[1][\scad]{\sca[#1]{icip}\XS}
\newcommand{\icassp}[1][\scad]{\sca[#1]{icassp}\XS}
\newcommand{\spie}[1][\scad]{\sca[#1]{spie}\XS}
\newcommand{\gdr}[1][\scad]{\sca[#1]{gdr--prc~isis}\XS}
\newcommand{\rfia}[1][\scad]{\sca[#1]{rfia}\XS}
\newcommand{\eusipco}[1][\scad]{\sca[#1]{eusipco}\XS}

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% Ajout de famouss.tex
% modif le 25/5/2000 Phiphi
% GPInuts - Famous people - VIP
%
\newcommand{\beaunom}[4][0]{{\taille[#1]#2 \uppercase{#3}}%
{\sca[#1]{#4}}\XS}

%\newcommand{\MA}[1][\scad]{{\taille[#1]Marc \uppercase{A}}% 
%{\sca[#1]{llain}\XS}}
\newcommand{\MA}[1][\scad]{\beaunom[#1]{Marc}{A}{llain}}
\newcommand{\AB}[1][\scad]{\beaunom[#1]{Amandine}{B}{lanc}}
\newcommand{\ABF}[1][\scad]{\beaunom[#1]{Amandine}{B}{lanc~--~Fusco}}
\newcommand{\AH}[1][\scad]{\beaunom[#1]{Alain}{H}{erment}}
\newcommand{\AMD}[1][\scad]{\beaunom[#1]{Ali}{M}{ohammad--Djafari}}
\newcommand{\AR}[1][\scad]{\beaunom[#1]{Andrea}{R}{idolfi}}
\newcommand{\CH}[1][\scad]{\beaunom[#1]{Christian}{H}{einrich}}
\newcommand{\CS}[1][\scad]{\beaunom[#1]{Charles}{S}{oussen}}
\newcommand{\FC}[1][\scad]{\beaunom[#1]{Fr\'ed\'eric}{C}{hampagnat}}
\newcommand{\FDU}[1][\scad]{\beaunom[#1]{Fr\'ed\'eric}{D}{ublanchet}}
\newcommand{\GDE}[1][\scad]{\beaunom[#1]{Guy}{D}{emoment}}
\newcommand{\GLB}[1][\scad]{\beaunom[#1]{Guy Le}{B}{esnerais}}
\newcommand{\HC}[1][\scad]{\beaunom[#1]{Herv\'e}{C}{arfantan}}
\newcommand{\JFB}[1][\scad]{\beaunom[#1]{Jean--Fran\c{c}ois}{B}{ercher}}
\newcommand{\JFG}[1][\scad]{\beaunom[#1]{Jean--Fran\c{c}ois}{G}{iovannelli}}
\newcommand{\JI}[1][\scad]{\beaunom[#1]{J\'er\^ome}{I}{dier}}
\newcommand{\MF}[1][\scad]{\beaunom[#1]{Marc}{F}{ayolle}}
\newcommand{\Mimi}[1][\scad]{\beaunom[#1]{Myriam}{F}{iani}}
\newcommand{\MN}[1][\scad]{\beaunom[#1]{Mila}{N}{ikolova}}
\newcommand{\PC}[1][\scad]{\beaunom[#1]{Philippe}{C}{iuciu}}
\newcommand{\SB}[1][\scad]{\beaunom[#1]{St\'ephane}{B}{rette}}
\newcommand{\RB}[1][\scad]{\beaunom[#1]{Redha}{B}{oubertakh}}
\newcommand{\TM}[1][\scad]{\beaunom[#1]{Thierry}{M}{artin}}
\newcommand{\VS}[1][\scad]{\beaunom[#1]{Vincent}{S}{amson}}
\newcommand{\YG}[1][\scad]{\beaunom[#1]{Yves}{G}{oussard}}
\newcommand{\KS}[1][\scad]{\beaunom[#1]{Ken}{S}{auer}}
\newcommand{\HS}[1][\scad]{\beaunom[#1]{Hichem}{S}{noussi}}

% Acronymes et raccourcis
%
\def\ATER{Attach\'e temporaire d'enseignement et de recherche}
\newcommand{\ater}[1][\scad]{\sca[#1]{ater}\xspace}
\def\ADDOC{Agir pour les Doctorants et les jeunes Docteurs\XS}
\newcommand{\addoc}[1][\scad]{\sca[#1]{addoc}\XS}

\def\CdO{Centre d'Orsay\XS}

\def\CEA{Commissariat � l'\'Enegie Atomique\XS}
\newcommand{\cea}[1][\scad]{\sca[#1]{cea}\XS}
\newcommand{\dam}[1][\scad]{\sca[#1]{dam}\XS}
\newcommand{\leti}[1][\scad]{\sca[#1]{leti}\XS}
\def\SHFJ{Service Hospitalier Fr�d�ric Joliot\XS}
\newcommand{\shfj}[1][\scad]{\sca[#1]{shfj}\XS}

\newcommand{\chu}[1][\scad]{\sca[#1]{chu}\XS}
\def\CHU{Centre Hospitalier Universitaire\XS}
\def\salpe{Salp�tri�re\XS}
\def\pisalpe{Piti�-Salp�tri�re\XS}
\def\gpisalpe{Groupe Piti�-Salp�tri�re\XS}

\def\CNRS{Centre National de la Recherche Scientifique\XS}
\newcommand{\cnrs}[1][\scad]{\sca[#1]{cnrs}\XS}

\newcommand{\edf}[1][\scad]{\sca[#1]{edf}\XS}
\def\EDF{\'Electricit� de France\XS}

\def\ENSTA{\'Ecole Nationale Sup\'erieure de Techniques Avanc\'ees\XS}
\newcommand{\ensta}[1][\scad]{\sca[#1]{ensta}\XS}
\def\ESIEA{\'Ecole Sup�rieure d'Informatique \'Electronique Automatique\XS}
\newcommand{\esiea}[1][\scad]{\sca[#1]{esiea}\XS}
\newcommand{\ensea}[1][\scad]{\sca[#1]{ensea}\XS}


\def\GPI{Groupe Probl\`emes Inverses\XS}
\newcommand{\gpi}[1][\scad]{\sca[#1]{gpi}\XS}

\def\INSERM{Institut National pour la Sant\'e et la Recherche M\'edicale\XS}
\newcommand{\inserm}[1][\scad]{\sca[#1]{inserm}\XS}

\def\IAS{Institut d'Astrophysique Spatiale\XS}
\newcommand{\ias}[1][\scad]{\sca[#1]{ias}\XS}

\newcommand{\ipsn}[1][\scad]{\sca[#1]{ipsn}\XS}
\def\IPSN{Institut de Protection et de Suret� Nucl�aire\XS}

\def\LSS{Laboratoire des Signaux et Syst\`emes\XS}
\newcommand{\lss}[1][\scad]{\sca[#1]{l{\footnotesize 2}s}\XS}

\def\LIMSI{Laboratoire d'Informatique pour la M�canique et les Sciences de l'Ing�nieur\XS}
\newcommand{\limsi}[1][\scad]{\sca[#1]{limsi}\XS}

\def\ONERA{Office National d'\'Etudes et de Recherches A\'erospatiales\XS}
\newcommand{\onera}[1][\scad]{\sca[#1]{on\'era}\XS}

\def\MDC{Ma\^\i tre de conf\'erences\XS}
\def\MDCs{Ma\^\i tres de conf\'erences\XS}
\newcommand{\mdc}[1][\scad]{\sca[#1]{mdc}\XS}

\def\UPS{Universit\'e de Paris-Sud\XS}
\newcommand{\ups}[1][\scad]{\sca[#1]{ups}\XS}
\newcommand{\upxi}[1][\scad]{\sca[#1]{upxi}\XS}

% Unites Thomson
\newcommand{\thomair}[1][\scad]{\sca[#1]{Thomson-csf Airsys}\XS}
\newcommand{\oscar}[1][\scad]{\sca[#1]{Oscar}\XS}
\newcommand{\thomrcm}[1][\scad]{\sca[#1]{Thomson-rcm}\XS}

% Unite propre de recherche
\newcommand{\upr}[1][\scad]{\sca[#1]{upr}\XS}
% LENA, lena upr640,...
\newcommand{\uprlen}[1][\scad]{\sca[#1]{upr}{\footnotesize 640}\XS}
\def\LENA{Laboratoire d'Enc�phalographie et de Neurophysiologie Appliqu�e\XS}
\newcommand{\lena}[1][\scad]{\sca[#1]{lena}\XS}
% Nous
\def\GsY{Gif--sur--Yvette\XS}
\newcommand{\ESE}[1][\scad]{\sca[#1]{sup\'elec}\XS}

\def\adrshort{\ESE, Plateau de Moulon, 91192 Gif--sur--Yvette Cedex\XS}
%\def\adresse{\ESE, Plateau de Moulon, 91192 Gif--sur--Yvette Cedex, France\XS}
\newcommand{\adresse}[1][\scad]{\ESE[#1], Plateau de Moulon, 91192 Gif--sur--Yvette Cedex, France\XS}
\newcommand{\tutelle}[1][\scad]{\cnrs[#1]~--~\ESE[#1]~--~\ups[#1]}
%\newcommand{\tutelle}[1][\scad]{\sca[#1]{Cnrs~--~Sup\'elec~--~Ups}}
\def\ALSS{\LSS~(\tutelle)\\\adresse}
\def\ALSSun{\LSS (\tutelle) \adresse}
%


% ---------------------------------------------------- 
% Modele d'en-tete par A. Mohammad-Djafari et J. Idier 
% ---------------------------------------------------- 
%
\def\helvnormal#1{{\fontfamily{phv}\fontsize{11}{14pt}\selectfont #1}}
\def\helvfoot  #1{{\fontfamily{phv}\fontsize{9}{14pt} \selectfont #1}}
\def\helvscript#1{{\fontfamily{phv}\fontsize{8}{14pt} \selectfont #1}}
%
\newcommand{\entete}[5][Nom.Pr\'enom]{
%%%%%%\def\entete#1#2#3#4#5{%
\helvscript{\renewcommand{\arraystretch}{.8}
\begin{tabular}{c} Sup�lec
\\--------\end{tabular} \hfill
\begin{tabular}{c}Centre national de la recherche scientifique          \\
--------\end{tabular} \hfill
\begin{tabular}{c}\UPS  \\--------\end{tabular}} \\[-.4cm]%
%
\begin{center}\helvnormal{\textbf{LABORATOIRE DES SIGNAUX \& SYST\`EMES}}\\%
\helvfoot{Unit\'e mixte de recherche n$^{\rm o}$ 8506\\%
Sup\'elec, plateau de Moulon, 
3 rue Joliot-Curie, 91192 \sca[-3]{Gif--sur--Yvette} Cedex (France)\\%
T\'el\'ephone\,: 01 69 85 17 12 --- 
T\'el\'ecopie\,: 01 69 85 17 65 --- 
Courriel\,: {#1}@lss.supelec.fr}%
\\[0.5cm]\end{center}% 
\begin{tabular}{p{8.5cm}l}%
\helvscript{V/R\'ef.\,: #2} & \helvscript{\sca[-3]{Gif}, le} #5\\%
\helvscript{N/R\'ef.\,: #3} \\%
\helvscript{Objet\,: #4}%
\end{tabular}
} 

\newcommand{\inhead}[5][Name.Surname]{
%%%%%%\def\entete#1#2#3#4#5{%
\helvscript{\renewcommand{\arraystretch}{.8}
\begin{tabular}{c} Sup�lec
\\--------\end{tabular} \hfill
\begin{tabular}{c}Centre national de la recherche scientifique  \\
--------\end{tabular} \hfill
\begin{tabular}{c}\UPS  \\--------\end{tabular}} \\[-.4cm]%
%
\begin{center}\helvnormal{\textbf{LABORATOIRE DES SIGNAUX \& SYST\`EMES}}\\%
\helvfoot{Unit\'e mixte de recherche n$^{\rm o}$ 8506\\%
Sup\'elec, plateau de Moulon, 
3 rue Joliot-Curie, 91192 \sca[-3]{Gif--sur--Yvette} Cedex (France)\\%
Telephone: 01 69 85 17 12 --- 
Fax: 01 69 85 17 65 --- 
E-mail: {#1}@lss.supelec.fr}%
\\[0.5cm]\end{center}% 
\begin{tabular}{p{8.5cm}l}%
\helvscript{Y/Ref.: #2} & \helvscript{\sca[-3]{Gif}, } #5\\%
\helvscript{O/Ref.: #3} \\%
\helvscript{Subject: #4}%
\end{tabular}
} 


alphabet.tex
% Document LaTeX        Auteur : Guy Le Besnerais       Date : 25/05/94
% Modification : J.-F. Gio      (15/09/95)
%                               J.I. (format Latex 2e, taille adaptative 22/03/99)
%                                    (fonte blackboard 27/08/99)

% Objet : Definitions compactes d'alphabets pour ecriture en Math
% b pour bold italique : \Ab
% v pour bold droit : \Av
% D pour droit : \AD
% c pour cursif : \Ac
% bb pour blackboard : \Abb
% Capitale pour capitale : \Lambda

\typeout{\space}
\typeout{\space\space\space\space Fichier 'alphabet.tex' -- LSS 1994/95/99}
\typeout{\space\space\space\space (necessite les packages xspace, amsmath, bbm)}
\typeout{\space}

\RequirePackage{amsmath}
\RequirePackage{xspace}
\RequirePackage{bbm}

%-------------------------------------------

% DEFINITIONS PRELIMINAIRES 

\def\XS{\xspace}
\DeclareMathAlphabet{\mathb}{OML}{cmm}{b}{it}
\def\sbm#1{\ensuremath{\mathb{#1}}}                % Style gras italique (necessite amsmath)           
%
\def\sbmm#1{\ensuremath{\boldsymbol{#1}}}          % Style gras italique (necessite amsmath)           
\def\sdm#1{\ensuremath{\mathrm{#1}}}               % Style droit en math
\def\sbv#1{\ensuremath{\mathbf{#1}}}               % Style gras droit
\def\scu#1{\ensuremath{\mathcal{#1\XS}}}           % Style cursif
\def\scb#1{\ensuremath{\boldsymbol{\mathcal{#1}}}} % Style gras cursif
\def\sbl#1{\ensuremath{\mathbbm{#1}}}              % Style blackboard (necessite bbm)

%-------------------------------------------

% ALPHABET GRAS ITALIQUE, taille adaptative

\def\Ab{{\sbm{A}}\XS}  \def\ab{{\sbm{a}}\XS}
\def\Bb{{\sbm{B}}\XS}  \def\bb{{\sbm{b}}\XS}
\def\Cb{{\sbm{C}}\XS}  \def\cb{{\sbm{c}}\XS}
\def\Db{{\sbm{D}}\XS}  \def\db{{\sbm{d}}\XS}
\def\Eb{{\sbm{E}}\XS}  \def\eb{{\sbm{e}}\XS}
\def\Fb{{\sbm{F}}\XS}  \def\fb{{\sbm{f}}\XS}
\def\Gb{{\sbm{G}}\XS}  \def\gb{{\sbm{g}}\XS}
\def\Hb{{\sbm{H}}\XS}  \def\hb{{\sbm{h}}\XS}
\def\Ib{{\sbm{I}}\XS}  \def\ib{{\sbm{i}}\XS}
\def\Jb{{\sbm{J}}\XS}  \def\jb{{\sbm{j}}\XS}
\def\Kb{{\sbm{K}}\XS}  \def\kb{{\sbm{k}}\XS}
\def\Lb{{\sbm{L}}\XS}  \def\lb{{\sbm{l}}\XS}
\def\Mb{{\sbm{M}}\XS}  \def\mb{{\sbm{m}}\XS}
\def\Nb{{\sbm{N}}\XS}  \def\nb{{\sbm{n}}\XS}
\def\Ob{{\sbm{O}}\XS}  \def\ob{{\sbm{o}}\XS}
\def\Pb{{\sbm{P}}\XS}  \def\pb{{\sbm{p}}\XS}
\def\Qb{{\sbm{Q}}\XS}  \def\qb{{\sbm{q}}\XS}
\def\Rb{{\sbm{R}}\XS}  \def\rb{{\sbm{r}}\XS}
\def\Sb{{\sbm{S}}\XS}  \def\sb{{\sbm{s}}\XS}
\def\Tb{{\sbm{T}}\XS}  \def\tb{{\sbm{t}}\XS}
\def\Ub{{\sbm{U}}\XS}  \def\ub{{\sbm{u}}\XS}
\def\Vb{{\sbm{V}}\XS}  \def\vb{{\sbm{v}}\XS}
\def\Wb{{\sbm{W}}\XS}  \def\wb{{\sbm{w}}\XS}
\def\Xb{{\sbm{X}}\XS}  \def\xb{{\sbm{x}}\XS}
\def\Yb{{\sbm{Y}}\XS}  \def\yb{{\sbm{y}}\XS}
\def\Zb{{\sbm{Z}}\XS}  \def\zb{{\sbm{z}}\XS}


% ALPHABET CURSIF MAJUSCULE, taille adaptative
% Cursif,               % Cursif gras

\def\Ac{{\scu{A}}\XS}   \def\Acb{{\scb{A}}\XS}
\def\Bc{{\scu{B}}\XS}   \def\Bcb{{\scb{B}}\XS}
\def\Cc{{\scu{C}}\XS}   \def\Ccb{{\scb{C}}\XS}
\def\Dc{{\scu{D}}\XS}   \def\Dcb{{\scb{D}}\XS}
\def\Ec{{\scu{E}}\XS}   \def\Ecb{{\scb{E}}\XS}
\def\Fc{{\scu{F}}\XS}   \def\Fcb{{\scb{F}}\XS}
\def\Gc{{\scu{G}}\XS}   \def\Gcb{{\scb{G}}\XS}
\def\Hc{{\scu{H}}\XS}   \def\Hcb{{\scb{H}}\XS}
\def\Ic{{\scu{I}}\XS}   \def\Icb{{\scb{I}}\XS}
\def\Jc{{\scu{J}}\XS}   \def\Jcb{{\scb{J}}\XS}
\def\Kc{{\scu{K}}\XS}   \def\Kcb{{\scb{K}}\XS}
\def\Lc{{\scu{L}}\XS}   \def\Lcb{{\scb{L}}\XS}
\def\Mc{{\scu{M}}\XS}   \def\Mcb{{\scb{M}}\XS}
\def\Nc{{\scu{N}}\XS}   \def\Ncb{{\scb{N}}\XS}
\def\Oc{{\scu{O}}\XS}   \def\Ocb{{\scb{O}}\XS}
\def\Pc{{\scu{P}}\XS}   \def\Pcb{{\scb{P}}\XS}
\def\Qc{{\scu{Q}}\XS}   \def\Qcb{{\scb{Q}}\XS}
\def\Rc{{\scu{R}}\XS}   \def\Rcb{{\scb{R}}\XS}
\def\Sc{{\scu{S}}\XS}   \def\Scb{{\scb{S}}\XS}
\def\Tc{{\scu{T}}\XS}   \def\Tcb{{\scb{T}}\XS}
\def\Uc{{\scu{U}}\XS}   \def\Ucb{{\scb{U}}\XS}
\def\Vc{{\scu{V}}\XS}   \def\Vcb{{\scb{V}}\XS}
\def\Wc{{\scu{W}}\XS}   \def\Wcb{{\scb{W}}\XS}
\def\Xc{{\scu{X}}\XS}   \def\Xcb{{\scb{X}}\XS}
\def\Yc{{\scu{Y}}\XS}   \def\Ycb{{\scb{Y}}\XS}
\def\Zc{{\scu{Z}}\XS}   \def\Zcb{{\scb{Z}}\XS}


% ALPHABET DROIT EN MATHS, taille adaptative

\def\AD{{\sdm{A}}\XS}  \def\aD{{\sdm{a}}\XS}
\def\BD{{\sdm{B}}\XS}  \def\bD{{\sdm{b}}\XS}
\def\CD{{\sdm{C}}\XS}  \def\cD{{\sdm{c}}\XS}
\def\DD{{\sdm{D}}\XS}  \def\dD{{\sdm{d}}\XS}
\def\ED{{\sdm{E}}\XS}  \def\eD{{\sdm{e}}\XS}
\def\FD{{\sdm{F}}\XS}  \def\fD{{\sdm{f}}\XS}
\def\GD{{\sdm{G}}\XS}  \def\gD{{\sdm{g}}\XS}
\def\HD{{\sdm{H}}\XS}  \def\hD{{\sdm{h}}\XS}
\def\ID{{\sdm{I}}\XS}  \def\iD{{\sdm{i}}\XS}
\def\JD{{\sdm{J}}\XS}  \def\jD{{\sdm{j}}\XS}
\def\KD{{\sdm{K}}\XS}  \def\kD{{\sdm{k}}\XS}
\def\LD{{\sdm{L}}\XS}  \def\lD{{\sdm{l}}\XS}
\def\MD{{\sdm{M}}\XS}  \def\mD{{\sdm{m}}\XS}
\def\ND{{\sdm{N}}\XS}  \def\nD{{\sdm{n}}\XS}
\def\OD{{\sdm{O}}\XS}  \def\oD{{\sdm{o}}\XS}
\def\PD{{\sdm{P}}\XS}  \def\pD{{\sdm{p}}\XS}
\def\QD{{\sdm{Q}}\XS}  \def\qD{{\sdm{q}}\XS}
\def\RD{{\sdm{R}}\XS}  \def\rD{{\sdm{r}}\XS}
\def\SD{{\sdm{S}}\XS}  \def\sD{{\sdm{s}}\XS}
\def\TD{{\sdm{T}}\XS}  \def\tD{{\sdm{t}}\XS}
\def\UD{{\sdm{U}}\XS}  \def\uD{{\sdm{u}}\XS}
\def\VD{{\sdm{V}}\XS}  \def\vD{{\sdm{D}}\XS}
\def\WD{{\sdm{W}}\XS}  \def\wD{{\sdm{w}}\XS}
\def\XD{{\sdm{X}}\XS}  \def\xD{{\sdm{x}}\XS}
\def\YD{{\sdm{Y}}\XS}  \def\yD{{\sdm{y}}\XS}
\def\ZD{{\sdm{Z}}\XS}  \def\zD{{\sdm{z}}\XS}


% ALPHABET GRAS DROIT, taille adaptative

\def\Av{{\sbv{A}}\XS}  \def\av{{\sbv{a}}\XS}
\def\Bv{{\sbv{B}}\XS}  \def\bv{{\sbv{b}}\XS}
\def\Cv{{\sbv{C}}\XS}  \def\cv{{\sbv{c}}\XS}
\def\Dv{{\sbv{D}}\XS}  \def\dv{{\sbv{d}}\XS}
\def\Ev{{\sbv{E}}\XS}  \def\ev{{\sbv{e}}\XS}
\def\Fv{{\sbv{F}}\XS}  \def\fv{{\sbv{f}}\XS}
\def\Gv{{\sbv{G}}\XS}  \def\gv{{\sbv{g}}\XS}
\def\Hv{{\sbv{H}}\XS}  \def\hv{{\sbv{h}}\XS}
\def\Iv{{\sbv{I}}\XS}  \def\iv{{\sbv{i}}\XS}
\def\Jv{{\sbv{J}}\XS}  \def\jv{{\sbv{j}}\XS}
\def\Kv{{\sbv{K}}\XS}  \def\kv{{\sbv{k}}\XS}
\def\Lv{{\sbv{L}}\XS}  \def\lv{{\sbv{l}}\XS}
\def\Mv{{\sbv{M}}\XS}  \def\mv{{\sbv{m}}\XS}
\def\Nv{{\sbv{N}}\XS}  \def\nv{{\sbv{n}}\XS}
\def\Ov{{\sbv{O}}\XS}  \def\ov{{\sbv{o}}\XS}
\def\Pv{{\sbv{P}}\XS}  \def\pv{{\sbv{p}}\XS}
\def\Qv{{\sbv{Q}}\XS}  \def\qv{{\sbv{q}}\XS}
\def\Rv{{\sbv{R}}\XS}  \def\rv{{\sbv{r}}\XS}
\def\Sv{{\sbv{S}}\XS}  \def\sv{{\sbv{s}}\XS}
\def\Tv{{\sbv{T}}\XS}  \def\tv{{\sbv{t}}\XS}
\def\Uv{{\sbv{U}}\XS}  \def\uv{{\sbv{u}}\XS}
\def\Vv{{\sbv{V}}\XS}  \def\vv{{\sbv{v}}\XS}
\def\Wv{{\sbv{W}}\XS}  \def\wv{{\sbv{w}}\XS}
\def\Xv{{\sbv{X}}\XS}  \def\xv{{\sbv{x}}\XS}
\def\Yv{{\sbv{Y}}\XS}  \def\yv{{\sbv{y}}\XS}
\def\Zv{{\sbv{Z}}\XS}  \def\zv{{\sbv{z}}\XS}


% Nouveau (08/99) : ALPHABET BLACKBOARD, taille adaptative

\def\Abb{{\sbl{A}}\XS}  \def\abb{{\sbl{a}}\XS}
\def\Bbb{{\sbl{B}}\XS}  \def\bbb{{\sbl{b}}\XS}
\def\Cbb{{\sbl{C}}\XS}  \def\cbb{{\sbl{c}}\XS}
\def\Dbb{{\sbl{D}}\XS}  \def\dbb{{\sbl{d}}\XS}
\def\Ebb{{\sbl{E}}\XS}  \def\ebb{{\sbl{e}}\XS}
\def\Fbb{{\sbl{F}}\XS}  \def\fbb{{\sbl{f}}\XS}
\def\Gbb{{\sbl{G}}\XS}  \def\gbb{{\sbl{g}}\XS}
\def\Hbb{{\sbl{H}}\XS}  \def\hbb{{\sbl{h}}\XS}
\def\Ibb{{\sbl{I}}\XS}  \def\ibb{{\sbl{i}}\XS}
\def\Jbb{{\sbl{J}}\XS}  \def\jbb{{\sbl{j}}\XS}
\def\Kbb{{\sbl{K}}\XS}  \def\kbb{{\sbl{k}}\XS}
\def\Lbb{{\sbl{L}}\XS}  \def\lbb{{\sbl{l}}\XS}
\def\Mbb{{\sbl{M}}\XS}  \def\mbb{{\sbl{m}}\XS}
\def\Nbb{{\sbl{N}}\XS}  \def\nbb{{\sbl{n}}\XS}
\def\Obb{{\sbl{O}}\XS}  \def\obb{{\sbl{o}}\XS}
\def\Pbb{{\sbl{P}}\XS}  \def\pbb{{\sbl{p}}\XS}
\def\Qbb{{\sbl{Q}}\XS}  \def\qbb{{\sbl{q}}\XS}
\def\Rbb{{\sbl{R}}\XS}  \def\rbb{{\sbl{r}}\XS}
\def\Sbb{{\sbl{S}}\XS}  \def\sbb{{\sbl{s}}\XS}
\def\Tbb{{\sbl{T}}\XS}  \def\tbb{{\sbl{t}}\XS}
\def\Ubb{{\sbl{U}}\XS}  \def\ubb{{\sbl{u}}\XS}
\def\Vbb{{\sbl{V}}\XS}  \def\vbb{{\sbl{v}}\XS}
\def\Wbb{{\sbl{W}}\XS}  \def\wbb{{\sbl{w}}\XS}
\def\Xbb{{\sbl{X}}\XS}  \def\xbb{{\sbl{x}}\XS}
\def\Ybb{{\sbl{Y}}\XS}  \def\ybb{{\sbl{y}}\XS}
\def\Zbb{{\sbl{Z}}\XS}  \def\zbb{{\sbl{z}}\XS}


% ALPHABET GREC GRAS, taille adaptative

\def\alphab      {{\sbmm{\alpha}}\XS}
\def\betab       {{\sbmm{\beta}}\XS}
\def\gammab      {{\sbmm{\gamma}}\XS}      \def\Gammab    {{\sbmm{\Gamma}}\XS}
\def\deltab      {{\sbmm{\delta}}\XS}      \def\Deltab    {{\sbmm{\Delta}}\XS}
\def\epsilonb    {{\sbmm{\epsilon}}\XS}
\def\varepsilonb {{\sbmm{\varepsilon}}\XS}
\def\etab        {{\sbmm{\eta}}\XS}
\def\thetab      {{\sbmm{\theta}}\XS}      \def\Thetab    {{\sbmm{\Theta}}\XS}
\def\varthetab   {{\sbmm{\vartheta}}\XS}
\def\iotab       {{\sbmm{\iota}}\XS}
\def\kappab      {{\sbmm{\kappa}}\XS}
\def\lambdab     {{\sbmm{\lambda}}\XS}    \def\Lambdab   {{\sbmm{\Lambda}}\XS}
\def\mub         {{\sbmm{\mu}}\XS}
\def\nub         {{\sbmm{\nu}}\XS}
\def\xib         {{\sbmm{\xi}}\XS}                 \def\Xib        {{\sbmm{\Xi}}\XS}
\def\pib         {{\sbmm{\pi}}\XS}                 \def\Pib        {{\sbmm{\Pi}}\XS}
\def\varpib      {{\sbmm{\varpi}}\XS}
\def\rhob        {{\sbmm{\rho}}\XS}
\def\varrhob     {{\sbmm{\varrho}}\XS}
\def\sigmab      {{\sbmm{\sigma}}\XS}      \def\Sigmab    {{\sbmm{\Sigma}}\XS}
\def\varsigmab   {{\sbmm{\varsigma}}\XS}  \def\Varsigmab {{\sbmm{\Varsigma}}\XS}
\def\phib        {{\sbmm{\phi}}\XS}        \def\Phib       {{\sbmm{\Phi}}\XS}
\def\varphib     {{\sbmm{\varphi}}\XS}
\def\chib        {{\sbmm{\chi}}\XS}
\def\psib        {{\sbmm{\psi}}\XS}        \def\Psib       {{\sbmm{\Psi}}\XS}
\def\omegab      {{\sbmm{\omega}}\XS}      \def\Omegab    {{\sbmm{\Omega}}\XS}
\def\taub        {{\sbmm{\tau}}\XS}
\def\upsilonb    {{\sbmm{\upsilon}}\XS}   \def\Upsilonb  {{\sbmm{\Upsilon}}\XS}
\def\zetab       {{\sbmm{\zeta}}\XS}


% CHIFFRES
% Gras et Blackboard, taille adaptative

\def\zerob   {{\sbv{0}}\XS}        \def\zerobb   {{\sbl{O}}\XS}
\def\unb     {{\sbv{1}}\XS}        \def\unbb     {{\sbl{1}}\XS}
\def\deuxb   {{\sbv{2}}\XS}        \def\deuxbb   {{\sbl{2}}\XS}
\def\troisb  {{\sbv{3}}\XS}        \def\troisbb  {{\sbl{3}}\XS}
\def\quatreb {{\sbv{4}}\XS}        \def\quatrebb {{\sbl{4}}\XS}
\def\cinqb   {{\sbv{5}}\XS}        \def\cinqbb   {{\sbl{5}}\XS}
\def\sixb    {{\sbv{6}}\XS}        \def\sixbb    {{\sbl{6}}\XS}
\def\septb   {{\sbv{7}}\XS}        \def\septbb   {{\sbl{7}}\XS}
\def\huitb   {{\sbv{8}}\XS}        \def\huitbb   {{\sbl{8}}\XS}
\def\neufb   {{\sbv{9}}\XS}        \def\neufbb   {{\sbl{9}}\XS}


% QUELQUES ENSEMBLES

\def\eC{\Cbb}
\def\eE{\Ebb}
\def\eN{\Nbb}
\def\eQ{\Qbb}
\def\eR{\Rbb}
\def\eZ{\Zbb}


snoussi.tex
\documentclass[nomathfonts]{aipproc}
\layoutstyle{6s}

\usepackage{xspace,epic,eepic,graphicx}
\graphicspath{{figures/}}
\def\figures{figures/}

\title[Bayesian source separation]
      {Bayesian source separation with mixture of Gaussians prior 
       for sources and Gaussian prior for mixture coefficients}  
\author{Hichem Snoussi}
  {
  address = {Laboratoire des Signaux et Syst\`emes (L2S), \linebreak   
             Sup\'elec, Plateau de Moulon,  
				 91192 Gif-sur-Yvette Cedex, France}
  ,email   = snoussi@lss.supelec.fr
  }
\author{Ali Mohammad-Djafari}
  {
  address = {Laboratoire des Signaux et Syst\`emes (L2S), \linebreak  
             Sup\'elec, Plateau de Moulon, 
				 91192 Gif-sur-Yvette Cedex, France}
  ,email   = djafari@lss.supelec.fr
  }

%----------- MACROS et Definitions Speciales pour CE document ----------%%
\input{inputs/alphabet}       % Definitions de macros pour les alphabets
\input{inputs/abrege}         % Abbreviations latines, anglaises et fran�aises
\input{inputs/abrmath}        % Abbreviations math�matiques        
\input{inputs/beginend}       % Abbreviations des \begin{} et des \end{}

%----------- MACROS et Definitions Speciales pour CE document ----------%%
\def\d#1{\,{\rm d}#1} % Derivees
\def\cste{\mbox{cste}} % Terme constant
\def\pasr{\mbox{pas}}
\def\plg{polygone }
\def\as{angle solide }
\def\plh{poly\`edre }
\def\DP{d\'eriv\'ee partielle }
\def\DPs{d\'eriv\'ees partielles }
\def\ttr{t\'etra\`edre }
\def\mod{\text{mod}}
\def\ku{1)}
\def\kuu{(\ku}
\def\kd{2-a)}
\def\kdd{(\kd}
\def\kt{2-b)}
\def\ktt{(\kt}
\def\pxta{p\left( \xb(1),\ldots,\xb(T) | \Ab \right)}
\def\pxas{p\left( \xb_{1..T} | \Ab, \sb_{1..T} \right)}
\def\pasx{p\left( \Ab,\sb_{1..T} | \xb_{1..T} \right)}
\def\pax{p\left( \Ab | \xb_{1..T} \right)}
\def\psx{p\left( \sb_{1..T} | \xb_{1..T} \right)}
\def\sbh{\wh{\sb}}
\def\abh{\wh{\ab}}
\def\Sbh{\wh{\Sb}}
\def\Abh{\wh{\Ab}}
\def\Bbh{\wh{\Bb}}
\def\zbh{\wh{\zb}}
\def\xbh{\wh{\xb}}
\def\hsp{\hspace*{.5cm}}
\def\hspu{\hspace*{.1cm}}
\def\hspd{\hspace*{.2cm}}
\def\hspt{\hspace*{.3cm}}
\def\hspq{\hspace*{.4cm}}
\def\hspn{\hspace*{-.6cm}}
\def\hspnd{\hspace*{-1.2cm}}
\def\vspq{\vspace*{.4cm}}
\def\vspd{\vspace*{.2cm}}
\def\vspt{\vspace*{.3cm}}
\def\vspu{\vspace*{.1cm}}
\def\vspuu{\vspace*{.05cm}}
\def\vspn{\vspace*{-.1cm}}
\def\vspnd{\vspace*{-.2cm}}
\def\vspnq{\vspace*{-.4cm}}


\bdoc
\begin{abstract}
In this contribution, we present new algorithms to source separation for the case of  noisy instantaneous linear mixture, within the Bayesian statistical framework. The source distribution prior is modeled by a mixture of Gaussians
\cite{Moulines97} and the mixing matrix elements distributions by a Gaussian  \cite{Djafari99a}. We model the mixture of Gaussians hierarchically by  mean of hidden variables representing the labels of the mixture. Then, we consider the joint a posteriori distribution of sources, mixing matrix elements, labels of the mixture and other parameters of the mixture with appropriate prior probability laws to eliminate degeneracy of the likelihood function of variance parameters and we propose two iterative algorithms to estimate jointly sources,  mixing matrix and hyperparameters: Joint MAP (Maximum \apost) algorithm and penalized EM algorithm. The illustrative example is taken in \cite{Macchi99} to compare with other algorithms proposed in literature.
\keywords{Source separation, Gaussian mixture, classification, JMAP algorithm, Penalized EM algorithm}
\end{abstract}

\maketitle
%\thispagestyle{empty}

\section{Problem description}
We consider a linear instantaneous mixture of $n$ sources. Observations could be corrupted by an additive noise. This noise may represent measurement errors or model incertainty:

\beq\label{model}
\xb(t)=\Ab\sb(t)+\epsilon(t),\qquad t=1,..,T
\eeq

\noindent where $\xb(t)$ is the  ($m \times 1$) measurement vector,  $\sb(t)$ is the ($n \times 1$) source vector which components have to be separated, $\Ab$ is the mixing matrix of dimension ($m \times n$) and  $\epsilon(t)$ represents noise affecting the measurements. We assume that the ($m \times T$) noise matrix $\epsilon(t)$ is statistically independant of sources, centered, white and Gaussian with known variance $\sigma_{\epsilon}^{2}\,\Ib$. We note $\sb_{1..T}$ the matrix $n \times T$ of sources and $\xb_{1..T}$ the matrix $m \times T$ of data.
\\
\\
Source separation problem consists of two sub-problems: Sources restoration and mixing matrix identification. Therefore, three directions can be followed:
\begin{enumerate}
\item {\em Supervised learning}: Identify $\Ab$ knowing a training sequence of sources $\sb$, then use it to reconstruct the sources.
\item {\em Unsupervised learning}: Identify $\Ab$ directly from a part or the whole observations and then use it to recover $\sb$.
\item {\em Unsupervised joint estimation}: Estimate jointly $\sb$ and $\Ab$
\end{enumerate}
In the following, we investigate the third direction. This choice is motivated by  practical cases where sources and mixing matrix are unknown.

% ajouter biblio It's well known that there are two fundamental indeterminations inherent to %source separation \cite{Bermond2000}: Scale and order indeterminations.

This paper is organised as follows: We begin in section II by proposing a Bayesian approach to source separation. We set up the notations, present the prior laws of the sources and the mixing matrix elements and present the joint MAP estimation algorithm assuming known hyperparameters. We introduce, in section III, a hierarchical modelisation of the sources by  mean of  hidden variables representing the labels of the mixture of Gaussians in the prior modeling and present a version of JMAP using the estimation of these hidden variables (classification) as an intermediate step. In both algorithms, we assumed known the hyperparameters which is not realistic in applications. That is why, in section IV, we present an original method for the estimation of  hyperparameters which takes advantages of using this hierarchical modeling. Finally, since EM algorithm has been used extensively in source separation \cite{Bermond00}, we considered this  algorithm and  propose, in section V, a penalized version of the EM algorithm for source separation. This penalization of the likelihood function is necessary to eliminate its degeneracy when some variances of Gaussian mixture approche zero \cite{Ridolfi99}. Each section  is supported by one typical simulation result and partial conclusion. At the end, we compare the two last algorithms.

\section{Bayesian approach to source separation}
Given the observations $\xb_{1..T}$, the joint \apost distribution of unknown variables $\sb_{1..T}$ and $\Ab$ is:
\beq\label{distribution}
\pasx \propto \pxas\,p(\Ab)\,p(\sb_{1..T})
\eeq
\noindent where $p(\Ab)$ and $p(\sb_{1..T})$ are the prior distributions through which we modelise our \aprio information about sources $\sb$ and mixing matrix $\Ab$. $\pxas$ is the joint likelihood distribution. We have, now, three directions:
\begin{enumerate}
\item First, integrate (\ref{distribution}) with respect to $\sb_{1..T}$ to obtain the marginal in $\Ab$ and then estimate  $\Ab$ by:
\beq \label{JMAP2}
\Abh=\argmax_{\Ab}{\{J(\Ab)=\ln \,\pax \}}
\eeq
\item Second, integrate (\ref{distribution})  with respect to $\Ab$ to obtain the marginal in $\sb_{1..T}$ and then estimate $\sb_{1..T}$ by:
\beq \label{JMAP3}
\sbh_{1..T}=\argmax_{\sb_{1..T}}{\{J(\sb_{1..T})=\ln \,\psx \}}
\eeq
\item Third, estimate jointly $\sb_{1..T}$ and  $\Ab$:
\beq \label{JMAP1}
(\Abh,\sbh_{1..T})=\argmax_{(\Ab,\sb_{1..T})}{\{J(\Ab,\sb_{1..T})=\ln \,\pasx\}}
\eeq
\end{enumerate}   

\subsection{Choice of \aprio distributions}
The \aprio distribution  reflects our knowledge concerning the parameter to be estimated. Therefore, it must be neither very specific to a particular problem nor too general (uniform) and non informative. A parametric model for these distributions seems to fit this goal: Its stucture expresses the particularity of the problem and  its parameters allow a certain flexibility.

\noindent{\bf{Sources \aprio:}} For sources $\sb$, we choose a mixture of Gaussians \cite{Moulines97}: 
\beq
p(s_{j})=\sum_{i=1}^{q_{j}}\alpha_{ji}\Nc(m_{ji}, \sigma_{ji}^{2}),\qquad j=1..n
\eeq
Hyperparameters $q_{j}$ are supposed to be known. 

\noindent This choice was motivated by the following points:
\bit
\item It represents a  general class of distributions and is convenient in many digital communications and image processing applications.
\item For a Gaussian likelihood $p\left( \xb_{1..T} | \sb_{1..T}, \Ab \right)$ (considered as a function of $\sb_{1..T}$), the \apost law remains in the same class (conjugate prior). We then have only to update the parameters of the mixture with the data.
\eit
\noindent {\bf{Mixing matrix \aprio:}} To account for some model uncertainty, we assign a Gaussian prior law to each element of the mixing matrix $\Ab$: 
\beq
p(\Ab_{ij})=\Nc(\Mb_{ji}, \sigma_{a,ij}^{2})
\eeq
which can  be interpreted as knowing  every element ($\Mb_{ji}$) with some uncertainty ($\sigma_{a,ij}^{2}$). We underline here the advantage of estimating  the mixing matrix $\Ab$ and not a separating matrix $\Bb$ (inverse of $\Ab$) which is the case of almost all the existing methods for source separation (see for example \cite{Cardoso96}). This approach has at least two advantages: ({\bf{i}}) $\Ab$ does not need to be invertible ($n \neq m$), ({\bf{ii}}) naturally, we have some \aprio information on the mixing matrix not on its inverse which may not exist. 
\vspd

\subsection{JMAP algorithm}
We propose an alternating iterative algorithm to estimate jointly $\sb_{1..T}$ and $\Ab$ by extremizing the log-posterior distribution:
\beq 
\left\{\barr{ccc}
\sbh_{1..T}^{(k)} &=& \argmax_{\sb_{1..T}}{ \ln\,p\left( \Abh^{(k-1)},\sb_{1..T} | \xb_{1..T} \right) }
\\ 
 \Abh^{(k)} &=& \argmax_{\Ab}{ \ln\,p\left( \Ab,\sbh_{1..T}^{(k)} | \xb_{1..T} \right) }
\earr\right.
\eeq
In the following, we suppose that sources are white and spatially independant. This assumption is not necessary in our approach but we start from here to be able to compare later with  other classical methods in which this hypothesis is fundamental.

\noindent With this hypothesis, in step $(k+1)$, the criterion to optimize with respect to $\sb_{1..T}$ is:
\beq
J(\sb_{1..T})=\sum_{t=1}^{T} \left[ \ln\,p\left( \xb(t) | \Abh^{(k)},\sb(t) \right) +\sum_{j=1}^{n}\ln\,p_j\left( \sb_j(t) \right) \right]
\eeq
Therefore, the optimisation is done independantly at each time $t$:
\beq
\sbh(t)^{(k+1)} = \argmax_{\sb(t)}{ \{ \ln\,p\left( \xb(t) | \Abh^{(k)},\sb \right) + \sum_{j=1}^{n} \ln\,p_{j} \left( \sb_{j}(t) \right) \}}
\eeq
The \apost distribution of $\sb$ is a mixture of $\prod_{j=1}^{n}q_j$ Gaussians. This leads to a high computational cost. To obtain a more reasonable algorithm, we propose an iterative scalar algorithm. The first step consists in estimating each source component knowing the other components estimated in the previous iteration:
\beq
\sbh_{j}(t)^{(k+1)} =\argmax_{\sb_{j}(t)}{ \{\ln\,p\left( \sb_{j}(t) | \xb(t),\Abh^{(k)},\sbh_{l \neq  j}(t)^{(k)} \right) \}}    
\eeq  
The \apost distribution of $\sb_j$ is a mixture of $q_j$ Gaussians: $\sum_{z=1}^{q_{j}}\alpha_{jz}^{'}\Nc(m_{jz}^{'}, {\sigma_{jz}^{'}}^{2})$, with:
\beq\left\{\barr{c}\label{zapost}
\displaystyle{m_{jz}^{'}=\frac{\sigma_j^{2}m_{jz}+\sigma_{jz}^{2}m_{j}}{\sigma_j^{2}+\sigma_{jz}^{2}}}\\~\\
\displaystyle{{\sigma_{jz}^{'}}^{2}=\frac{\sigma_j^{2}\,\sigma_{jz}^{2}}{\sigma_j^{2}+\sigma_{jz}^{2}}}\\~\\
\displaystyle{\alpha_{jz}^{'}= \alpha_{jz}\,\sqrt{\frac{1}{\sigma_{jz}^{2}+\sigma_j^{2}}}\Exp{\frac{-1}{2}\frac{1}{\sigma_{jz}^{2}+\sigma_j^{2}}(m_j-m_{jz})^2}}
\earr\right.
\eeq
where
\beq\label{treize}
\left\{\barr{c}
\displaystyle{\sigma_j^2=\frac{\sigma_{\epsilon}^2}{\sum_{i=1}^{m}\Ab_{ij}^2}}
\\~\\
\displaystyle{m_j = \frac{\sum_{i=1}^n \Ab_{ij}\,(\xb_i-\xbh_i)}{\sum_{i=1}^{m}\Ab_{ij}^2}}
\\~\\
\displaystyle{\xbh_i =\sum_{l \neq j}\Ab_{il}\,\sb_l}
\earr\right.
\eeq
\\

If the means $m_{jz}^{'}$ aren't close to each other, we are in the case of a multi-modal distribution. The algorithm to estimate $\sb_j$ is to first compute $\xbh_i$, $i=1, \dots, m$, $m_j$ and $\sigma_j^2$ by (\ref{treize}) and then $\alpha_{jz}^{'}$, ${\sigma_{jz}^{'}}^{2}$ and $m_{jz}^{'}$ by (\ref{zapost}), and select the $m_{jz}^{'}$ for which the ratio $\frac{\alpha_{jz}^{'}}{\sigma_{jz}^{'}}$ is the greatest one.

After a full update of all sources $\sb_{1..T}$, the estimate of $\Ab$ is obtained by optimizing:
\beq 
\barr{lcl} 
J(\Ab)=
\sum_{t=1}^{T}\ln\,p\left( \xb(t) | \Ab,\sbh^{k+1}(t) \right) +\ln\,p\left( \Ab(t) \right) + cte
\earr
\eeq 
which is quadratic in elements of $\Ab$. The  gradient has then a simple expression:
\beq
\frac{\partial{J(\Ab)}}{\partial{\Ab_{i,j}}}=\sum_{t=1}^{T}\frac{1}{\sigma_{\epsilon}^2}\sbh_{j}^{k+1}(t) \left( \xb_{i}(t)-\left[ \Ab\,\sbh^{k+1}(t) \right]_{i}\right)-\frac{1}{\sigma_{a;i,j}^2} \left( \Ab_{i,j}-\Mb_{i,j} \right)
\eeq
Cancelling the gradient to zero and defining $\Lambdab_{i,j}=\frac{\sigma_{\epsilon}^{2}}{\sigma_{a;i,j}^2}$, we obtain the following relation:
\beq\label{matrice}
\left[\sum_{t=1}^{T} \left( \xb(t)-\Ab \sbh^{k+1}(t) \right) \sbh^{k+1}(t)^{T}\right]_{i,j}-\Lambdab_{i,j}\left( \Ab_{i,j}-\Mb_{i,j} \right)=0
\eeq
We define the operator {\bf{Vect}} transforming a matrix to a vector by the concatenation of the transposed rows. Operator {\bf{Mat}} is the inverse of {\bf{Vect}}. Applying operator {\bf{Vect}} to relation (\ref{matrice}), we obtain the following expression:
\beq
Vect\left( \xb_{1..T}(\sbh_{1..T}^{k+1})^{T}\right)+\mu Vect(M)=\left( \mu+S^{*} \right) Vect\Ab
\eeq
\noindent where $\mu$ is a diagonal matrix $(nm \times nm)$ which diagonal vector is $Vect((\Lambda_{i,j})_{i=1..m,j=1..n})$ and $\Sb^{*}$ the matrix $(nm \times nm)$ with block diagonals  $\sbh_{1..T}\sbh_{1..T}^{T}$ estimated at  iteration $(k+1)$. We have finally the explicit estimation of $\Ab$:
\beq
\Abh^{k+1}=Mat\left(\left[ \mu+ S^{*} \right]^{-1}\left[ \mu Vect(\Mb)+Vect\left( \xb_{1..T}(\sbh_{1..T}^{k+1})^{T}\right) \right] \right)
\eeq

To show the faisability of this algorithm, we consider in the following a telecommunication example. For this, we simulated synthetic data with sources described by a mixture of $4$ Gaussians centered at $-3$, $-1$, $1$ and $3$, with the same variance $0.01$ and weighted by 0.3, 0.1, 0.4 and 0.2. The unknown mixing matrix is $\Ab=\left( \barr{cc}1&-0.6\\0.6&1 \earr \right)$. We fixed the \aprio parameters of $\Ab$ to: $\Mb=\left( \barr{cc}1&0\\0&1 \earr \right)$ and $\Lambdab=\left( \barr{cc}150&0.009\\0.009&150 \earr \right)$, meaning that we are nearly sure of diagonal values but we are very uncertain about the other elements of $\Ab$. Noise of variance  $\sigma_{\epsilon}^2=1$ was added to the data. The  figure $1$ illustrates the ability of the algorithm to perform the separation. However, we note that estimated sources are very centered arround the means. This is because  we fixed very low values for the \aprio variances of Gaussian mixture. Thus, the algorithm is sensitive to the \aprio parameters and  exploitation of data is useful. We will see in section IV how to deal with this issue.

\bigskip
\noindent\btabu{@{}c@{}}
\includegraphics[width=\textwidth,height=45mm]{figh1b}\\
\hspace*{0.5cm}(a)\hspace*{4.5cm}(b)\hspace*{4.5cm}(c)\\~\\
Figure 1- Results of separation with QAM-$16$ (Quadratic Amplitude Modulation) \\
using JMAP algorithm: (a) phase space distribution of sources,\\
(b) mixed signals, and (c) separated sources
\etabu
\\
\\
Now, we are going to re-examine closely the expression for the  \apost distribution of sources. It's a multi-modal distribution if the Gaussian means aren't too close. The maximum of this distribution doesn't correspond, in general, to the maximum of the most probable Gaussian. So, we intend to estimate first, at each time $t$, the \aprio Gaussian law according to which the source $\sb(t)$ is generated (classification) and then estimate $\sb(t)$ as the mean of the \apost Gaussian.  This leads us to the introduction of hidden variables and hierarchical modelization. 

\section{Hidden variables}
The \aprio distribution of the component $s_j$ is  $p(s_{j})=\sum_{i=1}^{q_{j}}\alpha_{ji}\Nc(m_{ji}, \sigma_{ji}^{2})$. We consider now the hidden variable $z_j$ taking its values in the discrete set $\Zc_j=\left( 1, \dots, q_j\right)$  so each source can belong to one of the $q_j$ sources, with $\alpha_{ji}=p\left( z_j=i \right)$. Given $z_j=i$, $\sb_j$  is normal $\Nc(m_{ji},\sigma_{ji}^{2})$. We can extend this notion to vectorial case by considering the vector $\zb=[z_1,\dots,z_{n}]$ taking its values in the set $\Zc=\Pi_{j=1}^{n}\Zc_{j}$. The $\sb$ distribution given $\zb$ is a normal law $p(\sb | \zb)=\Nc(\mb_z,\bf{\Gamma}_z)$ with:
\beq 
\mb_z=[m_{1z_1}, m_{2z_2},\dots, m_{nz_n}]
\eeq
\beq
\bf{\Gamma}_{z}=diag(\sigma_{1z_1}^{2}, \sigma_{2z_2}^{2},\dots, \sigma_{nz_n}^{2})
\eeq
The marginal \aprio law of $\sb$ is the mixture of $\Pi_{j=1}^{n}q_j$ Gaussians:
\beq
p(\sb)=\sum_{\zb \in \Zc} p(\zb)p(\sb | \zb)
\eeq
We can re-interpret this mixture by considering it as a discrete set of couples $\left( \Nc_z, p(\zb) \right)$ (see Figure $2$). Sources which belong to this class of distributions are generated as follows: First, generate the hidden variable $\zb \in \Zc$ according $p(\zb)$ and then, given this $\zb$, generate  $\sb$ according to $\Nc_z$. This model can be extended to include continuous values of $\zb$ (also continuous distribution $f(\zb)$) and then to take account of infinity of distributions in only one class  (see Figure $2$).
\\
\bcc 
\fbox{\input{\figures cache2.eepic}}\\~\\
Figure 2- Hierarchical modelization with hidden variables
\ecc

\subsection{\apost distribution of sources}
In the following, we suppose that mixing matrix is known. The joint law of $\sb$, $\zb$ and $\xb$ can be  factorized in two forms: $p(\sb,\zb,\xb)=p(\xb|\sb)p(\sb|\zb)p(\zb)$ or $p(\sb,\zb,\xb)=p(\sb | \xb,\zb)p(\zb|\xb)p(\xb)$. Thus, the marginal \apost law has two forms:
\beq
p(\sb | \xb)=\sum_{\zb \in \Zc}\frac{p(\zb)\,p(\xb|\sb)\,p(\sb|\zb)}{p(\xb)}
\eeq
or
\beq
p(\sb | \xb)=\sum_{\zb \in \Zc}p(\zb|\xb)\,p(\sb | \xb,\zb)
\eeq
We note in the second form that the \apost is in the same class that of the \aprio (same expressions but conditionally to $\xb$). This is due to the fact that mixture of Gaussians is a conjugate prior for Gaussian likelihood. Our strategy of estimation is based on this remark: The  sources are modeled hierarchically, we estimate them hierarchically; we begin by estimating the hidden variable using $p(\zb|\xb)$  and then estimate sources using $p(\sb | \xb,\zb)$ which is Gaussian of mean $\thetab_{xz}$:
\beq\label{moyenne}
\thetab_{xz}=\mb_{z}+\bf{\Gamma}_{z}\Ab^{t}R_z(\xb-\Ab \mb_z)
\eeq
and variance $\Vb_{xz}$:
\beq\label{variance}
\Vb_{xz}=\bf{\Gamma}_z-\bf{\Gamma}_z\Ab^{t}R_z\Ab \bf{\Gamma}_z
\eeq
where,
\beq
\Rb_{z}=(\Ab\bf{\Gamma}_z\Ab^{t}+ \Rb_n)^{-1}
\eeq
and $\Rb_n$ represent the noise covariance.

Now we have to estimate $\zb$ by using $p(\zb | \xb)$ which is  obtained by integrating the joint \apost of $\zb$ and $\sb$ with respect to $\sb$:
\beq
p(\zb | \xb) = \int p(\zb , \sb | \xb) d\sb \,\propto\, p(\zb)\int p(\xb|\sb)\,p(\sb | \zb) d\sb
\eeq


The expression to integrate is Gaussian in $\sb $. The result is immediate:

\beq
p(\zb | \xb) \propto p(\zb)\mid \Gammab_z \mid ^{-\frac{1}{2}}\mid \Vb_{xz} \mid ^{\frac{1}{2}} \Exp{K_{zx}}
\eeq
where:
\beq
\left\{\barr{ccc}
\Kb_{zx} &=& -\frac{1}{2}(\Ab \mb_z-\xb)^{t}\Qb_{xz}(\Ab \mb_z-\xb)
\\
\Qb_{xz} &=& (\Ib-\Rb_z\Ab \bf{\Gamma}_z \Ab^{t})\Rb_{n}^{-1}(\Ib-\Ab \bf{\Gamma}_z\Ab^{t}R_z)+R_z\Ab\bf{\Gamma}_z\Ab^{t}R_z
\earr\right\}
\eeq

\noindent If now we consider the whole observations, the law of $\zb_{1..T}$ is:
\beq
p(\zb_{1..T} | \xb_{1..T})\propto p(\zb_{1..T}) \int p(\xb_{1..T}|\sb_{1..T})\, p(\sb_{1..T}|\zb_{1..T}) \,d\sb_{1..T}
\eeq
Supposing that $\zb(t)$ are \aprio independant, the last relation becomes:
\beq
p(\zb_{1..T} | \xb_{1..T})\propto \Pi_{t=1}^{T}\left\{ p(\zb(t)) \int p(\xb(t)|\sb(t)) \,p(\sb(t)|\zb(t))\, d\sb(t) \right\}
\eeq  
Estimation of $\zb_{1..T}$ is then performed observation by observation:
\beq
\argmax_{\zb_{1..T}}{p(\zb_{1..T} | \xb_{1..T})}=\left( \argmax_{\zb(t)}{p(\zb(t) | \xb(t))} \right)_{t=1..T}
\eeq

\subsection{Hierarchical JMAP algorithm}
Taking into account of this hierarchical model, the JMAP algorithm is implemented in three steps. At iteration $(k)$:
\begin{enumerate}
\item First, estimate the hidden variable $\zbh_{MAP}$ (combinatary estimation) given observations and mixing matrix estimated in the previous iteration:
\beq
\zbh_{MAP}^{(k)}(t)=\argmax_{\zb(t)}{\{p\left(\zb(t)|\xb(t), \Abh^{(k-1)}\right)\}}
\eeq
\item Second, given the estimated $\zbh_{MAP}^{(k)}$, source vector $\sb$ follows Gaussian law \\
$\Nc(\thetab_{x\zbh_{MAP}^{(k)}},\Vb_{x\zbh_{MAP}^{(k)}})$ and then the source estimate is $\thetab_{x\zbh_{MAP}^{(k)}}$.
\item Third, given the estimated sources $\sbh^{k}$, mixing matrix is evaluated as in the algorithm of section II.
\end{enumerate}
We evaluated this algorithm using the same  synthetic data as in  section $2$. Separation was robust as shown in Figure $3$:  

\bigskip
\noindent\btabu{@{}c@{}}
\includegraphics[width=\textwidth,height=45mm]{fig3}\\
\hspace*{0.5cm}(a)\hspace*{4.5cm}(b)\hspace*{4.5cm}(c)\\~\\
Figure 3- Results of separation with QAM-$16$ \\
using Hierarchical JMAP algorithm: (a) phase space distribution of sources,\\
(b) mixed signals, and (c) separated sources
\etabu 
\\
\\  
The Bayesian approach allows us to express our \aprio information via parametric prior models. However, in general, we may not know the parameters of the \aprio distributions. This is the task of the next section where we estimate the unknown hyperparameters always in a Bayesian framework.   

\section{Hyperparameters estimation}
The hyperparameters considered here are the means and the variances of Gaussian mixture prior of sources: $
\sb_j \sim \sum_{z=1}^{q_j}\Pi_{jz} \Nc \left( m_{jz}, \frac{1}{\psi_{jz}} \right)$, $j=1,\dots, n$. 
We develop, in the following, a novel method to extract the hyperparameters from the observations $\xb_{1..T}$. The main idea is:  
conditioned on the hidden variables $(z_{j})_{1..T}=[z_j(1),\dots, z_j(T)]$,  hyperparameters $m_{jz}$ and $\psi_{jz}$ for $z \in \Zc_j=(1,\dots,q_{j})$ are means and variances of a Gaussian distribution. Thus, given the vector $(z_{j})_{1..T}=[z_j(1),\dots, z_j(T)]$, we can perform a partition of the set $\Tc=[1,\dots, T]$ into sub-sets $\Tc_z$ as:
\beq
\Tc_z=\left \{\, t \,| \,z_j(t) = z \right \} \, ,\,  \, z \in \Zc_j
\eeq
This is the classification step.

Suppose now  that mixing matrix $\Ab$ and components $\sb_{l \neq j}$ are fixed and we are interested in the estimation of $m_{jz}$ and $\psi_{jz}$. Let $\theta_{jz}=( m_{jz}\,,\, \psi_{jz} )$.\\

The joint \apost law of $\sb_j$ and $\theta_{jz}$  given $z_j$ at time  $t$ is:
\beq
p( \sb_j,\,\theta_{jz} \,|\, \xb,\,z_j) \propto p(\xb \,|\, \sb)\,p(\sb_j\,|\, \theta_{jz},\,z_j)\,p(\theta_{jz} \,|\,z_j)
\eeq
$p(\sb_j\,|\, \theta_{jz},\,z_j)$ is Gaussian of mean $m_{jz}$ and inverted variance $\psi_{jz}$.\\
$p(\theta_{jz} \,|\,z_j)\,=\,p(\theta_{jz})\,=\,p(m_{jz})\,p(\psi_{jz})$ is hyperparameters \aprio.
The marginal \apost distribution of $\theta_{jz}$ is obtained from previous relation by integration over $\sb_j$:
\beq
p(\theta_{jz} \,|\, \xb,\,z_j) \propto p(\theta_{jz})\int_{\sb_j}p(\xb \,|\, \sb)\,p(\sb_j\,|\, \theta_{jz},\,z_j)\,d\sb_j.
\eeq
The expression inside the integral is proportional to the joint \apost distribution of $(\sb_j\,,\,z_j)$ given $\xb$ and $\theta_{jz}$, thus:
\beq
p(\theta_{jz} \,|\, \xb,\,z_j) \propto p(\theta_{jz})\,p(z_j\,|\,\xb,\,\theta_{jz}) .
\eeq
where $p(z_j\,|\,\xb,\,\theta_{jz})$ is proportional to $\alpha_{jz}^{'}$ as defined in  expression (\ref{zapost}). Noting $\phi_j=1\,/\,\sigma_j^{2}$ and $\psi_{jz}=1\,/\,\sigma_{jz}^{2}$, we have:
\beq\label{vraisembance}
p(\theta_{jz} \,|\, \xb,\,z_j) \propto p(\theta_{jz})\sqrt{\frac{\phi_j\,\psi_{jz}}{\phi_j\,+\,\psi_{jz} }}\,\Exp{-\frac{1}{2}\,\frac{\phi_j\,\psi_{jz}}{\phi_j+\psi_{jz} }\,(m_{jz}-m_j)^2}
\eeq
Note that the likelihood is normal for means $m_{jz}$ and Gamma  for $\lambda_{jz}=\left(\phi_j\psi_{jz}\right)/\left(\phi_j+
\psi_{jz}\right) $. 

\noindent Choosing a uniform \aprio for the means, the estimate of $m_{jz}$ is:
\beq\label{mean}
\wh{m}_{jz}^{MAP} = \frac{\sum_{t \in \Tc_z}m_j(t)}{T_z} 
\eeq
For variances, we can choose ({\bf i}) an inverted Gamma prior $\Gc\left( \alpha, \beta \right)$ after developing the expression for $\lambda_{jz}$ knowing the relative order of $\psi_{jz}$ and $\phi_j$ (to make $\lambda_{jz}$ linear in $\psi_{jz}$) or ({\bf ii}) an a prior which is  Gamma in $\lambda_{jz}$. These choices are motivated by two points: First, it is a proper prior which eliminate degenaracy of some variances at zero (It is shown in \cite{Ridolfi99} that hyperparameter likelihood (noiseless case without mixing) is unbounded causing a variance degeneracy at zero). Second, it is a conjugate prior so estimation expressions remain simple to implement. The estimate of inverted variance (first choice when $\psi_{jz}$ is the same order of $\phi_j$) is:
\beq\label{variance2}
\wh{\psi}_{jz}^{MAP} = \frac{\alpha^{posteriori}-1}{\beta^{posteriori}}
\eeq
with $\alpha^{posteriori}=\alpha+\frac{\Tc_z}{2}$  and $\beta^{posteriori}=\beta+\frac{\sum_{t \in \Tc_z}(m_j(t)-\wh{m}_{jz}^{MAP})^2}{4}$. 
\vspd
\vspq

\subsection{Hierarchical JMAP  including estimation of hyperparameters}
Including the estimation of hyperparameters, the proposed hierarchical JMAP algorithm is composed of five steps:
\begin{enumerate}
\item Estimate hidden variables $(\wh{z}_j)_{1..T}^{MAP}$ by:
\beq
(\wh{z}_j)_{1..T}^{MAP} = 
(\argmax_{z_j}{p(z_j\,|\,\xb(t),\,m_{jz}\,,\,\psi_{jz},\,\Ab,\,\sb_{l \neq j})})_{1..T}   
\eeq
which permits to estimate  partitions:
\beq
\wh{\Tc_z} = \left \{ t \, |\, (\wh{z}_j)^{MAP}(t)=z \right\}
\eeq
This corresponds to the classification step in the previous algorithm
\item Given the estimate of partitions, hyperparameters $\wh{\psi}_{jz}^{MAP}$ and $\wh{m}_{jz}^{MAP}$ are updated according to equations (\ref{mean}) and (\ref{variance}).
The following steps are the same as those in the previous proposed algorithm
\item Re-estimation of hidden variables $(\wh{z}_j)_{1..T}^{MAP}$ given the estimated hyperparameters.
\item Estimation of sources  $(\wh{\sb})_{1..T}^{MAP}$.
\item Estimation of mixing matrix $(\wh{\Ab})^{MAP}$.
\end{enumerate}
\vspq

\subsection{Simulation results}
To be able to compare the results obtained by this algorithm and the Penalized likelihood algorithm developed in the next section with the results obtained by some other classical methods, we generated data according to the example described in \cite{Macchi99}.
 
\noindent{\bf{Data generation}}: $2$-D sources, every component \aprio is mixture of two Gaussians ($\pm 1$), $\psi=100$ for all Gaussians. Original sources are mixed with mixing matrix $\Ab=\left( \barr{cc} 1&-0.6\\0.4&1 \earr \right)$. A noise of variance $\sigma_{\epsilon}^2 = 0.03$ is added ($SNR=15\,dB$). Number of observations is $1000$.

\noindent {\bf{Parameters}}: $\Mb=\left( \barr{cc} 1&0\\0&1 \earr \right)$,  $\Lambdab=\left( \barr{cc} 150&0.009\\0.009&150 \earr \right)$, $\Pi=\left( \barr{cc} 0.5&0.5\\0.5&0.5 \earr \right)$, $\alpha=200$ and $\beta=2$.  

\noindent {\bf{Initial conditions}}: $\Ab^{(0)}=\left(\barr{cc}1 &0 \\0 &1\earr\right)$, $\psi^{(0)}=\left(\barr{cc}1 &1 \\1 &1\earr\right)$, $m^{(0)}=\left(\barr{cc}0 &0 \\0 &0\earr\right)$ and $\sb^{(0)}$ generated according to  $\sb_j^{(0)} \sim \sum_{z=1}^{q_j}\Pi_{jz} \Nc(m_{jz}^{(0)}, \frac{1}{\psi_{jz}^{(0)}})$.

Sources are recovered with negligible mean quadratic error: 
$MEQ(\sb_1)=0.0094$ and $MEQ(\sb_2)=0.0097$. The following figures illustrate separation results:

The non-negative performance index of \cite{Moreau96} is used to chacarterize mixing matrix identification achievement:
\[
ind(S=\Abh^{-1}\,\Ab) = \frac{1}{2} \left[ \sum_i \left( \sum_j\frac{|S_{ij}|^2}{max_{l}|S_{il}|^2}-1 \right)+\sum_j \left( \sum_i\frac{|S_{ij}|^2}{max_{l}|S_{lj}|^2}-1 \right) \right]
\]
Figure $7a$ represents the index evolution through iterations. Note the convergence of JMAP algorithm since iteration $30$ to a satisfactory value of $-45\,dB$. For the same SNR, algorithms PWS, NS \cite{Macchi99} and EASI \cite{Cardoso96} reach a value greater than $-35\,dB$ after $6000$ observations.
Figures $7b$ and $7c$  illustrate the identification of hyperparameters. We note the algorithm convergence to the original values ($-1$ for $m_{11}$ and $100$ for $\psi_{11}$).
In order to validate the idea of data classification before estimating hyperparameters, we can visualize the evolution of classification error (number of data badly classified). Figure $7d$ shows that this error converges to zero at iteration $15$. Then, after this iteration, hyperparameters identification is performed on the true classified data. Estimation of $m_{jz}$ and $\psi_{jz}$ takes into account only data which belong to this class and then it is not corrupted by other data which bring erroneous information on these hyperparameters. 

\bigskip
\bcc
\btabu[b]{@{}cc@{}}
\btabu[b]{@{}c@{}}
$\left\{\barr{@{}l@{}}
s_1(t) \\ ~\\~\\   
s_2(t)  
\earr\right.$  
~\\ ~\\ 
$\left\{\barr{@{}l@{}}
x_1(t) \\ ~\\~\\   
x_2(t)  
\earr\right.$ 
~\\ ~\\
$\left\{\barr{@{}l@{}}
\wh{s}_1(t) \\ ~\\ ~\\   
\wh{s}_2(t)     
\earr\right.$ 
~\\ ~\\
$\left\{\barr{@{}l@{}}
\wh{s}_1(t)-s_1(t) \\ ~\\ ~\\ 
\wh{s}_2(t)-s_2(t)     
\earr\right.$ 
\\ ~\\ 
\etabu
&
\btabu[b]{@{}c@{}}
\includegraphics[width=100mm,height=100mm]{fig1zh15sd}\\
\etabu
\etabu\\~\\
Figure 4- Separation results with $SNR=15\,dB$
\ecc

\bigskip
\noindent\btabu{@{}c@{}}
\includegraphics[width=\textwidth,height=45mm]{fig3zh15} \\
Figure 5- Separation results with $SNR=15\,dB$: Phase space distribution of sources,\\
mixed signals and separated sources.
\etabu

\bigskip
\noindent\btabu{@{}c@{}}
\includegraphics[width=\textwidth,height=60mm]{fig2zh15}\\~\\
Figure 6- Separation results with $SNR=15\,dB$: Histograms of sources,\\
mixed signals and separated sources.
\etabu

\bigskip
\noindent\btabu{@{}c@{}c@{}}
\includegraphics[width=\textwidth/2,height=50mm]{fig6zh15}&
\includegraphics[width=\textwidth/2,height=50mm]{fig4zh15}
\\~\\
{\small Figure 7-a- Evolution of index through iterations}& 
{\small Figure 7-b- Identification of $m_{11}$}
\\~\\
\includegraphics[width=\textwidth/2,height=50mm]{fig5zh15}&
\includegraphics[width=\textwidth/2,height=50mm]{fig7zh15}\\~\\
{\small Figure 7-c- Identification of $\psi_{11}$}& 
{\small Figure 7-d- Evolution of classification error}
\etabu

\bigskip
Thus, a joint estimation of sources, mixing matrix and hyperparameters is performed successfully with a JMAP algorithm. The EM algorithm was used in \cite{Bermond00} to solve source separation problem in a maximum likelihood context. We  now  use the EM algorithm in a Bayesian approach  to take into account of our \aprio information on the mixing matrix.

\section{\textit{Penalized EM}} 
The EM algorithm has been used extensively in data analysis to find the maximum likelihood estimation of a set of parameters from given data \cite{Redner84}. Considering both the mixing matrix $\Ab$ and hyperparameters $\thetab$, at the same level, being unknown parameters and complete data $\xb_{1..T}$ and $\sb_{1..T}$. Complete data means jointly observed data $\xb_{1..T}$ and unobserved data $\sb_{1..T}$. The EM algorithm is executed in two steps: ({\bf{i}}) E-step (expectation) consists in forming the logarithm of the joint distribution of observed data $\xb$ and hidden data $\sb$ conditionally to parameters $\Ab$ and $\thetab$ and then compute its expectation conditionally to $\xb$ and estimated parameters  $\Ab^{'}$ and $\thetab^{'}$ (evaluated in the previous iteration), ({\bf{ii}}) M-step (maximization) consists of the maximization of the obtained functional with respect to the parameters $\Ab$ and $\thetab$:
\begin{enumerate}
\item E-step : 
\beq
Q \left( \Ab,\,\thetab\,|\,\Ab',\,\thetab' \right) = E_{\xb,\sb}\left[  \log\,p(\xb,\,\sb\,|\,\Ab,\,\thetab) \,|\,\xb,\,\Ab',\,\thetab' \right]
\eeq
\item M-step : 
\beq
\left( \wh{\Ab},\,\wh{\thetab} \right) = \argmax_{\left( \Ab,\,\thetab\right)}{\{Q \left( \Ab,\,\thetab\,|\,\Ab',\,\thetab' \right)\}}
\eeq
\end{enumerate}
Recently, in \cite{Bermond00}, an EM algorithm has been used in source separation with mixture of Gaussians as sources prior. In this work, we show that:
\begin{enumerate}
\item This algorithm fails in estimating variances of Gaussian mixture. We proved that this is because the degeneracy of the estimated variance to zero. 
\item The computational cost of this algorithm is very high.
\item The algorithm is very sensitive to initial conditions.
\item In \cite{Bermond00}, there's neither an \aprio distribution on the mixing matrix $\Ab$ or on the hyperparameters $\thetab$.
\end{enumerate}  

Here, we propose to extend this algorithm in two ways by:
\begin{enumerate}
\item Introducing an \aprio distribution for $\thetab$ to eliminate degeneracy and an \aprio distribution for $\Ab$ to express our previous knowledge on the mixing matrix. 
\item Taking advantage of our hierarchical model and the idea of classification to reduce the computational cost.
\end{enumerate}
To distinguish the proposed algorithm from the one proposed in \cite{Bermond00}, we call this algorithm the \textit{Penalized EM}. The two steps become:

\begin{enumerate}
\item E-step : 
\beq 
Q \left( \Ab,\,\thetab\,|\,\Ab',\,\thetab' \right) = E_{\xb,\sb}\left[  \log\,p(\xb,\,\sb\,|\,\Ab,\,\thetab)+\log\,p(\Ab)+\log\,p(\thetab) \,|\,\xb,\,\Ab',\,\thetab' \right]
\eeq

\item M-step : 
\beq 
\left( \wh{\Ab},\,\wh{\thetab} \right) = \argmax_{\left( \Ab,\,\thetab\right)}{Q \left( \Ab,\,\thetab\,|\,\Ab',\,\thetab' \right)}
\eeq
\end{enumerate}   

The joint distribution is factorized as: $p(\xb,\,\sb,\,\Ab,\,\thetab)\,=\,p(\xb\,|\,\Ab,\,\sb)\,p(\Ab)\,p(\sb\,|\,\thetab)\,p(\thetab)$. We can remark that $p(\xb,\,\sb,\,\Ab,\,\thetab)$ as a function of $(\Ab,\,\thetab)$ is separable in $\Ab$ and $\thetab$. Consequently, the functional is  separated into two factors: one representing an $\Ab$ functional and the other representing a $\thetab$ functional:
\beq
Q \left( \Ab,\,\thetab\,|\,\Ab',\,\thetab' \right) = Q_a \left( \Ab\,|\,\Ab',\,\thetab' \right)+Q_h \left( \thetab\,|\,\Ab',\,\thetab' \right)
\eeq
with:
\beq \left\{ \barr{ccc}
Q_a \left( \Ab\,|\,\Ab',\,\thetab' \right) &= & E\left[  \log\,p(\xb\,|\,\Ab,\,\sb)+\log\,p(\Ab) \,|\,\xb,\,\Ab',\,\thetab' \right] \\
Q_h \left( \thetab\,|\,\Ab',\,\thetab' \right) &= & E\left[  \log\,p(\sb\,|\,\thetab)+\log\,p(\thetab) \,|\,\xb,\,\Ab',\,\thetab' \right] \\
\earr \right.
\eeq 

\subsection*{- Maximisation with respect to $\Ab$}
The functional $Q_a$ is:
\beq
Q_a=\frac{-1}{2\,\sigma_{\epsilon}^2} \sum_{t=1}^{T} E\left[ \left( \xb(t)-\Ab\,\sb(t)\right)^{T}\left( \xb(t)-\Ab\,\sb(t)\right)|\,\xb,\,\Ab',\,\thetab' \right]+ \log\,p(\Ab).
\eeq
The gradient of this expression with respect to the elements of $\Ab$ is:

\beq \label{gradient}
 \frac{\partial{Q_a}}{\partial{\Ab_{i,j}}}=\frac{T}{\sigma_{\epsilon}^2}\left( \wh{\Rb}_{xs}-\Ab\,\wh{\Rb}_{ss} \right)_{i,j}-\frac{1}{\sigma_{aij}^2} \left( \Ab_{i,j}-\Mb_{i,j} \right).
\eeq

where:

\beq \left\{ \barr{ccc}
\wh{\Rb}_{xs} &=& \frac{1}{T}\sum_{t=1}^{T} E\left[  \xb(t)\, \sb(t)^{T}|\,\xb,\,\Ab',\,\thetab' \right]\\
\wh{\Rb}_{ss} &=& \frac{1}{T}\sum_{t=1}^{T} E\left[  \sb(t)\, \sb(t)^{T}|\,\xb,\,\Ab',\,\thetab' \right]
\earr \right.
\eeq

Evaluation  of $\wh{\Rb}_{xs}$ and $\wh{\Rb}_{ss}$ requires the computation of the expectations of $\xb(t)\, \sb(t)^{T}$ and $\sb(t)\, \sb(t)^{T}$. The main computational cost is due to the fact that the expectation of any function $ f\left( \sb \right)$ is given by:
\beq\label{somme}
E\left[\,  f\left( \sb \right)\,|\,\xb,\,\Ab',\,\thetab' \right]=\sum_{\zb' \,\in\, \prod_{i=1}^{n}\Zc_i}E\left[\,  f\left( \sb \right)\,|\,\xb,\zb=\zb',\,\Ab',\,\thetab' \right]\, p(\zb'\,|\,\xb,\,\Ab',\,\thetab').
\eeq

\noindent which involves a sum of $\prod_{j=1}^{n}q\left(j\right)$ terms corresponding to the whole combinations of labels. One way to obtain an approximate but fast estimate of this expression is to limit the summation to only one term corresponding to the MAP estimate of $\zb$:
\beq
E\left[\,  f\left( \sb \right)\,|\,\xb,\,\Ab',\,\thetab' \right]=E\left[\,  f\left( \sb \right)\,|\,\xb,\zb=\wh{\zb}^{MAP},\,\Ab',\,\thetab' \right].
\eeq
Then, given estimated labels $\zb_{1..T}$, the source $\sb(t)$  \apost law is Normal with  mean $\thetab_{xz}$ and variance $\Vb_{xz}$ given by (\ref{moyenne}) and (\ref{variance2}).

The source estimate is then $\thetab_{xz}$. $\wh{R}_{xs}$ and $\wh{R}_{ss}$ become:
\beq
\wh{R}_{xs} = \frac{1}{T}\sum_{t=1}^{T} \xb(t)\, \wh{\sb}(t)^{T}
\eeq
and
\beq
\wh{R}_{ss} = \frac{1}{T}\sum_{t=1}^{T} \wh{\sb}(t)\, \wh{\sb}(t)^{T}+\frac{1}{T}\sum_{t=1}^{T} (\Ab^{t}\Rb_n^{-1}\Ab+ \bf{\Gamma}_z^{-1})^{-1} 
\eeq

When $\Sb_{1..T}$ estimated and  using the matrix operations defined in section II and cancelling the gradient (\ref{gradient}) to zero, we obtain the expression of the estimate of $\Ab$:
\beq\label{melange}
\Abh^{k+1}=Mat\left(\left[ \Lambdab+ T\,\wh{R}_{ss}^{*} \right]^{-1}\left[ \Lambdab Vect(M)+T\,Vect\left( \wh{R}_{xs}\right) \right] \right)
\eeq

\subsection*{- Maximisation with respect to $\thetab$}

With a uniform \aprio for the means, maximisation of $Q_h$ with respect to $m_{jz}$ gives :
\beq\label{hyper1}
\wh{m}_{jz}=\frac{\sum_{t=1}^{T}\theta_{jz}(t)\,p(z(t)\,|\,\xb,\,\Ab',\,\thetab')}{\sum_{t=1}^{T}p(z(t)\,|\,\xb,\,\Ab',\,\thetab')}
\eeq

With an Inverted Gamma prior $\Gc\left(\alpha,\,\beta\right)$ ($\alpha > 0$ et $\beta > 1$) for the variances, the maximisation of $Q_h$ with respect to $\sigma_{jz}$ gives:

\beq\label{hyper2}
\wh{\sigma}_{jz}=\frac{2\,\beta+\sum_{t=1}^{T}\left(V_{jz}+\theta_{jz}^2-2\,\wh{m}_{jz}\theta_{jz}+\wh{m}_{jz}^2\right)p(z(t)\,|\,\xb,\,\Ab',\,\thetab')}{\sum_{t=1}^{T}p(z(t)\,|\,\xb,\,\Ab',\,\thetab')+2\,\left(\alpha-1\right)}
\eeq

\subsection*{Summary of the Penalized EM algorithm}
Based on the preceeding equations, we propose the following algorithm to estimate sources and parameters using the following five steps:
\begin{enumerate}
\item Estimate the hyperparameters according to (\ref{hyper1}) and (\ref{hyper2}).
\item Update of data classification by estimating $\wh{\zb}_{1..T}^{MAP}$.
\item Given this classification, sources estimate is the mean of the Gaussian \apost law (\ref{mean}).
\item Update of data classification.
\item Estimate the mixing matrix $\Ab$ according to the re-estimation equation (\ref{melange}).
\end{enumerate}

\section{Comparison with JMAP algorithm and its sensitivity to initial conditions}
The Penalized EM algorithm has an optimization cost  approximately $2$ times higher, per sample, than the JMAP algorithm. However, both algorithms have a reasonable computational  complexity, linearly increasing with the number of samples.
 Sensitivity to initial conditions is inherent to the EM-algorithm even to the penalized version. In order to illustrate this fact, we simulated the algorithm with the same parameters as in section IV. Note that initial conditions for hyperparameters are $\psi^{(0)}=\left(\barr{cc}1 & 1 \\ 1 & 1\earr\right)$ and  $m^{(0)}=\left(\barr{cc}0 & 0 \\ 0 & 0\earr\right)$. However, the Penalized EM algorithm fails in separating sources (see figure 11). We note then that JMAP algorithm is more robust to initial conditions.

\bigskip
\noindent\btabu{@{}c@{}}
\includegraphics[width=\textwidth,height=40mm]{fig3em150}\\
\hspace*{0.3cm}(a)\hspace*{4.5cm}(b)\hspace*{4.5cm}(c)\\~\\
Figure 11- Results of separation with the Penalized EM algorithm: \\
 (a) Phase space distribution of sources,\\ 
 (b) mixed signals and (c) separated sources
\etabu

We modified the initial condition to have means:  
$m^{(0)}=\left(\barr{cc}
-0.5 & 0.5 \\ 
-0.5 & 0.5
\earr\right)$. 
We noted, in this case,  the convergence of the Penalized EM algorithm 
to the correct solution.  
Figures $12$-$16$ illustrate the separation results:

\bigskip
\noindent\btabu{@{}c@{}}
\includegraphics[width=\textwidth,height=40mm]{fig3em15b}\\
\hspace*{0.3cm}(a)\hspace*{4.5cm}(b)\hspace*{4.5cm}(c)\\~\\
Figure 12- Results of separation with the Penalized EM algorithm: \\
 (a) Phase space distribution of sources,\\ 
 (b) mixed signals and (c) separated sources
\etabu

\bigskip
\noindent\btabu{@{}c@{}c@{}}
\includegraphics[width=\textwidth/2,height=40mm]{fig7em15b}&
\includegraphics[width=\textwidth/2,height=40mm]{fig6em15b}\\
Figure 13- Evolution of classification error&
Figure 14- Evolution of index
\etabu

\bigskip
\noindent\btabu{@{}c@{}c@{}}
\includegraphics[width=\textwidth/2,height=40mm]{fig4em15b}&
\includegraphics[width=\textwidth/2,height=40mm]{fig5em15b}\\
Figure 15- Identification of $m_{11}$&
Figure 16- Identification of $\psi_{11}$
\etabu

\section{Conclusion}
We have proposed solutions to source separation problem using a Bayesian framework. Specific aspects of the described approach include:
\bit
\item Taking account of errors on model and measurements.
\item Introduction of \aprio distribution for the mixing matrix and hyperparameters. This was motivated by two different reasons: Mixing matrix prior should exploit  previous information and variances prior should  regularize the log-posterior objective function.
\eit
We then consider the problem in terms of a  mixture of Gaussian priors to develop a hierarchical strategy for source estimation. This same interpretation leads us to classify data before estimating hyperparameters and to reduce computational cost in the case of the proposed Penalized EM algorithm.
{\small 

\bibliographystyle{maxent95}
%\bibliography{biben,revuedef,revueabr,baseAJ,baseKZ,gpipubli}

\begin{thebibliography}{1}

\bibitem{Moulines97}
E.~Moulines, J.~Cardoso, and E.~Gassiat,
\newblock ``Maximum likelihood for blind separation and deconvolution of noisy
  signals using mixture models'',
\newblock in {\em ICASSP-97}, M{\"u}nchen, Germany, {A}pril 1997.

\bibitem{Djafari99a}
A.~Mohammad-Djafari,
\newblock ``A {B}ayesian approach to source separation'',
\newblock in {\em MaxEnt99 Proceedings}. 1999, {K}luwer.

\bibitem{Macchi99}
O.~Macchi and E.~Moreau,
\newblock ``Adaptative unsupervised separation of discrete sources'',
\newblock in {\em Signal Processing 73}, 1999, pp. 49--66.

\bibitem{Bermond00}
O.~Bermond,
\newblock {\em M\'ethodes statistiques pour la s\'eparation de sources},
\newblock PhD thesis, Ecole Nationale Sup\'erieure des T\'el\'ecommunications,
  {J}anuary 2000.

\bibitem{Ridolfi99}
A.~Ridolfi and J.~Idier,
\newblock ``Penalized maximum likelihood estimation for univariate normal
  mixture distributions'',
\newblock in {\em {A}ctes du 17$^{\mathrm{e}}$ colloque {GRETSI}}, Vannes,
  France, {S}eptember 1999, pp. 259--262.

\bibitem{Cardoso96}
J.~Cardoso and B.~Labeld,
\newblock ``Equivariant adaptative source separation'',
\newblock {\em {S}ignal {P}rocessing}, vol. 44, pp. 3017--3030, 1996.

\bibitem{Moreau96}
E.~Moreau and O.~Macchi,
\newblock ``High-order contrasts for self-adaptative source separation'',
\newblock in {\em Adaptative Control Signal Process. 10}, 1996, pp. 19--46.

\bibitem{Redner84}
R.~A. Redner and H.~F. Walker,
\newblock ``Mixture densities, maximum likelihood and the \textsc{em}
  algorithm'',
\newblock {\em SIAM Rev.}, vol. 26, no. 2, pp. 195--239, {A}pril 1984.

\end{thebibliography}

}
\edoc


physics0111018.extracted_bib
\begin{thebibliography}{999}} \def\ebib{\end{thebibliography}