5-Clasific-LDA.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 125,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "The autoreload extension is already loaded. To reload it, use:\n",
      "  %reload_ext autoreload\n"
     ]
    }
   ],
   "source": [
    "# Estos dos comandos evitan que haya que hacer reload cada vez que se modifica un paquete\n",
    "%load_ext autoreload\n",
    "%autoreload 2\n",
    "%config InlineBackend.figure_format = 'retina'"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Anotaciones:\n",
    "Desde el lado de Estadistica querriamos ver que tan bueno se aproxima a una gaussiana. Desde ML que tan bien generaliza"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 126,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pandas as pd\n",
    "import numpy as np\n",
    "from matplotlib import pyplot as plt\n",
    "import pandas as pd\n",
    "from collections import Counter\n",
    "from mpl_toolkits.mplot3d import Axes3D"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 127,
   "metadata": {},
   "outputs": [],
   "source": [
    "N = -1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 128,
   "metadata": {},
   "outputs": [],
   "source": [
    "data = pd.read_csv('data/alturas-pesos-mils-train.csv')[:N]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 129,
   "metadata": {},
   "outputs": [],
   "source": [
    "test = pd.read_csv('data/alturas-pesos-mils-test.csv')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 130,
   "metadata": {},
   "outputs": [],
   "source": [
    "from BiVariateJoint import BiVariateJoint"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 131,
   "metadata": {},
   "outputs": [],
   "source": [
    "data_hombres = data.loc[data['Genero']=='Hombre'][['Peso','Altura']].values\n",
    "data_mujeres = data.loc[data['Genero']=='Mujer'][['Peso','Altura']].values"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Grafico los datos"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 101,
   "metadata": {},
   "outputs": [],
   "source": [
    "%matplotlib qt\n",
    "f = plt.figure(figsize=(20,10))\n",
    "plt.scatter(data_hombres[:,0], data_hombres[:,1], color='b', s=2, label='Hombres')\n",
    "plt.scatter(data_mujeres[:,0], data_mujeres[:,1], color='r', s=2, label='Mujeres')\n",
    "plt.xlabel('Pesos [cms]')\n",
    "plt.ylabel('Alturas [cms]')\n",
    "plt.legend()\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Regla de Bayes (Continuo)\n",
    "\\begin{equation}\n",
    "P(Hombre | X_1, X_2) = \\frac{p(X_1,X_2 | Hombre) P(Hombre)}{p(X_1,X_2)} \\\\ \\\\\n",
    "\\end{equation}\n",
    "\n",
    "\\begin{equation}\n",
    "P(Mujer | X_1, X_2) = \\frac{p(X_1,X_2 | Mujer) P(Mujer)}{p(X_1,X_2)}\n",
    "\\end{equation}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Que tenemos?\n",
    "\\begin{equation}\n",
    "P(Hombre) = \\frac{N_{H}}{N_H + N_M} \\quad\n",
    "\\end{equation}\n",
    "\n",
    "\\begin{equation}\n",
    "P(Mujer) = \\frac{N_{M}}{N_H + N_M} \\quad\n",
    "\\end{equation}\n",
    "\n",
    "\n",
    "\\begin{equation}\n",
    "\\large\n",
    "p(X_1,X_2) = p(X_1,X_2 | Mujer) P(Mujer) + p(X_1,X_2 | Hombre) P(Hombre)\n",
    "\\end{equation}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Con estimar:\n",
    "\\begin{equation}\n",
    "P(X_1,X_2 | Mujer) \\quad y \\quad\n",
    "P(X_1,X_2 | Hombre) \n",
    "\\end{equation}\n",
    "\\begin{equation}\n",
    "o\n",
    "\\end{equation}\n",
    "\\begin{equation}\n",
    "p(X_1,X_2 | Mujer) \\quad y \\quad\n",
    "p(X_1,X_2 | Hombre) \n",
    "\\end{equation}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Bayes Gaussiano (Quadratic Discriminant Analisys)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\\begin{equation}\n",
    "P(Hombre | X_1, X_2) = \\frac{p(X_1,X_2 | Hombre) P(Hombre)}{p(X_1,X_2)} \\\\\n",
    "P(Mujer | X_1, X_2) = \\frac{p(X_1,X_2 | Mujer) P(Mujer)}{p(X_1,X_2)} \\\\\n",
    "p(X_1,X_2) = p(X_1,X_2 | Hombre) P(Hombre) + p(X_1,X_2 | Mujer) P(Mujer)\n",
    "\\end{equation}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\\begin{equation}\n",
    "P(C_0 | X_1, X_2)  \\rightarrow p \\\\\n",
    "P(C_1 | X_1, X_2) \\rightarrow (1-p) \\\\\n",
    "\\end{equation}\n",
    "### Log-odds\n",
    "\\begin{equation}\n",
    "\\frac{p}{1-p} = \\frac{P(X_1,X_2 | C_0) P(C_0)}{\n",
    " P(X_1,X_2 | C_1) P(C_1)}\n",
    "\\end{equation}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 132,
   "metadata": {},
   "outputs": [],
   "source": [
    "peso_altura_hombres = data.loc[data['Genero']=='Hombre'][['Peso', 'Altura']].values\n",
    "peso_altura_mujeres = data.loc[data['Genero']=='Mujer'][['Peso', 'Altura']].values"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 133,
   "metadata": {},
   "outputs": [],
   "source": [
    "from scipy.stats import multivariate_normal"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 134,
   "metadata": {},
   "outputs": [],
   "source": [
    "def get_gauss_prob(data, mean_hombres, mean_mujeres, cov_hombres, cov_mujeres):\n",
    "    data_np = data\n",
    "    likelihood_class_1 = multivariate_normal.pdf(data_np, mean_hombres, cov_hombres)\n",
    "    likelihood_class_2 = multivariate_normal.pdf(data_np, mean_mujeres, cov_mujeres)\n",
    "    N_class_1 = len(peso_altura_hombres)\n",
    "    N_class_2 = len(peso_altura_mujeres)\n",
    "    prior_1 = N_class_1/(N_class_1 + N_class_2)\n",
    "    prior_2 = N_class_2/(N_class_1 + N_class_2)\n",
    "    total = likelihood_class_1 * prior_1 + likelihood_class_2 * prior_2\n",
    "    p_class_1 = likelihood_class_1 * prior_1/total\n",
    "    p_class_2 = likelihood_class_2 * prior_2/total\n",
    "    return p_class_1, p_class_2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 135,
   "metadata": {},
   "outputs": [],
   "source": [
    "def get_acc_gauss(data, mean_hombres, mean_mujeres, cov_hombres, cov_mujeres):\n",
    "    p_class_1, p_class_2 = get_gauss_prob(data[['Peso', 'Altura']].values, mean_hombres, mean_mujeres, cov_hombres, cov_mujeres)\n",
    "    return ((p_class_1>p_class_2)==(data['Genero']=='Hombre')).sum()/len(p_class_1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# LDA"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Supone $\\Sigma_0=\\Sigma_1$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 136,
   "metadata": {},
   "outputs": [],
   "source": [
    "mean_hombres = peso_altura_hombres.mean(axis=0)\n",
    "mean_mujeres = peso_altura_mujeres.mean(axis=0)\n",
    "cov_LDA = np.cov(np.vstack([peso_altura_hombres - mean_hombres, peso_altura_mujeres - mean_mujeres]).T)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 137,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[77.94597261 53.52970001]\n",
      " [53.52970001 50.2964779 ]]\n",
      "[[80.71054209 56.29362173]\n",
      " [56.29362173 53.03283768]]\n",
      "[[75.19882927 50.77709625]\n",
      " [50.77709625 47.57064786]]\n"
     ]
    }
   ],
   "source": [
    "print(cov_LDA)\n",
    "print(np.cov(peso_altura_hombres.T)) \n",
    "print(np.cov(peso_altura_mujeres.T))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 138,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.9167395924490561\n"
     ]
    }
   ],
   "source": [
    "acc_train_LDA =  get_acc_gauss(data, mean_hombres, mean_mujeres, cov_LDA, cov_LDA)\n",
    "print(acc_train_LDA)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 139,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.9225\n"
     ]
    }
   ],
   "source": [
    "acc_test_LDA = get_acc_gauss(test, mean_hombres, mean_mujeres, cov_LDA, cov_LDA)\n",
    "print(acc_test_LDA)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 140,
   "metadata": {},
   "outputs": [],
   "source": [
    "N = 300\n",
    "X = np.linspace(data.min()['Peso'], data.max()['Peso'], N)\n",
    "Y = np.linspace(data.min()['Altura'], data.max()['Altura'], N)\n",
    "X, Y = np.meshgrid(X, Y)\n",
    "\n",
    "# Pack X and Y into a single 3-dimensional array\n",
    "pos = np.empty(X.shape + (2,))\n",
    "pos[:, :, 0] = X\n",
    "pos[:, :, 1] = Y\n",
    "\n",
    "# The distribution on the variables X, Y packed into pos.\n",
    "Z_H_LDA = multivariate_normal.pdf(pos, mean_hombres, cov_LDA)\n",
    "Z_M_LDA = multivariate_normal.pdf(pos, mean_mujeres, cov_LDA)\n",
    "Z_LDA, _ = get_gauss_prob(pos, mean_hombres, mean_mujeres, cov_LDA, cov_LDA)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 142,
   "metadata": {},
   "outputs": [],
   "source": [
    "%matplotlib qt\n",
    "from matplotlib import cm\n",
    "from mpl_toolkits.mplot3d import Axes3D\n",
    "from matplotlib.colors import ListedColormap\n",
    "fig = plt.figure(figsize=(20,10))\n",
    "ax = fig.gca(projection='3d')\n",
    "cm = plt.cm.RdBu\n",
    "#cf = ax.contourf(X, Y, Z, 256, alpha=.8, vmin=0., vmax=1., cmap=cm)\n",
    "#plt.colorbar(cf, ax=ax)\n",
    "ax.contourf(X, Y, Z_H_LDA, 256)\n",
    "ax.contourf(X, Y, Z_M_LDA, 256)\n",
    "ax.view_init(70, -90)\n",
    "ax.set_ylabel('Alturas [cms]')\n",
    "ax.set_xlabel('Pesos [kgs]')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Log-Odds"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Partiendo de Bayes (QDA)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\\begin{equation}\n",
    "\\large\n",
    "-\\frac{1}{2}(\\vec{x}-\\vec{\\mu}_0)^T\\Sigma_0^{-1}(\\vec{x}-\\vec{\\mu}_0) + \\frac{1}{2}(\\vec{x}-\\vec{\\mu}_1)^T\\Sigma_1^{-1}(\\vec{x}-\\vec{\\mu}_1)  - ln(|\\Sigma_0|) + ln(|\\Sigma_1|) + log(\\frac{P(C_0)}{P(C_1)}) = 0\n",
    "\\end{equation}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Suponiendo $\\Sigma_0=\\Sigma_1=\\Sigma$\n",
    "### Como las matrices de covarianza son diagonales:  $\\quad \\vec{x}^T\\Sigma^{-1}\\vec{\\mu}_i = \\vec{\\mu}_i^T\\Sigma^{-1}\\vec{x}$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\\begin{equation}\n",
    "\\large\n",
    "-\\frac{1}{2}(\\vec{x}-\\vec{\\mu}_0)^T\\Sigma^{-1}(\\vec{x}-\\vec{\\mu}_0) + \\frac{1}{2}(\\vec{x}-\\vec{\\mu}_1)^T\\Sigma^{-1}(\\vec{x}-\\vec{\\mu}_1) + log(\\frac{P(C_0)}{P(C_1)}) = 0\n",
    "\\end{equation}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\\begin{equation}\n",
    "\\large\n",
    "-\\frac{1}{2}(\\vec{x}^T\\Sigma^{-1}\\vec{x}+\\vec{\\mu}_0^T\\Sigma^{-1}\\vec{\\mu}_0-\\vec{\\mu}_0^T\\Sigma^{-1}\\vec{x}-\\vec{x}^T\\Sigma^{-1}\\vec{\\mu}_0) + \\frac{1}{2}(\\vec{x}^T\\Sigma^{-1}\\vec{x}+\\vec{\\mu}_1^T\\Sigma^{-1}\\vec{\\mu}_1-\\vec{\\mu}_1^T\\Sigma^{-1}\\vec{x}-\\vec{x}^T\\Sigma^{-1}\\vec{\\mu}_1) + log(\\frac{P(C_0)}{P(C_1)}) = 0\n",
    "\\end{equation}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\\begin{equation}\n",
    "\\large\n",
    "-\\frac{1}{2}(\\vec{\\mu}_0^T\\Sigma^{-1}\\vec{\\mu}_0-2\\vec{\\mu}_0^T\\Sigma^{-1}\\vec{x}) + \\frac{1}{2}(\\vec{\\mu}_1^T\\Sigma^{-1}\\vec{\\mu}_1-2\\vec{\\mu}_1^T\\Sigma^{-1}\\vec{x}) + log(\\frac{P(C_0)}{P(C_1)}) = 0\n",
    "\\end{equation}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\\begin{equation}\n",
    "\\large\n",
    "-\\frac{1}{2}(\\vec{\\mu}_0+\\vec{\\mu}_1)^T\\Sigma^{-1}(\\vec{\\mu}_0-\\vec{\\mu}_1)+\\vec{x} ^T\\Sigma^{-1}(\\vec{\\mu}_0- \\vec{\\mu}_1) + log(\\frac{P(C_0)}{P(C_1)}) = 0\n",
    "\\end{equation}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "\\begin{equation}\n",
    "\\large\n",
    "X^TW + C = 0\n",
    "\\end{equation}\n",
    "\n",
    "\\begin{equation}\n",
    "\\large\n",
    "W = \\Sigma^{-1}(\\vec \\mu_0- \\vec \\mu_1)\n",
    "\\end{equation}\n",
    "\n",
    "\\begin{equation}\n",
    "\\large\n",
    "c =  log(\\frac{P(C_0)}{P(C_1)}) - \\frac{1}{2}(\\vec{\\mu}_0+\\vec{\\mu}_1)^T\\Sigma^{-1}(\\vec{\\mu}_0-\\vec{\\mu}_1) = log(\\frac{P(C_0)}{P(C_1)}) - \\frac{1}{2}(\\vec{\\mu}_0+\\vec{\\mu}_1)^T W\n",
    "\\end{equation}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 147,
   "metadata": {},
   "outputs": [],
   "source": [
    "%matplotlib qt\n",
    "fig = plt.figure(figsize=(20,10))\n",
    "ax = fig.gca(projection='3d')\n",
    "cm = plt.cm.RdBu\n",
    "ax.contourf(X, Y, Z_LDA, 256, vmin=0., vmax=1., cmap=cm)\n",
    "ax.set_title('TRAIN: ' + str(np.round(acc_train_LDA*100)/100) + ' - CV:' + str(np.round(acc_test_LDA*100)/100))\n",
    "ax.set_ylabel('Alturas [cms]')\n",
    "ax.set_xlabel('Pesos [kgs]')\n",
    "ax.view_init(70, -90)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 143,
   "metadata": {},
   "outputs": [],
   "source": [
    "W = np.linalg.inv(cov_LDA).dot(mean_hombres - mean_mujeres)\n",
    "C = np.log(peso_altura_hombres.shape[0]/peso_altura_mujeres.shape[0]) - 0.5*(mean_hombres + mean_mujeres).T.dot(np.linalg.inv(cov_LDA)).dot(mean_hombres - mean_mujeres) "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 144,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Verificación ecuación\n",
    "z_teorica = np.zeros([N, N])\n",
    "for row in range(pos.shape[0]):\n",
    "    for col in range(pos.shape[1]):\n",
    "        z_teorica[row, col] = pos[row, col].T.dot(W) + C"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 108,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Verificación\n",
    "%matplotlib qt\n",
    "fig = plt.figure(figsize=(20,10))\n",
    "ax = fig.gca()\n",
    "cm = plt.cm.RdBu\n",
    "cf = ax.contourf(X, Y, Z_LDA, 256, alpha=.8, vmin=0., vmax=1., cmap=cm)\n",
    "plt.colorbar(cf, ax=ax)\n",
    "z_levels = np.logspace(-5,-2,10)/4\n",
    "#ax.contour(X, Y, Z_H, z_levels)\n",
    "#ax.contour(X, Y, Z_M, z_levels)\n",
    "ax.scatter(data_hombres[:,0], data_hombres[:,1], color='b', s=2, label='Hombres')\n",
    "ax.scatter(data_mujeres[:,0], data_mujeres[:,1], color='r', s=2, label='Mujeres')\n",
    "ax.contour(X, Y, Z_LDA, (0.5,), colors='k', linewidths=0.25)\n",
    "#ax.contour(X, Y, z_teorica, (0,), colors='k', linewidths=0.25)\n",
    "ax.set_ylabel('Alturas [cms]')\n",
    "ax.set_xlabel('Pesos [kgs]')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 146,
   "metadata": {},
   "outputs": [],
   "source": [
    "%matplotlib qt\n",
    "fig = plt.figure(figsize=(20,10))\n",
    "ax = fig.gca()\n",
    "cm = plt.cm.RdBu\n",
    "cf = ax.contourf(X, Y, Z_LDA, 256, alpha=.8, vmin=0., vmax=1., cmap=cm)\n",
    "plt.colorbar(cf, ax=ax)\n",
    "z_levels = np.logspace(-5,-2,5)/4\n",
    "ax.contour(X, Y, Z_H_LDA, z_levels, alpha=0.5)\n",
    "ax.contour(X, Y, Z_M_LDA, z_levels, alpha=0.5)\n",
    "ax.contour(X, Y, Z_LDA, (0.5,), colors='k', linewidths=1)\n",
    "ax.scatter(data_hombres[:,0], data_hombres[:,1], color='b', s=1, label='Hombres')\n",
    "ax.scatter(data_mujeres[:,0], data_mujeres[:,1], color='r', s=1, label='Mujeres')\n",
    "ax.set_title('TRAIN: ' + str(np.round(acc_train_LDA*100)/100) + ' - CV:' + str(np.round(acc_test_LDA*100)/100))\n",
    "ax.set_ylabel('Alturas [cms]')\n",
    "ax.set_xlabel('Pesos [kgs]')\n",
    "plt.legend()\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "|Modelo| tipo| Train Acc| CV Acc|Comentarios\n",
    "|-| -| -| -|-|\n",
    "|Histogram - step 0.25|Bayes|0.97|0.68\n",
    "|Histogram - step 0.25|Naive Bayes|0.89|0.88\n",
    "|Histogram - step 0.5|Bayes|0.94|0.82\n",
    "|Histogram - step 0.5|Naive Bayes|0.89|0.88\n",
    "|Histogram - step 1|Bayes|0.92|0.90\n",
    "|Histogram - step 1|Naive Bayes|0.89|0.88\n",
    "|Histogram - step 4|Bayes|0.91|0.92| No generaliza por afuera de la zona de los puntos\n",
    "|Histogram - step 4|Naive Bayes|0.89|0.88\n",
    "|Histogram - step 8|Bayes|0.90|0.90\n",
    "|Histogram - step 8|Naive Bayes|0.87|0.86\n",
    "|Gaussian|Bayes|0.92|0.92|\n",
    "|Gaussian|Naive Bayes|0.89|0.89|\n",
    "|Regresion Logística|-|0.92|0.92|"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "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.6.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}