\section{La SVD}\label{sec:svd}

Nous commencons par quelques rappels d'alg\`ebre lin\'eaire \'el\'ementaire sur valeurs propres et vecteurs propres, qui vont nous \^etre indispensables pour expliquer la d\'ecomposition en valeurs singuli\`eres ou SVD.

\subsection{Quelques rappels d'alg\`ebre lin\'eaire}
On d\'efinit tout d'abord les notions de vecteur propre et de valeur propre associ\'ee.

\begin{dfn}
Soit $A$ une matrice $n\times n$. On dit que le vecteur $x\in\mathbb{R}^n\setminus\{0\}$ est un vecteur propre de $X$ associ\'e \`a la valeur propre $\lambda$ si:
\[
Ax=\lambda x
\]
\end{dfn}

\begin{dfn}
Une matrice $A$ de taille $n\times n$ pour laquelle il existe une base $\mathcal{B}$ de vecteurs propres est dite diagonalisable. On a  alors la d\'ecomposition matricielle:
\[
A=PDP^{-1}
\]
avec $P$ la matrice de passage de la base canonique de $\mathbb{R}^n$ \`a $\mathcal{B}$ et $D$ la matrice diagonale, avec comme coefficients diagonaux les valeurs propres correspondant aux vecteurs propres.
\end{dfn}

\noindent Rappelons que:

\begin{dfn}
Une matrice $A$ de taille $n\times n$ est dite sym\'etrique si elle est \'egale \`a sa transpos\'ee : $A=A^T$. Elle est dite sym\'etrique positive si, pour tout $x\in\mathbb{R}^n$, $x^T\,A x\geq 0$.
\end{dfn}

On a  alors:

\begin{prop}
Toute matrice sym\'etrique est diagonalisable et admet une base de vecteurs propres qui est orthonorm\'ee. Si elle est sym\'etrique positive, ses valeurs propres sont toutes positives ou nulles.
\end{prop}

\subsection{Approximation de rang faible}
La SVD va \^etre bas\'ee sur les notions de valeurs singuli\`eres et vecteurs singuliers d'une matrice rectangulaire. Rappelons tout d'abord le r\'esultat suivant, qui d\'ecoule du fait que $A^TA$ et $A A^T$ sont toutes deux sym\'etriques

\begin{prop}
Soit $A$ une matrice de taille $p\times n$. Les matrices $A^T A$ et $A A^T$ sont toutes deux sym\'etrique positives; elles sont donc toutes les deux diagonalisables \`a spectre positif ou nul. Leurs spectres sont identiques aux valeurs propres nulles pr\`es. On appelle alors valeurs singuli\`eres de la matrice $A$ les racines carr\'ees des valeurs propres positives de $AA^T$ et $A^TA$. Si $A$ est de rang $r$, elle a $r$ valeurs singuli\`eres. 
\end{prop}

On peut alors donner la d\'efinition suivante: 

\begin{dfn}
Soit $A$ une matrice de taille $p\times n$. Les vecteurs singuliers \`a droite de $A$ sont les vecteurs propres de $A^T\,A$, tandis que les vecteurs propres \`a gauche sont les vecteurs propres de $A A^T$.
\end{dfn}

\pagebreak[4]
On a le lien suivant entre vecteurs singuliers \`a gauche et \`a droite:

\begin{prop}
Si $v$ est un vecteur singulier \`a droite, $u=Av$ est un vecteur singulier \`a gauche de $A$, tandis que si $u$ est un vecteur singulier \`a droite, $v=A^T u$ est un vecteur singulier \`a droite de $A$. 
\end{prop}

\noindent On peut maintenant d\'efinir la SVD d'une matrice:

\begin{thm}
Soit $A$ une matrice de taille $p\times n$ et $s_1\geq\cdots\geq s_r>0$ ses valeurs singuli\`eres rang\'ees dans l'ordre d\'ecroissant. On note $U$ et $V$ les matrices de taille $p\times n$ et $p\times n$ dont les colonnes sont les vecteurs propres de $AA^T$ et $A^T\,A$ associ\'es aux valeurs propres $s_1^2,\cdots,s_r^2$. Notons $\Delta_r=\begin{pmatrix}s_1&&0\\&\ddots&\\0&&s_r\end{pmatrix}$. On a alors:
\begin{itemize}
\item Si $p< n$:
\[
A=U\Sigma V\mbox{ avec }\Sigma=\begin{pmatrix}\Delta_r &0\end{pmatrix}
\]
\item Si $p=n$:
\[
A=U\Sigma V\mbox{ avec }\Sigma=\Delta_r
\]  
\item Si $p>n$:
\[
A=U\Sigma V\mbox{ avec }\Sigma=\begin{pmatrix}\Delta_r\\0\end{pmatrix}
\]
\end{itemize}
C'est ce qu'on appelle la d\'ecomposition en valeurs singuli\`eres de la matrice $A$ (SVD en anglais). Elle est unique \`a permutation des vecteurs singuliers pr\`es.
\end{thm}

\noindent On peut maintenant faire le lien entre SVD et meilleure approximation par une matrice de rang faible. Pour cela, notons $\Delta_m=\begin{pmatrix}s_1&&0\\&\ddots&\\0&&s_m\end{pmatrix}$. \`A toute matrice $A$ et pour tout $m\leq r=rang(A)$, on associe la matrice $A_m$ suivante:
\[
A_m=U\Sigma_m V\mbox{ o\`u }\Sigma_m=\begin{pmatrix}\Delta_m &0\end{pmatrix},\,\mbox{  resp }\Delta_m\mbox{ ou }\begin{pmatrix}\Delta_m\\0\end{pmatrix}
\]
suivant que $p<n$ (resp $p=n$ ou $n<p$). Cette d\'ecomposition s'appelle la SVD tronqu\'ee de la matrice $A$ \`a l'ordre $m$.

\noindent Rappelons que, pour toute matrice $A=(a_{i,j})$ sa norme de Frobenius est d\'efinie par: 
%
\[
\|A\|_2=\sqrt{\sum_{i,j}a_{ij}^2}
\] 
%
Nous pouvons maintenant \'enoncer le r\'esultat suivant:

\begin{prop}
Soit $A_m$ la matrice d\'efinie pr\'ec\'edemment. Alors:
\[
\|A-A_m\|_2=\min_{M,\,rang(M)=m}\|A-M\|_2
\]
\end{prop}

\subsection{SVD et analyse en composantes principales}

Nous allons maintenant rappeler ce qu'est l'analyse en composantes principales et expliciter son lien \`a la SVD. On suppose que l'on a $n$ observations $x_1,\cdots,x_n$ qui sont la r\'ealisation d'un vecteur al\'eatoire de $\mathbb{R}^p$. On les stocke dans une matrice $X$ de taille $p\times n$. Pour simplifier, on suppose que cette derni\`ere est centr\'ee i.e. que la moyenne de chaque colonne a \'et\'e soustraite et que donc toutes les colonnes sont de moyenne nulle. 

Traditionellement, pour faire l'analyse en composantes principales de ces donn\'ees, on consid\`ere la matrice de covariance $C$ empirique, d\'efinie par $C=X^T\, X/(n-1)$. Cette derni\`ere \'etant sym\'etrique, on peut la diagonaliser et on a donc une \'ecriture de la forme $C=VLV^T$ o\`u $V$ est la matrice des vecteurs propres et $L$ une matrice diagonale o\`u on a rang\'e les valeurs propres $\lambda_i$ par ordre d\'ecroissant. Rappelons que les vecteurs propres sont ce qu'on appelle les directions ou axes principaux des donn\'ees. On projette alors les donn\'ees sur ces axes principaux et on obtient les composantes dites principales. La $j$-\`eme composante principale est donn\'ee par la $j$-\`eme colonne $XV$. 

On peut maintenant comparer cette information avec celle que nous apporte la SVD de $X$. On \'ecrit  $X=USV^T$ o\`u $S$ est une matrice diagonale contenant les valeurs singuli\`eres $s_i$. On a alors: 
\[
C=VSU^TUSV^T/(n-1)=V\times[S^2/(n-1)]\times V^T
\] 
Cela signifie que les vecteurs singuliers \`a droite de $V$ sont les directions principales et que les valeurs singuli\`eres sont reli\'ees aux valeurs propres de la matrice de covariance par la relation $\lambda_i=s^2_i/(n-1)$. Les composantes principales sont donn\'ees par $XV=USV^TV=US$.

Tout ceci n'est valable que si $X$ est centr\'ee. Ce n'est que dans ce cas que la matrice de covariance est \'egale \`a $X^TX/(n-1)$. Si on veut faire l'analyse en composantes principales \`a partir d'une matrice de corr\'elation (au lieu de la matrice de covariance), les colonnes de $X$ ne doivent pas uniquement \^etre centr\'ees mais aussi normalis\'ees.
\subsection{Algorithme pour d\'eterminer la SVD d'une matrice}
De nombreux algorithmes ont \'et\'e propos\'es pour le calcul num\'erique de la SVD d'une matrice $A$. Nous en pr\'esentons les plus simples et renvoyons \`a la bibliographie pour des versions plus \'elabor\'ees. \`A chaque fois, une \'etape pr\'eliminaire de la m\'ethode consiste \`a bidiagonaliser la matrice $A$ sous la forme $A=UBV$ o\`u $U,V$ sont deux matrices orthogonales et $B$ une matrice bidiagonale en utilisant ce qu'on appelle des r\'eflexions de Householder. Nous commencons par expliquer cette premi\`ere \'etape, puis nous pr\'esentons diff\'erentes alternatives pour en d\'eduire num\'eriquement une approximation de la SVD d'une matrice.
\subsubsection{Notions pr\'eliminaires sur les matrices orthogonales}
Nous commencons par rappeler ce qu'est une matrice orthogonale.

\begin{dfn}
Soit $A$ une matrice de taille $n\times n$. Elle est dite orthogonale si $A\,A^T=A^T\,A=I_n$ o\`u $I_n$ est la matrice identit\'e de taille $n$. Cette d\'efinition est aussi \'equivalente \`a dire que les colonnes de $A$ (ou ses lignes) forment une base orthonorm\'ee de $\mathbb{R}^n$.
\end{dfn}
Donnons deux exemples fondamentaux de transformations orthogonales, qui nous seront utiles dans la suite: \\
\begin{dfn}
Soit $\mathcal{H}$ un hyperplan. La sym\'etrie orthogonale par rapport \`a $\mathcal{H}$ est la transformation qui associe \`a tout vecteur $\mathbf{x}$:
\[
s(x)=x_{\mathcal{H}}-x_{\mathcal{H}^\perp}
\]
o\`u $(x_{\mathcal{H}},x_{\mathcal{H}^\perp})$ sont respectivement les projections orthogonales de $x$ sur $\mathcal{H}$ et $\mathcal{H}^\perp$.
\end{dfn}

\begin{dfn}
La matrice de la forme: 
\[
H_{\mathbf{a}}=Id-2\frac{\mathbf{a} \mathbf{a}^T}{\|\mathbf{a}\|^2_2}\mbox{ o\`u }\mathbf{a}\in\mathbb{R}^n\setminus\{0\}
\]
est dite de Householder associ\'ee au vecteur $\mathbf{a}$. C'est la matrice de la sym\'etrie orthogonale par rapport \`a $Vect(\mathbf{a})^{\perp}$.
\end{dfn}
 
\subsubsection{Premi\`ere \'etape de la SVD : la bidiagonalisation de Golub-Kahan}
Nous donnons maintenant un algorithme fournissant la d\'ecomposition d'une matrice. La m\'ethode de Householder consiste \`a multiplier $A$ successivement par des matrices $Q_1,\cdots,Q_k$ de la forme $Q_k=\begin{pmatrix}I&0\\0&\Gamma\end{pmatrix}$ o\`u $\Gamma$ est une matrice orthogonale, afin d'introduire des z\'eros en dessous de la diagonale de $A$ et au-dessus de la premi\`ere superdiagonale au fur et \`a mesure.  Remarquons que, si $A$ est diagonale par blocs de la forme $\begin{pmatrix}A_{11}&A_{12}\\0&A_{22}\end{pmatrix}$, on a: 
\[
\begin{pmatrix}I&0\\0&\Gamma\end{pmatrix}\times\begin{pmatrix}A_{11}&A_{12}\\0&A_{22}\end{pmatrix}=\begin{pmatrix}A_{11}&A_{12}\\0&\Gamma A_{22}\end{pmatrix}
\]
Si on veut cr\'eer une nouvelle colonne de z\'eros sur la diagonale principale et au-dessus de la premi\`ere superdiagonale, on voit qu'il suffit que la matrice $\Gamma A_{22}$ ait sa premi\`ere colonne de la forme $\begin{pmatrix}\star\\\star\\0\\\vdots\\0\end{pmatrix}$. Si  on appelle $\mathbf{a}=\begin{pmatrix}a_1\\\vdots\\a_p\end{pmatrix}$ la premi\`ere colonne de $A_{22}$ et si on choisit pour $\Gamma$ la matrice de reflexion orthogonale $S$ qui envoie $\mathbf{a}$ sur $\mathbf{a}'=\begin{pmatrix}a_1\\\sqrt{\sum_{i=2}^p a_i^2}\\0\\\vdots\\0\end{pmatrix}$ qui est un vecteur de m\^eme norme que $\mathbf{a}$, on voit qu'on a bien les propri\'et\'es voulues.\\

\noindent Ceci m\`ene ainsi \`a l'algorithme suivant, dit de r\'eduction de Householder d'une matrice sous forme bidiagonale. L'algorithme est sous sa forme g\'en\'erale (cas $p\geq n$, sinon on passe \`a la transpos\'ee) :
 
\restylealgo{algoruled}
\begin{algorithm}[H]
    \SetKwInOut{Entr}{Entrée}
    \SetKwInOut{Sort}{Sortie}
    \SetKwInOut{Init}{Initialisation}
\caption{R\'eduction de Householder d'une matrice $A$ sous forme bidiagonale}
\Entr{\begin{itemize}\item $p,n$, matrice $A$ de taille $p\times n$\end{itemize}}
       \textbf{Initialisation :}
       \begin{itemize}\item $U,V \leftarrow I_n,I_p$\end{itemize}
\Pour{$k=1,\cdots,n $}
{ \begin{itemize}\item Trouver la matrice orthogonale $Q$ qui transforme $\begin{pmatrix}\text{\Large 0} & b_{k-1,k} & \cdots & b_{p,k}\end{pmatrix}$ en $\begin{pmatrix}\text{\Large 0} & b_{k-1,k} & \sqrt{\sum_{j=k}^p b_{j,k}^2} & \text{\Large 0}\end{pmatrix}$
\item $B \leftarrow QB$
\end{itemize}
\Si{$k\leq n-2$ est pair}
{
\begin{itemize}
\item Trouver la matrice orthogonale $P$ qui transforme \\ $\begin{pmatrix}\text{\Large 0} &b_{k,k}&\cdots b_{k,n}\end{pmatrix}$ en $\begin{pmatrix}0&\cdots&b_{k,k}&\sqrt{\sum_{j=k}^n b_{k,j}^2}&\text{\Large 0}\end{pmatrix}$
\item $B \leftarrow BP$
\item $V\leftarrow PV$
\end{itemize}
}
}
\Sort{$B$ matrice bidiagonale sup\'erieure, $U$ et $V$ matrices orthogonales telles que $A=UBV$} 
\end{algorithm}

\subsubsection{Algorithme calculant la SVD d'une matrice bidiagonale}
Une deuxi\`eme \'etape est maintenant n\'ecessaire pour faire d\'ecro\^itre les valeurs situ\'ees sur la premi\`ere superdiagonale de $B$. Celle-ci est bas\'ee sur la d\'ecomposition $QR$.

\begin{prop}
Soit $B$ une matrice de taille $p\times n$. Une d\'ecomposition $QR$ de $B$, est une d\'ecomposition de la forme $B=QR$ o\`u $Q$ est orthogonale et $R$ triangulaire sup\'erieure.  Elle est unique si on impose les signes des coefficients situ\'es sur la diagonale de $R$.
\end{prop} 

Pour d\'eterminer num\'eriquement la d\'ecomposition QR d'une matrice $B$, on utilise l'algorithme suivant:

\restylealgo{algoruled}
\begin{algorithm}[H]
    \SetKwInOut{Entr}{Entrée}
    \SetKwInOut{Sort}{Sortie}
    \SetKwInOut{Init}{Initialisation}

\caption{D\'ecomposition QR d'une matrice $B$}
\Entr{\begin{itemize}\item $p,n$, matrice $B$ de taille $p\times n$\end{itemize}}
       \textbf{Initialisation :}
       \begin{itemize}\item $U,V \leftarrow I_n,I_p$
       \end{itemize}
\Pour{$k=1,\cdots,n $}
{
\begin{itemize}
\item Trouver la matrice orthogonale $Q$ qui transforme $\begin{pmatrix}\text{\Large 0} & b_{k-1,k} & \cdots & b_{p,k}\end{pmatrix}$ en $\left(\text{\Large 0}, b_{k-1,k}, \sqrt{\sum_{j=k}^p b_{j,k}^2},\text{\Large 0}\right)$
\item $B \leftarrow QB$
\item $U\leftarrow UQ$
\end{itemize}
\Si{$k\leq n-2$ est pair}
{
\begin{itemize}
\item Trouver la matrice orthogonale $P$ qui transforme \\ $\begin{pmatrix}\text{\Large 0} &b_{k,k}&\cdots b_{k,n}\end{pmatrix}$ en $\begin{pmatrix}0&\cdots&b_{k,k}&\sqrt{\sum_{j=k}^n b_{k,j}^2}&\text{\Large 0}\end{pmatrix}$
\item $B \leftarrow BP$
\item $V\leftarrow PV$
\end{itemize}
}
}
\Sort{$v_1,\cdots,v_n$ vecteurs non nuls et $R$ matrice orthogonale tels que $B=H_{v_n}\cdots H_{v_1}R$}
\end{algorithm}

\noindent Cette d\'ecomposition QR va nous permettre de calculer num\'eriquement une approximation de la SVD d'une matrice bidiagonale. Construisons par r\'ecurrence, pour une matrice $B$ bidiagonale donn\'ee, la suite suivante. Le premier terme $B_0$ est \'egal \`a $B$. Supposons ensuite avoir construit par r\'ecurrence $B_1,\cdots,B_k$. On pose alors:
\[
B_{k+1}=R_{k}Q_{k}\mbox{ o\`u }(Q_k,R_k)\mbox{ est la d\'ecomposition QR de }B_k\;.
\]
Remarquons que pour tout $k$:
\[
B=Q_0^T\cdots Q_k^T \, B_k Q_k\cdots Q_0
\]
et que les deux matrices $U_k=Q_0^T\cdots Q_k^T$ et $V_k=Q_k\cdots Q_0$ sont orthogonales comme produits de matrices orthogonales.

\begin{prop}
Soit $B$ une matrice bidiagonale. La suite $(B_k)$ de matrices d\'efinie pr\'ec\'edemment converge vers une matrice $\Sigma$ diagonale quand $k\to\infty$. 
\end{prop}

Cette proposition donne donc asymptotiquement la SVD recherch\'ee. En pratique, on consid\`ere $k$ tel que les valeurs de la super-diagonale de $B_k$ soient inf\'erieures \`a un certain seuil $\varepsilon>0$ fix\'e \`a l'avance et on obtient donc une approximation de la SVD de $B$ par :
\[
B\approx U_k B_k V_k
\]
$A$ s'exprimant en fonction de $B$ comme $A=UBV$, on en d\'eduit que: 
\[
A\approx U\Sigma V\mbox{ avec }U=PU_k,\,\Sigma_k= B_k\mbox{ et }V=V_k Q\;.
\] 

%\subsection{Extensions}
D'autres variantes ont \'et\'e propos\'ees dans la litt\'erature pour les matrices de grandes taille ou les matrices creuses. On renvoie \`a~\cite{berry2006parallel}.

On peut aussi proposer des versions parall\`eles de la SVD (voir par exemple le chapitre 4 de~\cite{berry2006parallel} dans \cite{kontoghiorghes2005handbook}).