31  Analyse en composantes principales kernelisée (Kernel PCA)

L’analyse en composantes principales est une technique de réduction de dimension qui vise à créer une nouvelle représentation des données en plus petite dimension par une transformation simple et linéaire. La version kernelisée permet de représenter des jeux de données plus complexes de manière plus compacte en appliquant une transformation non linéaire.

Cette transformation non linéaire est en fait construite à partir de la transformation linéaire classique de la PCA appliquée aux données préalablement projetées de manière non linéaire dans un nouvel espace de représentation. La kernelisation permet de bénéficier des techniques des méthodes à noyaux pour effectuer les calculs dans ce nouvel espace potentiellement en très grande dimension de manière efficace.

La méthode procède donc ainsi :

  1. Projeter les données dans un nouvel espace de représentation : \boldsymbol{X}=\begin{bmatrix}\boldsymbol{x}_1^\top\\\vdots\\\boldsymbol{x}_n^\top \end{bmatrix} \in\mathbb{R}^{n\times d} \qquad \mapsto \qquad \boldsymbol{\Phi }= \begin{bmatrix}\text{------}\ \phi(\boldsymbol{x}_1)^\top \text{------} \\\vdots\\\text{------}\ \phi(\boldsymbol{x}_n)^\top \text{------} \end{bmatrix} \in \mathbb{R}^{n\times d_{\phi}}
  2. Appliquer la PCA dans le nouvel espace : \boldsymbol{Z} = \boldsymbol{\Phi }\boldsymbol{U}^\top \qquad ou\qquad \boldsymbol{z}_i = \boldsymbol{U} \phi(\boldsymbol{x}_i) , \ i=1,\dots, n avec les lignes de \boldsymbol{U} contenant les vecteurs propres de \boldsymbol{\Phi}^\top \boldsymbol{\Phi} associés aux plus grandes valeurs propres.

Cela crée néanmoins deux problèmes :

Au final, comme démontré et détaillé ci-dessous, l’algorithme Kernel PCA s’implémente avec les grandes étapes suivantes :

  1. Choisir une fonction de noyau K et calculer la matrice de noyau \boldsymbol{K}.
  2. Calculer la matrice de noyau centrée \overline{\boldsymbol{K}}.
  3. Calculer les vecteurs propres \boldsymbol{w}_k de \overline{\boldsymbol{K}} associés aux p plus grandes valeurs propres \lambda_k.
  4. Normaliser les vecteurs propres : \boldsymbol{W} = \begin{bmatrix}\boldsymbol{w}_1^\top/\sqrt{\lambda_1}\\\vdots\\\boldsymbol{w}_p^\top /\sqrt{\lambda_p}\end{bmatrix}
  5. Calculer la nouvelle représentation des données : \boldsymbol{z}_i = \boldsymbol{U} \overline{\phi}(\boldsymbol{x}_i) = \boldsymbol{W} \boldsymbol{\Phi }\overline{\phi}(\boldsymbol{x}_i) = \boldsymbol{W} \begin{bmatrix} \overline{K}(\boldsymbol{x}_1,\boldsymbol{x}_i) \\\vdots\\\overline{K}(\boldsymbol{x}_n, \boldsymbol{x}_i) \end{bmatrix}
    Pour toutes les données d’un coup : \boldsymbol{Z}^\top = \boldsymbol{U} \boldsymbol{\Phi}^\top = \boldsymbol{W} \boldsymbol{\Phi }\boldsymbol{\Phi}^\top = \boldsymbol{W} \overline{\boldsymbol{K}}

31.1 Version duale de la PCA

La version duale de la PCA calcule les nouvelles composantes comme \boldsymbol{z}_i = \boldsymbol{W} \boldsymbol{\Phi }\overline{\phi}(\boldsymbol{x}_i),\qquad avec\quad \boldsymbol{W} = \begin{bmatrix}\boldsymbol{w}_1^\top/\sqrt{\lambda_1}\\\vdots\\\boldsymbol{w}_p^\top /\sqrt{\lambda_p}\end{bmatrix}\in\mathbb{R}^{p\times n} en calculant les vecteurs propres \boldsymbol{w}_k de \boldsymbol{\Phi }\boldsymbol{\Phi}^\top\in\mathbb{R}^{n\times n} au lieu des vecteurs propres \boldsymbol{u}_k de \boldsymbol{\Phi}^\top \boldsymbol{\Phi}\in\mathbb{R}^{d_{\phi}\times d_{\phi}}.

Cela est possible grâce au raisonnement suivant, qui montre que \boldsymbol{z}_i = \boldsymbol{W} \boldsymbol{\Phi }\overline{\phi}(\boldsymbol{x}_i) est équivalent à \boldsymbol{z}_i = \boldsymbol{U} \overline{\phi}(\boldsymbol{x}_i).

