Lab5_K_means_Clustering.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Lab5 - Machine Learning\n",
    "## K-means Clustering"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In this lab session we will implement the K-means algorithm. K-means is one of the most well-known algorithms for data clustering.\n",
    "\n",
    "Given data $\\mathbf{x} = \\{x_1, x_2,\\ldots,x_n\\}$ where each $x_n \\in \\mathbb{R}^D$ we want to group these datapoints in $K$ clusters ($K \\ll N$), where each cluster is represented by a mean vector $\\mu_k \\in \\mathbb{R}^D$.\n",
    "\n",
    "Our purpose is to find a way to divide the data in clusters, finding appropriate values for the means $\\mu_k, k=1,\\ldots,K$\n",
    "\n",
    "\n",
    "### Kmeans Algorithm steps\n",
    "<ol>\n",
    "  <li>From $ n =1 $ to $N$, assign each data instance in the cluster with the closest center $\\mu_k$</li>\n",
    "    \n",
    "  <br>\n",
    "  <li>Calculate the new cluster centers $$\\mu_k = \\frac{\\sum_{n=1}^N r_{nk} \\, x_n }{\\sum_{n=1}^N r_{nk}}$$ </li>\n",
    "  <br>\n",
    "  <li>Check for convergence and stop, otherwise go to step 1</li>\n",
    "</ol>\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$ r_{nk} \\in \\{0,1\\}, \\, \\sum_{k=1}^K r_{nk}=1$$\n",
    "so when $x_n$ belongs to cluster $k$, then $r_{nk}=1$ and $r_{nj}=0$ for each $j \\neq k$.\n",
    "\n",
    "The cost function \n",
    "$$ J(r_1,\\ldots,r_N, \\mu_1,\\ldots,\\mu_K))= \\sum_{n=1}^N \\sum_{k=1}^K r_{nk} \\lVert \\mathbf{x}_n - \\mathbf{\\mu}_k \\rVert ^2$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from __future__ import division\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "%matplotlib inline\n",
    "np.random.seed(0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "N = 2000\n",
    "X = np.zeros((N, 2))\n",
    "p = [0.3, 0.6, 1]\n",
    "Mtrue = np.array([[0,0], [-2,3], [2,3]])\n",
    "for i in range(N):\n",
    "    u = np.random.rand()\n",
    "    if u < p[0]:\n",
    "        k = 0\n",
    "    elif u>=p[0] and u < p[1]:\n",
    "        k = 1\n",
    "    else:\n",
    "        k = 2\n",
    "    X[i, :] = Mtrue[k]+ 0.5*np.random.randn(2)\n",
    "N, D  = X.shape\n",
    "K = 3\n",
    "Minit = np.array([np.mean(X, axis=0),]*K) + np.random.randn(K,2)\n",
    "print Minit  \n",
    "    "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Plot the data and the initial values of the centers\n",
    "plt.plot(X[:, 0], X[:, 1], 'o', color='lightgray', markersize = 1)\n",
    "plt.plot(Minit[:, 0], Minit[:, 1], 'b+', mew = 3, ms=25)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from matplotlib.colors import cnames as mcolors\n",
    "def plot_clusters(X, r, k, M):\n",
    "    if k >3:\n",
    "        colors = mcolors.keys()\n",
    "    else:\n",
    "        colors = ['r', 'g', 'b']\n",
    "    for k in range(K):\n",
    "        plt.plot(X[r[:,k]==1, 0], X[r[:,k]==1, 1], '.', color=colors[k], markersize = 3)\n",
    "    plt.plot(M[:, 0], M[:, 1], 'b+', mew = 3, ms=25)\n",
    "    plt.show()\n",
    "\n",
    "def plot_costs(costs):\n",
    "    x = range(1, len(costs)+1)\n",
    "    y = costs\n",
    "    plt.plot(x, y)\n",
    "    plt.ylabel('cost')\n",
    "    plt.xlabel('iterations')\n",
    "    plt.title(\"Cost Function =\")\n",
    "    plt.xticks(x)\n",
    "    plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def ml_kmeans(X, M):\n",
    "    N, D = X.shape\n",
    "    K = M.shape[0]\n",
    "    # Apply the two steps of K means until convergence\n",
    "    tol = 1e-6\n",
    "    Jold = np.inf\n",
    "    maxIters = 100    \n",
    "    M = np.copy(Minit)\n",
    "    costs = []\n",
    "    for it in range(maxIters):        \n",
    "        r = np.zeros((N, K))\n",
    "        # Step 1 -- Assgin data to clusters\n",
    "        # ******************************************************************\n",
    "        # **********************Your code here *****************************\n",
    "        # ******************************************************************\n",
    "                \n",
    "               \n",
    "        # Step 2 -- Update mean centers\n",
    "        # ******************************************************************\n",
    "        # **********************Your code here *****************************\n",
    "        # ******************************************************************    \n",
    "        \n",
    "        # Step 3 -- Calculate cost function and check for convergence\n",
    "        \n",
    "        \n",
    "        # ******************************************************************\n",
    "        # **********************Your code here *****************************\n",
    "        # ******************************************************************\n",
    "        print(\"Iteration #{}, Cost function value: {}\".format(it, J))\n",
    "        \n",
    "        plot_clusters(X, r, k, M)\n",
    "    plot_costs(costs)\n",
    "    return M, r\n",
    "        "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "M, r = ml_kmeans(X, Minit)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "x = range(1, len(costs)+1)\n",
    "y = costs\n",
    "plt.plot(x, y)\n",
    "plt.ylabel('cost')\n",
    "plt.xlabel('iterations')\n",
    "plt.title(\"Cost Function =\")\n",
    "plt.xticks(x)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Plot the clustering and the centers\n",
    "fig = plt.figure()\n",
    "f, ax = plt.subplots(1,2)\n",
    "f.set_figwidth(12)\n",
    "ax[0].plot(X[:, 0], X[:, 1], 'o', color='lightgray', markersize = 1)\n",
    "ax[0].plot(Minit[:, 0], Minit[:, 1], 'b+', mew = 3, ms=25)\n",
    "if k >3:\n",
    "    colors = mcolors.keys()\n",
    "else:\n",
    "    colors = ['r', 'g', 'b']\n",
    "for k in range(K):\n",
    "    ax[1].plot(X[r[:,k]==1, 0], X[r[:,k]==1, 1], '.', color=colors[k], markersize = 3)\n",
    "ax[1].plot(M[:, 0], M[:, 1], 'b+', mew = 3, ms=25)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print M\n",
    "print Minit"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 2",
   "language": "python",
   "name": "python2"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.12"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}