{
 "cells": [
  {
   "attachments": {},
   "cell_type": "markdown",
   "id": "88431645",
   "metadata": {},
   "source": [
    "# TP ACP\n",
    "Le but de se TP est de calculer le premier axe de l'analyse en composantes principales. Pour ce faire, nous allons charger un jeu de données, calculer la matrice de covariance et calculer sa décomposition en valeurs propres/vecteurs propres.\n",
    "\n",
    "\n",
    "\n",
    "### Compétences associées :\n",
    "\n",
    "| Numéro        | Compétence           | \n",
    "|:------------- |:--------------------:|\n",
    "|AS401|Connaître la fonction de minimisation associée à l’ACP|\n",
    "|AS404|Savoir interpréter le premier vecteur propre et plus grande valeur propre de la matrice de covariance|\n",
    "|PY304|Être capable de calculer la décomposition valeurs/vecteurs propres d’une matrice|\n",
    "|PY305|Être capable de calculer une projection au sens de l’ACP|\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e8d96222",
   "metadata": {},
   "source": [
    "1) Charger le jeu de données iris de `sklearn`\n",
    "> utilisez la fonction `load_iris` de scikit-learn"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "5db0db6e",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "510cef4d",
   "metadata": {
    "tags": [
     "correction"
    ]
   },
   "outputs": [],
   "source": [
    "from sklearn.datasets import load_iris\n",
    "X,_ = load_iris(return_X_y=True)\n",
    "n,p= X.shape"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "3f94ae5b",
   "metadata": {},
   "source": [
    "2) Calculer une matrice $X_r$ correspondant aux données centrées et réduites. Vérifiez que la normalisation s'est bien passée.\n",
    "\n",
    "> `np.mean(X,axis=0)` permet de calculer la moyenne pour chaque colonne de `X` \n",
    "\n",
    "> `np.std(X,axis=0)` permet de calculer l'écart type pour chaque colonne de `X`\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "e3c3942b",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "826fadde",
   "metadata": {
    "tags": [
     "correction"
    ]
   },
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "Xr = (X - np.mean(X,axis=0))/np.std(X,axis=0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "4c18b2e8",
   "metadata": {
    "tags": [
     "correction"
    ]
   },
   "outputs": [],
   "source": [
    "plt.boxplot(Xr);"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "e030ad0c",
   "metadata": {
    "tags": [
     "correction"
    ]
   },
   "outputs": [],
   "source": [
    "print(np.mean(Xr))# tout est à 0\n",
    "print(np.std(Xr)) # tout est à 1\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "bc8f21c7",
   "metadata": {},
   "source": [
    "3) Calculer la matrice de variance/covariance `cov` associée à $X_r$. Quels points devez vous vérifier ? \n",
    "\n",
    "\n",
    "**NB** : Comme ici nous avons $\\bar{X_r} = 0$, la matrice de variance/covariance est égale à $\\frac{1}{n} X_r^\\top X_r$. Elle est également égale à la matrice de corrélation étant donné que $\\sigma_{x_j} = 1$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "e29da0e9",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "8b5c4cb2",
   "metadata": {
    "tags": [
     "correction"
    ]
   },
   "outputs": [],
   "source": [
    "import matplotlib.pyplot as plt \n",
    "cov =  Xr.T@Xr / n\n",
    "plt.imshow(cov)\n",
    "plt.colorbar();\n",
    "cov"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "8699717a",
   "metadata": {},
   "source": [
    "4) Calculer la décomposition valeurs/vecteurs propres de la matrice `cov` en utilisant la fonction `eigh` du module `numpy.linalg`\n",
    "\n",
    "Question bonus : pourquoi `eigh` plutot que `eig`, `eigvalsh`, ... "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "fcb74ac1",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "5f0701de",
   "metadata": {
    "tags": [
     "correction"
    ]
   },
   "outputs": [],
   "source": [
    "from numpy.linalg import eigh\n",
    "eigenvalues,eigenvectors = eigh(cov)\n",
    "print(eigenvalues)\n",
    "print(eigenvectors)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "9f078f2e",
   "metadata": {},
   "source": [
    "5) Quelle est la valeur de la plus grande valeur propre ? Récupérez le vecteur propre `v` associé. On pourra regarder directement les valeurs ou utiliser la fonction `argmax` de `numpy`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "1dc8df16",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "19a23edf",
   "metadata": {
    "tags": [
     "correction"
    ]
   },
   "outputs": [],
   "source": [
    "i = np.argmax(eigenvalues)\n",
    "print(eigenvalues[i])\n",
    "v = eigenvectors[:,i]\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "191d90b1",
   "metadata": {},
   "source": [
    "6) Projetez les données `X` sur le vecteur $v$. Vous devriez obtenir une matrice de $n$ lignes"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "aacaa138",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "b49a4928",
   "metadata": {
    "tags": [
     "correction"
    ]
   },
   "outputs": [],
   "source": [
    "u = Xr@v\n",
    "u.shape\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "694a6aff",
   "metadata": {},
   "source": [
    "Ce vecteur `u` correspond à une projection des données de $X$ en 1D maximisant la variance. On verra en détail ce point au prochain cours. En attendant, mesurez la variance de `u`. Que remarquez vous ?\n",
    "Comparez à la variance totale.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "ab982659",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "83f19009",
   "metadata": {
    "tags": [
     "correction"
    ]
   },
   "outputs": [],
   "source": [
    "print(f\"Variance de u = {np.var(u):.2f}\")\n",
    "print(f\"Variance de l'ensemble des colonnes {np.var(Xr,axis=0).sum():.2f}\")\n",
    "# Ici il faut observer que la variance du vecteur u  \n",
    "# est inférieure à la somme des variances de chaque colonne (1 car centré réduit)\n",
    "# Ceci indique que le vecteur u explique \"2.92 quantité de variance (information)\" sur les 4 au total\n",
    "# D'autre part, il explique \"plus\" de variance qu'en prenant une colonne toute seule (= 1 car centré réduit)\n",
    "# on approfondira ce point au prochain CM/TD"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "625bbcae",
   "metadata": {},
   "source": [
    "### Avec `sklearn`\n",
    "Maintenant, nous allons faire la même chose avec `sklearn` et vérifier que tout va bien.\n",
    "1) Créer un objet de la classe `PCA` du module `sklearn.decomposition`.  Après avoir lu la documentation, réglez le paramètre `n_components` à 1 afin de reproduire notre code de `numpy`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "ad36b659",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "eedce1e2",
   "metadata": {
    "tags": [
     "correction"
    ]
   },
   "outputs": [],
   "source": [
    "from sklearn.decomposition import PCA\n",
    "pca = PCA(n_components=1)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "9c211ec0",
   "metadata": {},
   "source": [
    "Cette classe suit le modèle de beaucoup de classes de `sklearn`, c'est à dire le couple de méthodes `fit` et `transform`. La méthode `fit` permet de régler les paramètres du modèle, ici le vecteur de projection; et la méthode `transform` permet d'appliquer l'opération sur des données. \n",
    "\n",
    "2) Appliquez la méthode `fit` sur les **mêmes** données que l'exercice précedent et calculez la projection. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "f1a9443c",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "4f7814f4",
   "metadata": {
    "tags": [
     "correction"
    ]
   },
   "outputs": [],
   "source": [
    "pca.fit(Xr)\n",
    "u_sklearn = pca.transform(Xr)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "6fe212b9",
   "metadata": {},
   "source": [
    "3) Vérifiez que le vecteur de projection calculé par `sklearn` est similaire à celui que vous avez calculez. Quelle petite différence notez vous (pas obligatoire) ?\n",
    "\n",
    "Comment interprétez vous les valeurs du vecteur de projection en termes d'importance de variable ?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "83d6802f",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "49a4661a",
   "metadata": {
    "tags": [
     "correction"
    ]
   },
   "outputs": [],
   "source": [
    "print(pca.components_)\n",
    "print(v)\n",
    "\n",
    "#Le signe peut être différent (aka le sens du vecteur). Mais ça ne change rien.\n",
    "# Plus la valeur est importante, plus la variable contribue au premier axe factoriel"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "a05eed4e",
   "metadata": {},
   "source": [
    "Enfin, comparez que les projections calculées par `sklearn` sont équivalentes (au signe près) aux projections que vous avez calculées."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "2e9d0ac2",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "ffc86262",
   "metadata": {
    "tags": [
     "correction"
    ]
   },
   "outputs": [],
   "source": [
    "print(u[:10])\n",
    "print(u_sklearn[:10])"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "fa30ba41",
   "metadata": {},
   "source": [
    "Lors du prochain TP, nous allons considérer plusieurs axes factoriels et tenter de reconstruire les données."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "mlenv",
   "language": "python",
   "name": "mlenv"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.10.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}