Pour tout vecteur propre \boldsymbol{u}_k de \boldsymbol{\Phi}^\top \boldsymbol{\Phi} : \boldsymbol{\Phi}^\top \boldsymbol{\Phi }\boldsymbol{u}_k = \lambda_k \boldsymbol{u}_k et, en posant \boldsymbol{v}_k = \frac{1}{\lambda_k}\boldsymbol{\Phi }\boldsymbol{u}_k\in \mathbb{R}^n cela donne \boldsymbol{u}_k = \frac{1}{\lambda_k}\boldsymbol{\Phi}^\top \boldsymbol{\Phi }\boldsymbol{u}_k = \boldsymbol{\Phi}^\top \left(\frac{1}{\lambda_k}\boldsymbol{\Phi }\boldsymbol{u}_k\right) = \boldsymbol{\Phi}^\top \boldsymbol{v}_k et \boldsymbol{\Phi}\boldsymbol{\Phi}^\top\boldsymbol{v}_k = \frac{1}{\lambda_k}\boldsymbol{\Phi}\boldsymbol{\Phi}^\top \boldsymbol{\Phi}\boldsymbol{u}_k = \frac{1}{\lambda_k}\boldsymbol{\Phi }(\lambda_k\boldsymbol{u}_k) = \boldsymbol{\Phi}\boldsymbol{u}_k = \lambda_k\boldsymbol{v}_k (les deux dernières égalités sont dues au fait que \boldsymbol{u}_k est vecteur propre de \boldsymbol{\Phi}^\top \boldsymbol{\Phi}, et à la définition de \boldsymbol{v}_k).

Donc \boldsymbol{v}_k est bien un vecteur propre de \boldsymbol{\Phi}\boldsymbol{\Phi}^\top associé à la même valeur propre \lambda_k que u_k.

Cependant, si l’on calcule les vecteurs propres de \boldsymbol{\Phi}\boldsymbol{\Phi}^\top avec un ordinateur, nous n’obtiendrons pas exactement les \boldsymbol{v}_k, car ceux-ci ne sont définis qu’à un facteur multiplicatif près. Typiquement, un solveur retourne des vecteurs propres \boldsymbol{w}_k tels que \|\boldsymbol{w}_k\|=1, alors que \|\boldsymbol{v}_k\|^2 = \boldsymbol{v}_k^\top \boldsymbol{v}_k = \frac{1}{\lambda_k^2}\boldsymbol{u}_k^\top \boldsymbol{\Phi}^\top \boldsymbol{\Phi }\boldsymbol{u}_k = \frac{1}{\lambda_k^2} \boldsymbol{u}_k^\top (\lambda_k \boldsymbol{u}_k) = \frac{1}{\lambda_k} \boldsymbol{u}_k^\top \boldsymbol{u}_k = \frac{1}{\lambda_k} Mais il suffit de renormaliser les \boldsymbol{w}_k pour obtenir les \boldsymbol{v}_k = \boldsymbol{w}_k/\sqrt{\lambda_k} : \|\boldsymbol{w}_k\| = 1 \qquad \Rightarrow\qquad \|\boldsymbol{w}_k/\sqrt{\lambda_k}\|^2 = 1/\lambda_k

Au final, on obtient donc les \boldsymbol{u}_k à partir des \boldsymbol{w}_k : \boldsymbol{u}_k = \boldsymbol{\Phi}^\top \boldsymbol{v}_k = \frac{1}{\sqrt{\lambda_k}}\boldsymbol{\Phi}^\top \boldsymbol{w}_k \qquad \text{et}\quad \boldsymbol{U} = \boldsymbol{W} \boldsymbol{\Phi},\quad \boldsymbol{z}_i = \boldsymbol{U} \overline{\phi}(\boldsymbol{x}_i) = \boldsymbol{W} \boldsymbol{\Phi }\overline{\phi}(x_i)

31.1.1 Astuce du noyau

La forme duale de la PCA permet d’éviter le calcul des vecteurs propres en grande dimension. Cependant, il est toujours nécessaire de calculer explicitement la projection des données \boldsymbol{\Phi} et les produits scalaires entre les images \phi(\boldsymbol{x}_i) de tous les points inclus dans le produit \boldsymbol{\Phi}\boldsymbol{\Phi}^\top.

Pour palier à cette difficulté, les méthodes à noyaux calculent ces produits scalaires implicitement comme K(\boldsymbol{x},\boldsymbol{x}') = \left\langle \phi(\boldsymbol{x}), \phi(\boldsymbol{x}') \right\rangle à partir d’une fonction de noyau K valide, par exemple le noyau gaussien : K(\boldsymbol{x},\boldsymbol{x}') = \exp(-\|\boldsymbol{x}-\boldsymbol{x}'\|^2/2\sigma^2).

