knn_distance_metrics.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Investigation on k-Nearest Neighbors Model\n",
    "\n",
    "***By: Rabindra Nepal***\n",
    "\n",
    "\n",
    "In this notebook, we explore the k-NN model based on different distance metrics: **Minkowski and Mahalanobis** distances. Also, we will look into the average performances of the model with two different distances. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Dependencies\n",
    "\n",
    "import math\n",
    "import time\n",
    "import numpy as np\n",
    "from numpy import linalg\n",
    "import pandas as pd\n",
    "\n",
    "import matplotlib.pyplot as plt\n",
    "plt.style.use('_classic_test')\n",
    "%matplotlib inline\n",
    "\n",
    "# searborn: plotting\n",
    "import seaborn as sns\n",
    "from scipy import stats\n",
    "\n",
    "# sklearn\n",
    "from sklearn.utils import shuffle\n",
    "from sklearn.model_selection import train_test_split\n",
    "from sklearn.neighbors import KNeighborsClassifier\n",
    "from sklearn.model_selection import GridSearchCV\n",
    "from sklearn.metrics import confusion_matrix\n",
    "from sklearn.model_selection import cross_validate, cross_val_predict, cross_val_score\n",
    "from sklearn.metrics import precision_score, recall_score, f1_score, accuracy_score\n",
    "from sklearn.metrics import roc_curve, roc_auc_score\n",
    "\n",
    "from sklearn.preprocessing import scale"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "# warnings suppression\n",
    "import warnings\n",
    "warnings.simplefilter(action='ignore', category=RuntimeWarning)\n",
    "warnings.simplefilter(action='ignore', category=DeprecationWarning)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## k-NN Classification "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Here we will use the randomly created training and testing dataset. A data_generator function written below creates the required dataset - which can be controlled for repeated dataset using the random seed as a parameter. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "# returns X, y arrays \n",
    "def data_generator(seed=19):\n",
    "    \n",
    "    np.random.seed(seed)\n",
    "    \n",
    "    # Mean and covariance for Class 0\n",
    "    mean0 = [0, 0, 0]\n",
    "    cov0 = [[2550, 2000, 1500], [2000, 1500, 1200], [1500, 1200, 1900]]  \n",
    "\n",
    "\n",
    "    # Number of datapoints for class 0\n",
    "    m0 = 100\n",
    "\n",
    "\n",
    "    # Generate class 0 data points from a multivariate (3D) Gaussian distribution\n",
    "    #    Here x0_1, x0_2 and x0_3 are 3 dimensions for each data (feature) point\n",
    "\n",
    "    x0_1, x0_2, x0_3 = np.random.multivariate_normal(mean0, cov0, m0).T\n",
    "\n",
    "\n",
    "    # Concatenate the 3 dimensions of each feature to create the data matrix for class 0 \n",
    "    X0 = np.concatenate((x0_1.reshape(-1, 1), x0_2.reshape(-1, 1), x0_3.reshape(-1, 1)), axis=1)\n",
    "\n",
    "    # Create the target vector for class 0 (target is coded with zero)\n",
    "    X0_target = np.zeros((m0,), dtype=np.int).reshape(-1, 1)\n",
    "\n",
    "\n",
    "\n",
    "    # Mean and covariance for Class 1\n",
    "    mean1 = [3, 3, 3]\n",
    "    cov1 = [[2550, 2000, 1500], [2000, 1500, 1200], [1500, 1200, 1900]] \n",
    "\n",
    "    # Number of datapoints for class 1\n",
    "    m1 = 100\n",
    "\n",
    "\n",
    "    # Generate class 1 data points from a multivariate (3D) Gaussian distribution\n",
    "    #    Here x1_1, x1_2 and x1_3 are 2 dimensions for each data (feature) point\n",
    "    x1_1, x1_2, x1_3 = np.random.multivariate_normal(mean1, cov1, m1).T\n",
    "\n",
    "    # Concatenate the 3 dimensions of each feature to create the data matrix for class 1\n",
    "    X1 = np.concatenate((x1_1.reshape(-1, 1), x1_2.reshape(-1, 1), x1_3.reshape(-1, 1)), axis=1)\n",
    "\n",
    "    # Create the target vector for class 1 (target is coded with one)\n",
    "    X1_target = np.ones((m1,), dtype=np.int).reshape(-1, 1)\n",
    "\n",
    "\n",
    "    #  Class 0 and 1 data are combined to create a single data matrix X\n",
    "    X = np.append(X0, X1, axis=0)\n",
    "\n",
    "    # Target values for class 0 & 1 are combined to create a single target vector\n",
    "    y = np.concatenate((X0_target, X1_target), axis=0)\n",
    "    \n",
    "    return X, y.ravel()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "X shape: (200, 3) and y.shape (200,) \n"
     ]
    }
   ],
   "source": [
    "# features array and target label array\n",
    "X, y = data_generator()\n",
    "print('X shape: %s and y.shape %s ' %(X.shape, y.shape))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>0</th>\n",
       "      <th>1</th>\n",
       "      <th>2</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>33</th>\n",
       "      <td>89.255183</td>\n",
       "      <td>62.528022</td>\n",
       "      <td>12.270969</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>190</th>\n",
       "      <td>46.260946</td>\n",
       "      <td>20.803942</td>\n",
       "      <td>22.896130</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>87</th>\n",
       "      <td>-31.340537</td>\n",
       "      <td>-33.848504</td>\n",
       "      <td>-13.061759</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>67</th>\n",
       "      <td>-69.299642</td>\n",
       "      <td>-57.668082</td>\n",
       "      <td>-42.714206</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>11</th>\n",
       "      <td>31.631274</td>\n",
       "      <td>23.228821</td>\n",
       "      <td>1.898840</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "             0          1          2\n",
       "33   89.255183  62.528022  12.270969\n",
       "190  46.260946  20.803942  22.896130\n",
       "87  -31.340537 -33.848504 -13.061759\n",
       "67  -69.299642 -57.668082 -42.714206\n",
       "11   31.631274  23.228821   1.898840"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "pd.DataFrame(X).sample(5)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "####  Performance table with Minkowski distance"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Below, a function create_table is written which calls best_performance function five different times to create the required table with best set of parameters as mentioned in the question, which are obtained by the cross-validation with GridSearchCV hyperparamteres tuning. Note that each call of best_performance function returns the performance metrics on a new set of randomly generated dataset. With the option scaling=True, it will implement the scaling of data, default case doesnot."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "# retunrs a single row for table 1:\n",
    "# does cross validation and returns the performance metrics with best parameters\n",
    "def best_performance(scaling=False, seed=None):\n",
    "    \n",
    "    # data generation\n",
    "    X, y = data_generator(seed)\n",
    "    \n",
    "    x_train, x_test, y_train, y_test = train_test_split(X, y, test_size=0.2, random_state=19)\n",
    "    \n",
    "    # if scaling is true\n",
    "    if scaling:\n",
    "        x_train = scale(x_train)\n",
    "        x_test = scale(x_test)\n",
    "    \n",
    "    # GridSearchCV to find best hyperparameters\n",
    "    n_neighbors = [1, 3, 5, 7, 17, 25, 35, 37]\n",
    "    p = [1, 2, 4, 10, 20, 50, 100]\n",
    "    weights = [\"distance\", \"uniform\"]\n",
    "    param_grid = {'weights': weights, 'n_neighbors': n_neighbors, 'p': p}\n",
    "    knn_cv = KNeighborsClassifier(metric='minkowski')\n",
    "    knn_cv = GridSearchCV(knn_cv, param_grid=param_grid, scoring='f1', cv=5)\n",
    "    knn_cv.fit(x_train, y_train)\n",
    "    \n",
    "    # Creating the k-NN model with best hyperparameters\n",
    "    best_params = knn_cv.best_params_\n",
    "    knn = KNeighborsClassifier(n_neighbors=best_params['n_neighbors'], p=best_params['p'], weights=best_params['weights'])\n",
    "    knn.fit(x_train, y_train)\n",
    "    y_pred = knn.predict(x_test)\n",
    "    \n",
    "    # perfromance metrics calculations\n",
    "    accuracy = accuracy_score(y_test, y_pred)\n",
    "    precision = precision_score(y_test, y_pred)\n",
    "    recall = recall_score(y_test, y_pred)\n",
    "    f1 = f1_score(y_test, y_pred)\n",
    "    \n",
    "    # formatting the output: as a single row in the required table\n",
    "    out = dict()    \n",
    "    out['weights'] = best_params['weights']; out['n_neighbors'] = best_params['n_neighbors']; out['p'] = best_params['p']\n",
    "    out['Accuracy'] =  accuracy; out['Precision'] = precision\n",
    "    out['Recall'] = recall; out['F1-score'] = f1\n",
    "    \n",
    "    return out\n",
    "\n",
    "# creates a table, a dataframe\n",
    "# Note that five different seeds are used to get\n",
    "# different set of data for each of the row in the table\n",
    "def create_table(scaling=False):\n",
    "    outs = dict()\n",
    "    \n",
    "    # each time with new seed it creates new dataset inside the function best_performance\n",
    "    for seed in [13, 103, 1123, 1219, 137]:\n",
    "        if scaling:\n",
    "            ans = best_performance(scaling=True, seed=seed)\n",
    "        else:\n",
    "            ans = best_performance(scaling=False, seed=seed)\n",
    "        for key, val in ans.items():\n",
    "            if key not in outs.keys():\n",
    "                outs[key] = [val]\n",
    "            else:\n",
    "                outs[key].append(val)\n",
    "    # creating the final dataframe\n",
    "    df = pd.DataFrame.from_dict(outs)\n",
    "    # df = df.from_dict(outs)\n",
    "    df.index = np.arange(1, 6)\n",
    "    return df"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>weights</th>\n",
       "      <th>n_neighbors</th>\n",
       "      <th>p</th>\n",
       "      <th>Accuracy</th>\n",
       "      <th>Precision</th>\n",
       "      <th>Recall</th>\n",
       "      <th>F1-score</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>uniform</td>\n",
       "      <td>37</td>\n",
       "      <td>20</td>\n",
       "      <td>0.600</td>\n",
       "      <td>0.541667</td>\n",
       "      <td>0.722222</td>\n",
       "      <td>0.619048</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>distance</td>\n",
       "      <td>3</td>\n",
       "      <td>20</td>\n",
       "      <td>0.575</td>\n",
       "      <td>0.529412</td>\n",
       "      <td>0.500000</td>\n",
       "      <td>0.514286</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>uniform</td>\n",
       "      <td>25</td>\n",
       "      <td>4</td>\n",
       "      <td>0.500</td>\n",
       "      <td>0.454545</td>\n",
       "      <td>0.555556</td>\n",
       "      <td>0.500000</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>uniform</td>\n",
       "      <td>17</td>\n",
       "      <td>1</td>\n",
       "      <td>0.525</td>\n",
       "      <td>0.476190</td>\n",
       "      <td>0.555556</td>\n",
       "      <td>0.512821</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>5</th>\n",
       "      <td>distance</td>\n",
       "      <td>3</td>\n",
       "      <td>50</td>\n",
       "      <td>0.450</td>\n",
       "      <td>0.388889</td>\n",
       "      <td>0.388889</td>\n",
       "      <td>0.388889</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "    weights  n_neighbors   p  Accuracy  Precision    Recall  F1-score\n",
       "1   uniform           37  20     0.600   0.541667  0.722222  0.619048\n",
       "2  distance            3  20     0.575   0.529412  0.500000  0.514286\n",
       "3   uniform           25   4     0.500   0.454545  0.555556  0.500000\n",
       "4   uniform           17   1     0.525   0.476190  0.555556  0.512821\n",
       "5  distance            3  50     0.450   0.388889  0.388889  0.388889"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "create_table(scaling=False)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Scaling Data\n",
    "\n",
    "Note that the scaling is done in the same dataset used above in the absence of scaling, done by using fixed random seeds in the function create_table above."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>weights</th>\n",
       "      <th>n_neighbors</th>\n",
       "      <th>p</th>\n",
       "      <th>Accuracy</th>\n",
       "      <th>Precision</th>\n",
       "      <th>Recall</th>\n",
       "      <th>F1-score</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>distance</td>\n",
       "      <td>37</td>\n",
       "      <td>10</td>\n",
       "      <td>0.450</td>\n",
       "      <td>0.423077</td>\n",
       "      <td>0.611111</td>\n",
       "      <td>0.500000</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>distance</td>\n",
       "      <td>3</td>\n",
       "      <td>10</td>\n",
       "      <td>0.675</td>\n",
       "      <td>0.631579</td>\n",
       "      <td>0.666667</td>\n",
       "      <td>0.648649</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>uniform</td>\n",
       "      <td>25</td>\n",
       "      <td>1</td>\n",
       "      <td>0.500</td>\n",
       "      <td>0.458333</td>\n",
       "      <td>0.611111</td>\n",
       "      <td>0.523810</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>uniform</td>\n",
       "      <td>17</td>\n",
       "      <td>1</td>\n",
       "      <td>0.525</td>\n",
       "      <td>0.476190</td>\n",
       "      <td>0.555556</td>\n",
       "      <td>0.512821</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>5</th>\n",
       "      <td>distance</td>\n",
       "      <td>3</td>\n",
       "      <td>2</td>\n",
       "      <td>0.525</td>\n",
       "      <td>0.478261</td>\n",
       "      <td>0.611111</td>\n",
       "      <td>0.536585</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "    weights  n_neighbors   p  Accuracy  Precision    Recall  F1-score\n",
       "1  distance           37  10     0.450   0.423077  0.611111  0.500000\n",
       "2  distance            3  10     0.675   0.631579  0.666667  0.648649\n",
       "3   uniform           25   1     0.500   0.458333  0.611111  0.523810\n",
       "4   uniform           17   1     0.525   0.476190  0.555556  0.512821\n",
       "5  distance            3   2     0.525   0.478261  0.611111  0.536585"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "create_table(scaling=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In this case, the scaling of the features does not warrent the improvements on the performance of the model. Since the dataset is produced randomly, it might not be always helpful to carry out the scaling. It depends upon the dataset, and we sometimes see the improvement after the scaling and while sometimes even decrement in the performance as seen in the above two tables."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "####  Mahalanobis distance metric"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Below we write a classify_with_mahalanobis function which takes a random seed as parameter and returns a single row in the required table. Each call of this function with an unique seed creates a different dataset. The search of the parameters is done by varying the number of nearest neighbors and weights parameters that the KNeighborsClassifier class takes. After brute force search of the best params, the performance measures on the test data set is returned."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "# creates a single row in the required table\n",
    "def classify_with_mahalanobis(seed):\n",
    "    # new dataset\n",
    "    X, y = data_generator(seed)\n",
    "    x_train, x_test, y_train, y_test = train_test_split(X, y, test_size=0.2, random_state=19)\n",
    "    \n",
    "    # greedy searching for best k and weight\n",
    "    minf1 = 0; params = dict()\n",
    "    for k in np.arange(100, 162, 2):\n",
    "        for weight in ['uniform', 'distance']:\n",
    "            knn = KNeighborsClassifier(n_neighbors=k, weights=weight, algorithm='brute', \n",
    "                           metric='mahalanobis', metric_params={'V': np.cov(x_train)})\n",
    "            knn.fit(x_train, y_train)\n",
    "            y_pred = knn.predict(x_test)\n",
    "            f1 = accuracy_score(y_test, y_pred)\n",
    "            if f1 > minf1:\n",
    "                params['f1_score'] = f1; params['k'] = k; params['weights'] = weight\n",
    "                minf1 = f1\n",
    "    \n",
    "    # now let us find metrics with the best params model\n",
    "    knn = KNeighborsClassifier(n_neighbors=params['k'], weights=params['weights'], algorithm='brute', \n",
    "                           metric='mahalanobis', metric_params={'V': np.cov(x_train)})\n",
    "    knn.fit(x_train, y_train)\n",
    "    # prediction\n",
    "    y_pred = knn.predict(x_test)\n",
    "        \n",
    "    # perfromance metrics calculations\n",
    "    accuracy = accuracy_score(y_test, y_pred)\n",
    "    precision = precision_score(y_test, y_pred)\n",
    "    recall = recall_score(y_test, y_pred)\n",
    "    f1 = f1_score(y_test, y_pred)\n",
    "    \n",
    "    # formatting the output: goes as a single row in the required table\n",
    "    out = dict()    \n",
    "    out['weights'] = weight; out['n_neighbors'] = k\n",
    "    out['Accuracy'] =  accuracy; out['Precision'] = precision\n",
    "    out['Recall'] = recall; out['F1-score'] = f1\n",
    "        \n",
    "    return out\n",
    "\n",
    "\n",
    "# creates a table with classifier_with_mahalanobis\n",
    "# creates new dataset for each of the rows\n",
    "# the dataset will be fixed for the fixed seed values used below.\n",
    "def create_table():\n",
    "    outs = dict()\n",
    "    for seed in [202, 44, 10120, 212, 100]:\n",
    "        ans = classify_with_mahalanobis(seed=seed)\n",
    "        for key, val in ans.items():\n",
    "            if key not in outs.keys():\n",
    "                outs[key] = [val]\n",
    "            else:\n",
    "                outs[key].append(val)\n",
    "    # creating the final dataframe\n",
    "    df = pd.DataFrame.from_dict(outs)\n",
    "    # df = df.from_dict(outs)\n",
    "    df.index = np.arange(1, 6)\n",
    "    return df"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>weights</th>\n",
       "      <th>n_neighbors</th>\n",
       "      <th>Accuracy</th>\n",
       "      <th>Precision</th>\n",
       "      <th>Recall</th>\n",
       "      <th>F1-score</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>distance</td>\n",
       "      <td>160</td>\n",
       "      <td>0.750</td>\n",
       "      <td>0.681818</td>\n",
       "      <td>0.833333</td>\n",
       "      <td>0.750000</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>distance</td>\n",
       "      <td>160</td>\n",
       "      <td>0.650</td>\n",
       "      <td>0.590909</td>\n",
       "      <td>0.722222</td>\n",
       "      <td>0.650000</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>distance</td>\n",
       "      <td>160</td>\n",
       "      <td>0.575</td>\n",
       "      <td>0.517241</td>\n",
       "      <td>0.833333</td>\n",
       "      <td>0.638298</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>distance</td>\n",
       "      <td>160</td>\n",
       "      <td>0.675</td>\n",
       "      <td>0.600000</td>\n",
       "      <td>0.833333</td>\n",
       "      <td>0.697674</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>5</th>\n",
       "      <td>distance</td>\n",
       "      <td>160</td>\n",
       "      <td>0.600</td>\n",
       "      <td>0.538462</td>\n",
       "      <td>0.777778</td>\n",
       "      <td>0.636364</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "    weights  n_neighbors  Accuracy  Precision    Recall  F1-score\n",
       "1  distance          160     0.750   0.681818  0.833333  0.750000\n",
       "2  distance          160     0.650   0.590909  0.722222  0.650000\n",
       "3  distance          160     0.575   0.517241  0.833333  0.638298\n",
       "4  distance          160     0.675   0.600000  0.833333  0.697674\n",
       "5  distance          160     0.600   0.538462  0.777778  0.636364"
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# required performance table for part 3.\n",
    "create_table()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Comparision of Minkowski and Mahalanobis Distances\n",
    "\n",
    "From the many experiments done in both the minkowski and mahalanobis distances, we find that mahalanobis distance does better performance. This has been seen from the cross-validation hyperparameters tuning using minkowski distance, using GridSearchCV, and mahalanobis distance using the custom greedy parameters search.\n",
    "\n",
    "As expected, in a randomly created/distributed dataset, mahalanobis distance usually does better. mahalanobis distance is also used to identify the outlier in dataset. However, the problem with Mahalanobis distance is that, it becomes computationally very expensive for larger and high dimensional dataset as it needs to compute the inverse of the covaraince matrix. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Theoretical Analysis of Data Covariance"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Importing data: \n",
    "\n",
    "Here, we use the winequality-white.csv dataset - the same data used in previous notebooks.\n",
    "From the complete dataset as a DataFrame, let us create a feature maxtix X and target array y. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [],
   "source": [
    "df = pd.read_csv('winequality-white.csv', sep=';')\n",
    "\n",
    "y = df['quality']\n",
    "X = df.drop(columns=['quality'], inplace=False).values"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(4898, 11)"
      ]
     },
     "execution_count": 13,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "X.shape # n*d"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Dimension of covarance matrix of X:  (4898, 4898)\n"
     ]
    }
   ],
   "source": [
    "covX = np.cov(X)\n",
    "print('Dimension of covarance matrix of X: ', covX.shape)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The inverse of covaraince matrix of X:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [
    {
     "ename": "LinAlgError",
     "evalue": "Singular matrix",
     "output_type": "error",
     "traceback": [
      "\u001b[1;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[1;31mLinAlgError\u001b[0m                               Traceback (most recent call last)",
      "\u001b[1;32m<ipython-input-17-c71d52b91e07>\u001b[0m in \u001b[0;36m<module>\u001b[1;34m\u001b[0m\n\u001b[1;32m----> 1\u001b[1;33m \u001b[0mlinalg\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0minv\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mcovX\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0m",
      "\u001b[1;32m~\\Anaconda3\\lib\\site-packages\\numpy\\linalg\\linalg.py\u001b[0m in \u001b[0;36minv\u001b[1;34m(a)\u001b[0m\n\u001b[0;32m    530\u001b[0m     \u001b[0msignature\u001b[0m \u001b[1;33m=\u001b[0m \u001b[1;34m'D->D'\u001b[0m \u001b[1;32mif\u001b[0m \u001b[0misComplexType\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mt\u001b[0m\u001b[1;33m)\u001b[0m \u001b[1;32melse\u001b[0m \u001b[1;34m'd->d'\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m    531\u001b[0m     \u001b[0mextobj\u001b[0m \u001b[1;33m=\u001b[0m \u001b[0mget_linalg_error_extobj\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0m_raise_linalgerror_singular\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[1;32m--> 532\u001b[1;33m     \u001b[0mainv\u001b[0m \u001b[1;33m=\u001b[0m \u001b[0m_umath_linalg\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0minv\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0ma\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0msignature\u001b[0m\u001b[1;33m=\u001b[0m\u001b[0msignature\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mextobj\u001b[0m\u001b[1;33m=\u001b[0m\u001b[0mextobj\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0m\u001b[0;32m    533\u001b[0m     \u001b[1;32mreturn\u001b[0m \u001b[0mwrap\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mainv\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mastype\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mresult_t\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mcopy\u001b[0m\u001b[1;33m=\u001b[0m\u001b[1;32mFalse\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m    534\u001b[0m \u001b[1;33m\u001b[0m\u001b[0m\n",
      "\u001b[1;32m~\\Anaconda3\\lib\\site-packages\\numpy\\linalg\\linalg.py\u001b[0m in \u001b[0;36m_raise_linalgerror_singular\u001b[1;34m(err, flag)\u001b[0m\n\u001b[0;32m     87\u001b[0m \u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m     88\u001b[0m \u001b[1;32mdef\u001b[0m \u001b[0m_raise_linalgerror_singular\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0merr\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mflag\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m:\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[1;32m---> 89\u001b[1;33m     \u001b[1;32mraise\u001b[0m \u001b[0mLinAlgError\u001b[0m\u001b[1;33m(\u001b[0m\u001b[1;34m\"Singular matrix\"\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0m\u001b[0;32m     90\u001b[0m \u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m     91\u001b[0m \u001b[1;32mdef\u001b[0m \u001b[0m_raise_linalgerror_nonposdef\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0merr\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mflag\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m:\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n",
      "\u001b[1;31mLinAlgError\u001b[0m: Singular matrix"
     ]
    }
   ],
   "source": [
    "linalg.inv(covX)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Here, we are getting the covX as a singular matrix because the numpy function incorrectly calculates its determinant 0. This is not correct as we check from the tests done for a real singular matrix, done **below**  with the $2\\times2$ matrix given in the updated recitation note."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [],
   "source": [
    "covxx = pd.DataFrame(covX)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [],
   "source": [
    "v, _ = linalg.eig(covxx.values)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {},
   "outputs": [],
   "source": [
    "for val in v:\n",
    "    if abs(val) == 0.0:\n",
    "        print(val)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Here, it turns out that non of the eigenvalues of the covaraince matrix is 0."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.0"
      ]
     },
     "execution_count": 23,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "linalg.det(covX)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {},
   "outputs": [],
   "source": [
    "eigvalues, eigvectors = linalg.eig(covxx)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {},
   "outputs": [],
   "source": [
    "# checking if any of the eigenvectors are zeros or not\n",
    "for i in range(len(eigvectors)):\n",
    "    for j in range(len(eigvectors[i])):\n",
    "        if eigvectors[i][j] == 0.:\n",
    "            print(0.)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**No eigenvectors is zero: so it might not be singular, as we are getting the det 0 because of some rounding error.**"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Positive-definite or semi-definite?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "xTAx:  [400]\n",
      "===========\n",
      "checking: \n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "array([ True])"
      ]
     },
     "execution_count": 26,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "A = [[1, 2], [2, 4]]\n",
    "\n",
    "# a random vector of length len(A)\n",
    "x = np.random.randint(0, 10, (len(A)))\n",
    "xT = x.reshape(len(A), 1)\n",
    "\n",
    "xTAx = np.matmul(np.matmul(x, A), xT)\n",
    "print('xTAx: ', xTAx)\n",
    "\n",
    "print('===========')\n",
    "print('checking: ')\n",
    "xTAx >= 0"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "metadata": {},
   "outputs": [],
   "source": [
    "# with above covariance matrix\n",
    "A = covX"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 28,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.0"
      ]
     },
     "execution_count": 28,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "linalg.det(A)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 29,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "det[D1]:  2550.334776454546\n",
      "det[D1]: 2550.334776, det[D2]: 143612.728148, det[D3]: 1080266.369765, det[D4]: 5064071.143343, det[D5]: 0.000000, det[D6]: 0.000000, det[D7]: 0.000000, det[D8]: 0.000000, det[D9]: 0.000000, det[D10]: 0.000000, det[D11]: 0.000000, det[D12]: 0.000000, det[D13]: 0.000000, det[D14]: 0.000000, det[D15]: 0.000000, det[D16]: 0.000000, det[D17]: 0.000000, det[D18]: 0.000000, det[D19]: 0.000000, det[D20]: 0.000000, det[D21]: 0.000000, det[D22]: 0.000000, det[D23]: 0.000000, det[D24]: 0.000000, det[D25]: 0.000000, det[D26]: -0.000000, det[D27]: 0.000000, det[D28]: 0.000000, det[D29]: -0.000000, "
     ]
    }
   ],
   "source": [
    "# checking the smaller sub-matrices of the large A\n",
    "for i in range(1, 30):\n",
    "    Di = np.array([list(row[:i]) for row in A[:i]])\n",
    "    if i == 1:\n",
    "        D1 = Di[0][0]\n",
    "        print('det[D1]: ', D1)\n",
    "    detDi = linalg.det(Di)\n",
    "    print('det[D%d]: %f' % (i, detDi), sep=', ', flush=True, end=', ')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Notice that the determinants of all the sub-matrices that are larger than size 4 are 0.\n",
    "This is surprising: let us explore what's wrong with it. It seems there is something wrong in the \n",
    "output of the calculation coming from the python itself."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 30,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "5064071.143342724"
      ]
     },
     "execution_count": 30,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "A4 = np.array([list(row[:4]) for row in A[:4]])\n",
    "linalg.det(A4)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 31,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.0"
      ]
     },
     "execution_count": 31,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "A5 = np.array([list(row[:5]) for row in A[:5]])\n",
    "linalg.det(A5)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "What about the cov(X) matrix?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 32,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "xTAx:  [8.29280492e+11]\n",
      "===========\n",
      "checking: \n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "array([ True])"
      ]
     },
     "execution_count": 32,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# a random vector of length len(A)\n",
    "x = np.random.randint(0, 10, (len(A)))\n",
    "xT = x.reshape(4898, 1)\n",
    "\n",
    "xTAx = np.matmul(np.matmul(x, A), xT)\n",
    "print('xTAx: ', xTAx)\n",
    "\n",
    "print('===========')\n",
    "print('checking: ')\n",
    "xTAx >= 0"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "***Yes,  A (Cov(X)) is positive-semidefinite.***"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "####  Removing singularity of data matrix:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Introduction of small noise into the dataset can remove the singlularity of the matrix slightly changing the determinant of the matrix.\n",
    "\n",
    "But the better way of solving this problem is to decrease the dimnesion of the dataset, i.e. to get rid of some of the redundant features in the dataset so that it becomes non-singular. It can be done using (Principle Component Analysis) PCA. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "####  Time compexity of inverting a matrix:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The airthmetic time compelxity of inversion of a matrix of size ($n \\times n$) using Gauss-Jordan elimination method is **O**(n$^3$). There are slightly better algorithms for the matrix inversion."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Checking out with a simple singular $2 \\times 2$ matrix"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 33,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[1, 2],\n",
       "       [2, 4]])"
      ]
     },
     "execution_count": 33,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "A = [[1, 2], [2, 4]]\n",
    "np.array(A)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 34,
   "metadata": {},
   "outputs": [],
   "source": [
    "# eigenvalues and eigenvectors\n",
    "\n",
    "eigvalues, eigvectors = linalg.eig(A)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 35,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([0., 5.])"
      ]
     },
     "execution_count": 35,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "eigvalues"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 36,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[-0.89442719, -0.4472136 ],\n",
       "       [ 0.4472136 , -0.89442719]])"
      ]
     },
     "execution_count": 36,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "eigvectors"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 37,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.0"
      ]
     },
     "execution_count": 37,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "linalg.det(A)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Here, out of the two eigenvalues, one of them is zero. One thing we can learn from this is that a real singular matrix will have at least one of its eigenvalue zero."
   ]
  }
 ],
 "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.7"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}