Il serait donc possible de calculer les vecteurs propres \boldsymbol{w}_k de
\boldsymbol{\Phi}\boldsymbol{\Phi}^\top = \begin{bmatrix} \left\langle \phi(\boldsymbol{x}_1), \phi(\boldsymbol{x}_1) \right\rangle& \dots & \left\langle \phi(\boldsymbol{x}_1), \phi(\boldsymbol{x}_n) \right\rangle\\\vdots& \ddots\\\left\langle \phi(\boldsymbol{x}_n), \phi(\boldsymbol{x}_1) \right\rangle & \dots & \left\langle \phi(\boldsymbol{x}_n), \phi(\boldsymbol{x}_n) \right\rangle\end{bmatrix} à partir de \boldsymbol{K} = \begin{bmatrix} K(\boldsymbol{x}_1, \boldsymbol{x}_1) & \dots & K(\boldsymbol{x}_1,\boldsymbol{x}_n) \\\vdots& \ddots\\K(\boldsymbol{x}_n,\boldsymbol{x}_1) & \dots & K(\boldsymbol{x}_n,\boldsymbol{x}_n) \end{bmatrix} Mais il faut en réalité calculer les vecteurs propres d’une version centrée de la matrice : \begin{align*} \boldsymbol{\Phi}\boldsymbol{\Phi}^\top &= \begin{bmatrix} \left\langle \phi(\boldsymbol{x}_1) - \phi_0, \phi(\boldsymbol{x}_1)- \phi_0 \right\rangle& \dots & \left\langle \phi(\boldsymbol{x}_1)- \phi_0, \phi(\boldsymbol{x}_n)- \phi_0 \right\rangle\\\vdots& \ddots\\\left\langle \phi(\boldsymbol{x}_n)- \phi_0, \phi(\boldsymbol{x}_1)- \phi_0 \right\rangle & \dots & \left\langle \phi(\boldsymbol{x}_n)- \phi_0, \phi(\boldsymbol{x}_n)- \phi_0 \right\rangle\end{bmatrix}\\ &= \begin{bmatrix} \overline{K}(\boldsymbol{x}_1,\boldsymbol{x}_1) & \dots & \overline{K}(\boldsymbol{x}_1,\boldsymbol{x}_n) \\\vdots& \ddots\\\overline{K}(\boldsymbol{x}_n,\boldsymbol{x}_1) & \dots &\overline{K}(\boldsymbol{x}_n,\boldsymbol{x}_n) \end{bmatrix} = \overline{\boldsymbol{K}} \end{align*} où \boldsymbol{\phi}_0 = \frac{1}{n} \sum_{i=1}^n \phi(\boldsymbol{x}_i).

Pour cela, calculons : \begin{align*} \overline{K}(\boldsymbol{x},\boldsymbol{x}') &= \left\langle \phi(\boldsymbol{x}) - \boldsymbol{\phi}_0, \phi(\boldsymbol{x}')- \boldsymbol{\phi}_0 \right\rangle = \left\langle \phi(\boldsymbol{x}), \phi(\boldsymbol{x}') \right\rangle - \left\langle \phi(\boldsymbol{x}) , \boldsymbol{\phi}_0 \right\rangle - \left\langle \boldsymbol{\phi}_0, \phi(\boldsymbol{x}') \right\rangle + \left\langle \boldsymbol{\phi}_0, \boldsymbol{\phi}_0 \right\rangle\\ &= K(\boldsymbol{x},\boldsymbol{x}') - \frac{1}{n} \sum_{i=1}^n K(\boldsymbol{x},\boldsymbol{x}_i) - \frac{1}{n} \sum_{i=1}^n K(\boldsymbol{x}_i,\boldsymbol{x}') + \frac{1}{n^2} \sum_{i=1}^n \sum_{j=1}^n K(\boldsymbol{x}_i,\boldsymbol{x}_j) \end{align*} soit, de manière matricielle : \overline{\boldsymbol{K}} = \boldsymbol{K} - \frac{1}{n} \boldsymbol{K}\boldsymbol{1}\boldsymbol{1}^\top - \frac{1}{n} \boldsymbol{1}\boldsymbol{1}^\top \boldsymbol{K} + \frac{\boldsymbol{1}^\top \boldsymbol{K}\boldsymbol{1}}{n^2}\boldsymbol{1}\boldsymbol{1}^\top Ainsi, tous les calculs nécessaires peuvent se faire efficacement à partir de la matrice de noyau \boldsymbol{K} standard.