Assignment_5.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np;\n",
    "from numpy import linalg as la;\n",
    "import matplotlib.pyplot as plt;\n",
    "from sklearn.datasets import load_wine;\n",
    "samples, labels = load_wine(return_X_y=True);\n",
    "samples = np.transpose(samples);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Standardization of Data"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### PCA on Centered Data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "samplesmean = np.repeat(np.reshape(np.asmatrix(samples.mean(1)),(13,1)),178,axis=1);\n",
    "samples1 = samples - samplesmean;"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "There are the two-dimensional features of the centered dataset.\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXwAAAD8CAYAAAB0IB+mAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAIABJREFUeJzt3X2UVPWd5/H3l0YajI0P24geEdFZzRlMDEqLGoTNBCcxjivzcNjjhsw6JqO7huwkyjxo2FVnk5xmEs/qmVnRYR5ykjPxRFF2dXLiovEYkES0G2yNympQgxh56szOwNotCPXbP+6tprr6VtWtqvtY9/M6pw9Vt27V/dbt5nt/9/dozjlERKTzTUo7ABERSYYSvohIQSjhi4gUhBK+iEhBKOGLiBSEEr6ISEEo4YuIFIQSvohIQSjhi4gUxOS0A6jU29vr5syZk3YYIiK5snXr1mHn3IxG+2Uq4c+ZM4fBwcG0wxARyRUz2xlmP1XpiIgUhBK+iEhBKOGLiBSEEr6ISEEo4YuIFIQSvohIQSjhi4gUhBK+iERueGSY+WvnMzwynHYoUkEJX0Qi17+5n6E9Q6zevDrtUKSCEr6IRGp4ZJj7Bu6j5EqsGVijUn6GKOGLSKT6N/dTciUASq6kUn6GZGouHRHJv007N2FmTJ08FYCNOzemHJGUKeGLSKQGbhhIOwSpQVU6IiIFoYQvIlIQSvgiIgWhhC8iUhBK+CIiBaGELyJSEEr4IiIFoYQvIlIQSvgiIgWhhF9U06eD2bGf6dPTjkhEYtZ2wjezM83saTPbbmavmNmX/e2nmNmTZvZz/9+T2w9XInPwYP3nItJxoijhHwFWOud+HbgUWGFmc4Fbgaecc+cCT/nPRUQkJW0nfOfcbufcNv/xQWA7cAawFPiOv9t3gN9u91gikm1a6SrbIq3DN7M5wIXAc8BM59xu8C4KwKk13nOjmQ2a2eD+/fujDEfq6emp/1ykBVrpKtsiS/hmdgLwCPAV59yBsO9zzq11zvU55/pmzJgRVTjSyIED4NyxnwOhf2UigbTSVfZFkvDN7Di8ZP8959x6f/NeMzvdf/10YF8UxxKRbNJKV9kXRS8dA/4O2O6c++8VLz0GXOc/vg54tN1jiUh2Va50ZWZa6SqDoljxaiHw+8DPzGzI3/ZVYDXwkJl9AXgbWBbBsUQko5Jc6Wp4ZJhP/8OneeB3H+Cz6z/Lhs9toPf43sSOn1dtJ3zn3GbAary8pN3PFxGpVm4cXr5++Vgj8V2fuivtsDJPI21FJBWtduGsbBzeunurGomboIQv2aHpHgql1S6clY3DZWokDkcJX7JD0z0URjtdODft3DTuuWFqJA5JCV+CqbQtMWqnC+fADQOsWLCCSTaJlZetpHRHiV0376LkSqrWaUAJX4IFlbZbvQBk6eKRpVgKrJ0unEF3BxrhG44559KOYUxfX58bHBxMOwwBLxk2EvZvJ+izgt4bdr92JHGMDlTuBpmF7o8rn1jJvc/fy6Gjh+ju6uYLF36Bbw99m9Ejo0ybPI23b3479RiTZmZbnXN9jfZTCV9EGspSCbr67uCR7Y+Mqx664+k7mu79s3cvLF7s3fAtXuw970Qq4ctE06eHazBVCb8QhkeGmX33bEaPjGIY21ds58O9H047rDEX/83FvLzv5bHnJ3afyP6R/dx86c2h++YvXgzPPgtHjsDkyXDZZbBpU+P3ZYVK+NK6MMm+mdk1w87MmcQMnpolNFC9PvGVDawOx/L1y4HslIoHbhhgdNUoo6tG2XXzLg4cOtB075+hIS/Zg/fv0FDwfnmf/lkJX8Lr6Wltds2wM3MmMYOnZgkNVK/Kprob5Lbd2xgeGWbZMq9UfPCg9++yDEye0mrvn3nzvJI9eP/Om1f787NStdUKJfxWFLGnh5Jjx2rUJ77cDbK7qxuAKV1TWL15dehScZJa7f2zbp1XjdPT4/27bt3EfTph+ucoJk8rnk4fINTTM/47qdqjowWViqvrvisTKcDGnRuZN298vXetUnGSWp3AbebMxnX2Yc5T1inhy0QqyefS3r1etcrQkJd8163zElkjQcm8WlAi3XvNxOO1GkPWDY8Mc//g/QB1z1PWqZdOK9TTQzIoCz1N0owhzrECK59YyT1b7mmq50+S1EsnTurpIRnUSp16vZ429V6r1VsliXr9WseOq0G1E+ruy5TwW6GeHpJBYXuaVKrX06bea7WSaysxNCvo2HEm5U5aulEJX4qtg3pchelpUq1eibzWa/WSaysxlIXp417r2HEm5U5aulEJX4qlOsF3UI+rck+TAwe8f8M0ltYrkdd6rV5ybSWGsjBVMrWO/dSbT3H46GG6u7qbSsphLjKPL3+cuTPmsuvmXYyuGk10KceoKeFLseQ4ocehXom81mtRlHir2wdefetXoapkah17yTlLMDO+tOBLjK4a5QfXDIQaBVyreqjyIpD3wVbjOOcy8zN//nzXtJ6eytp073nashhTs7L8HVqJrfo99X6kbXv2OLdokXfaFy3ynldatMi5yZO90z15snNnfGSHm/Tnkxx34ib9+SS3csPK0Mfa/95+N+3r0xx34qZ9fZrb/97+CZ+/aJG330V/fZHb/97+mu9zzrlbNtwyFkPlPlO/NtVdsOaCsf2yBBh0IXJs/kv4Wbwlz2JMzcryd2gltrDxq8dVJBpNu1DdPrB7x8xxVTVPvfXUhM98+WWYPr2EmWP69BIv+/OlBVXzvDBUmtD+UF1Sr37f7U/fzvn3ns/dz949dqdx+9O3j+1z+OhhXtr3Uq5L+vlP+FmT40a/QqqcH0g9riJTndB/8pPx1SrV7QMzf233uPd/fNbHJ3zmwoVw8KABxsGDxsKF3vagap4Tz3oLJn0w9vnnf/SDCVVG1e9bv309rw6/isMbU1NyJR5+9WEOHz3MlElTKOEl/jx3zVTCj1qWSsJR64CeLOOUp5DotO8VszCzZFZ3xyyVxpfyq9sHPvi93xm3/yPbH5nwmQcOOKA86NH85+NnyxxdNcrjyx/nV1cvgVnPwpQD9C34gI/90dcn3AVUvm/oPw6x7719445nZkyySZgZJ0096dh3yXHXzPwn/DwMgspiTI3Ui7nVi1pUXSCj+p1nudoqw8LMkrluHXR1jd9W2eWzujfPnFnTmDp56tjPmSeeOeEzpxz/PlAe0e785xP1b+7HfWgPfP7f0P1fT2XhHbex9cAP6zY0f3b9Z8dK9gCTmMT1864fm2p538ixi0Hl+/M2XXL+E34eBkFlMaZGKs9rVKJKsHn4nedAq/PZV1fXPPPMxPfPnAkf/3j4QVjVpfSgro//+k/+ALoPAA66D3jPA77T2hXLOfS1fdi3n4H3ZrJx58bAzy8n69eGX+OF3S+M+xyHY/329WN3BWXdXd2suHjFWHx568GT/4SfNXm44yiavFXVJDQYrNX57Cvr38uC3t/OIKwgr9z5IO79E3HOcO+fyCt3Pjhhn2XL4P1fXASHp9P1y8tZ8JOdNfvNl5P18vXLmdI1BfAS+k19N3Hh6Rdy2gmnYWaYX41k2ITSfbuje5O+Q1DCj1onlj6zehELmxiD7iR6erL7vRKqamp13ptyIq8U9P52BmE1q5w4g3rn1Nq/nKy37d4GMK7xdmjPEFeccwWjq0Yp3VHC3eEo3VEad/cRxejepO8QlPClsaguYlEn2HYS44EDqV2c01wasPLYZsfq2SurXILiqyyJlhP5okXxz5sTVjlxnnjWL0KvXFVO1lO6powN2GpmicR2B6ClMSmbEr4kpxPvflqQ5tKAlcceGYEPfWhilUtQfEEl0airbMqavSBWJs5fXf1J+hZ80DCmWsm6mVJ7UO+g8kUxTFVNGpOyaQEU6XzVK3jBsTUNenoSv/A0rEqJccWx6mMHXXer93lhqMRgRUn01stvpff43sBVoqKYk758wTly5NgFp96c+uMWWf/QHhbecRvPNpizvla9fpjFYOrFUb4oOtzY46D589NaUEUlfMmvWlVE1XX7ULvHUQpdMRtOIRzjnVCY6Yur9znxrLdCl0SjqJNutm0hytksw/QWClJdPbPm+TV1q2r6N/fz/pH3x6qSkpqUTQlf8qtWYoyq0bPV3jIN3hdXVUgYYY5dvU/vH3wxVEKNqk662Tn1W03SUaq8yzh89DCHjx4Ggi+QaS6ooiodKYbp05svKbd64WjwvjALZsclzLEn7rMh1GdHtcj3unUT18VtVZzLHlaqvMs4dOQQULuqJs3F0JXwpRjKSTfG+vG8azc5tlP/XSnKC2JlFVOcSbXWXUX5nA6PDI+d06jOUytUpSOdp14SV0+hmtqtfy9Xrey6eRdzZ8zl8eWPRxxhcwOVsrAWbdA5TbMKSglfguV56b+oknir4wYavC+L869EmRyDklxUYw+auSilvRZtFi441WJP+GZ2pZm9ZmY7zOzWuI8nEQmqh87TRSCKQV6t3g00eF8W519plBzDXqRqJbkoxh40m0DTXos27QtOkFgTvpl1AfcCnwHmAv/ezObGeczcykMyzdPsku1W3cT0+8hiqQ8aJ8f/8o9/xbZv3M0ZvSc0XDIwKMk16moZ5g6g2QSadu+dtC84QeJutF0A7HDOvQlgZt8HlgKvxnzc7Jo+fWKj4YED+UqmRRDT7yPNHhr11EuGwyPD/M2f/SbsuoTD7jh++lPHeecZzh3rRVOeJ6dWg+S8eccGUwV1tQwz2CrNxs5WlM9pueE2jjaNZsVdpXMGsKvi+Tv+tuLKS2LP6sRiOddsqS/NeXfK+jf3U9p9AbjjADh61MbKKNXVM7VK1Y36/4cZbDVwwwC/+MNRPvaPoxy+fZQXbhrg0kvTOSfNyFQVXpiFb1v9AZYBf1vx/PeBv6ra50ZgEBicPXt2nOv8ZkOtRbPzsJh2lhc2j1pCv4/qhbWrBS3GnbS+tX3OznrGMemwfypKkf8ZhP2eixZN/LWkcU7CqrVQetTIyCLm7wCVS9fMAt6tuuCsdc71Oef6ZsyYEXM4GZaHEnWRujQm9PtoVPprdQrjKA3cMMDu5y5n0cLj6OmB6dMt8lkyw4wA3rvXWxu3WhrnJKysNdzGnfAHgHPN7GwzmwJcCzwW8zGzrVYiKVIyzYMEfh9hGnCbnWYAiKXBuXJu+9dfT2dqiGXLvLVxq6U5LXMjmWu4DXMb0M4PcBXwOvAGsKrevvPnz4/ldkfqKFI1TcbcsuEW1/21bseduO6vdbuVG1ZO2GfPHq/KoqfH+3fPnhAfXFXnsYdTm/+MCo2qnaIQpkqn+k8VnLvkkua/TyciI1U6OOd+6Jw7zzn3a865b8R9PGlSXhqRO1CY0l8Uq0YtY11bfeCTaHQMU3VVfbezaBFs2dLeSlpZHAQXJ420lWTkYZxBwpLqJz7EvJbbAZIaNxCm6iqOWUYz1YMmAUr4kgzdSSSnqp1o3qSftdzImlSjY5hkHuUaucMjw3zsvo81nLe+0yjhd4pWS9B56B1US5bvGtKMrarBed27C1suGUfd6FhrXEGSC56DdyF7ad9Ldeet70Tm1fdnQ19fnxscHEw7jHwqr+xUKUO/21jiy/J3rhVbrZHWWRVxvIsXjx9xe9llya0NUB7x+sDvPsCFf30ho0dGAeju6sbM+MipH0ll8ZQomNlW51xfo/1Uwpdk5PlOIkp5q9oKEW8zo4HTHFdQrq9fvn75WDVVd1f32DKDeU32zVDCl2QkNc4gi9U7KYt7eoZmZsJsaVxBBCobn7ft3hZ53/i89PZRwu8URSxB1/uOaZecM/T7iGJq4noaldork2Fa6/n2b+4fq6+fZJNYcfGKSHtH5aW3jxJ+pyjiSN3K75w1tX4fKVwI2qpGCRFvo1J7ZTJMunG2fLF5YscTHHVHATjqjvLUW09FeowsTnkdRAlfJEkpXJjbqkYJEW+9UnvaybB8seme3E13Vzfg1dsvOXtJpMfI0nw59SjhS2fIUBVK1sRdjVKv1J5mMqyutwdimdMmc/Pl1BH3AigiychTFVbCXTPLCbkde/d6df9DQxMXPaknzUVLKi82U7qm8KUFX4plsZk89e5RwhdJWt66ZhJuRaogUSTDPF5sskoDryS78jZIKawsDxirUk62zzwzfnuSv4o0B2vlhQZeSf7lsCTcacol+0qt9J9vZyxAFhaBqZaXfvfVlPBFkpajBubKZFvWSsNvrbEAYS4EaQ3Wqicv/e6rKeFLvnTCKNocjZkImoO+lf7ztUrpYQaFpTVYq5a0u5q2QwlfsqtWybdDq3binAKh1SqIqJJtrVJ69YXgJwP/b0KMSQ/WaiRP/e6rKeFLbWlPP1wuCRdEnFMgtFoFEVWyrXXhqLwQWNdRSqduy3wCzVO/+2rqpSO1ZaU3SVbiiFlcnZKGR4aZffdsRo+MMm3yNN6++W16j+9t/4MjUO4F9MJQiZFTfkpp2e8x7aSDmYqxnvKUyxs+tyHVeNVLR1pTWarPihw1crYjrsbJLFdBlO8gbnz4TzjuD6+AE/ZlLsZ68tZ4q4Qv42WxfjxHjZztiKtxMg9VEHmIsVoeG2810lbCy1PJOoeDtqKYAiFIs6NdWx3Z2o48TU9QFnTnFMfUDVFSCV8ay2PJWoO2Whb3/PmdIo93JSrhy3g9PRNLxlIoWRzZmkV5vCtRCV/GK0h9udSWxZGtEg0lfOlMUffsSXJMQsrjH7I2slWio4QvnSnqO5Uk2wRSbn9odbBV3IulS/uU8EUkEmrszT4lfBGJRJEbe/MyXbISvkgYSY72zenI4iI39uZlxK0SvkgYSfZeauZYaU9wVyGNxt4slKzzNOJWCV8kzzI0wCyqmTWbafytLlmn0XCc5bmKqinhi0imhG38DSpZp9FwnKcRtxppKyKZErbxN6hkPTR0V+INx3kacasSvkie5bSBt56wjb9BJesiNxyHoRK+SAbt3QtLl8LgoNdue/HF8OijAfXirTQeZ3wm0XXrJs7WGSSoZL33mnDvLaq2Vrwys28B/xY4DLwBXO+c+2f/tduALwBHgT9yzm1o9Hla8UrEs3gxPPPM+G3lBcTbVpAVxIokqRWvngQ+4py7AHgduM0/+FzgWuB84EpgjZl1tXkskcIIqnsu0kAmiUdbCd8594Rzzm8iYQswy3+8FPi+c+6Qc+4tYAewoJ1jiRRJUN2z6qOlXVE22n4eeNx/fAawq+K1d/xtIhLCunVwySXQ1QWTJnmPI6uP7sCGXgmnYaOtmf0IOC3gpVXOuUf9fVYBR4Dvld8WsH9gJaGZ3QjcCDB79uwQIYt0vpkzYcuWmD48Qw20kqyGCd85d0W9183sOuBqYIk71gL8DnBmxW6zgHdrfP5aYC14jbYhYhapL+O9UETS0laVjpldCfwZcI1zbqTipceAa82s28zOBs4Fnm/nWCKhZWi6AZEsabcO/38APcCTZjZkZvcDOOdeAR4CXgX+N7DCOXe0zWOJiESqaIu2tNUPP2rqhy+R6KR+5qqeitXixd6cO0eOeCNzL7ssorEOCUuqH75I9nRSLxRVT8WqaIu2KOFLezI0H/uYJOeul1wr2tw7SvjSHpVAJcfSWLQlTUr4IlmWo+qpLDaANloRK6pFW/JCCV8ky3JUPZXG4iON5GWt2aQo4Ut7clQClXhlrQE0T2vNJkUJX9qToxJox0u5AT1rDaB5Wms2KUr4Ip0i5Qb0rDWA5mmt2aRoxSsRiUS5ATQr8rTWbFJUwhcRKQglfJFOoQZ0aUAJX6RTVDegQ/ZGQUcoi/3+s04JX6RTdfgo6Cz2+886JXwRyaWs9fvPAyV8EcmlrPX7zwMlfJFO1eGNuFnr958H6ocv0qk6fNRz1vr954FK+CIiBaGELyJSEEr4IiIFoYQvIlIQSvgiIgWhhC8iUhBK+CIiBaGELyJSEEr4IiIFoYQvIlIQSvgiNQyPDDN/7XyGR4bTDkUkEkr4IjX0b+5naM8QqzevTjsUkUgo4YsEGB4Z5r6B+yi5EmsG1uS+lK+7FQElfJFA/Zv7KbkSACVXyn0pP8m7FV1csksJXyTApp2bMDOmTp6KmbFx58a0Q2pZ0ncrqgrLLs2HLxJg4IaBtEOITNDdyl2fuiuWY1VfXG69/FZ6j++N5VjSPJXwRTpckncrnVYV1mlUwhfpcEnerVReXIBcV4V1IiV8EYlMJ1WFdSJV6YiIFEQkCd/M/tjMnJn1+s/NzP7SzHaY2UtmdlEUxxERkda1nfDN7EzgN4G3KzZ/BjjX/7kRuK/d44iISHuiKOHfDfwp4Cq2LQW+6zxbgJPM7PQIjiUiIi1qK+Gb2TXAL51zL1a9dAawq+L5O/62oM+40cwGzWxw//797YQjIiJ1NOylY2Y/Ak4LeGkV8FXgU0FvC9jmArbhnFsLrAXo6+sL3EdERNrXMOE7564I2m5mHwXOBl40M4BZwDYzW4BXoj+zYvdZwLttRysiIi1ruUrHOfcz59ypzrk5zrk5eEn+IufcHuAx4D/4vXUuBf7FObc7mpBFRKQVcQ28+iFwFbADGAGuj+k4IiISUmQJ3y/llx87YEVUny0iIu3TSFsRkYJQwhcRKQglfBGRglDCFxEpCCV8EZGCUMIXESkIJXwRkYJQwhcRKQglfBGRglDCFxEpCCV8EZGCUMIXESkIJXwRkYJQwhcRKQglfBGRglDCFxEpCCV8EZGCUMIXESkIJXwRkYJQwhcpuOGRYeavnc/wyHDaoUjMlPBFCq5/cz9De4ZYvXl12qFIzJTwRQpseGSY+wbuo+RKrBlYo1J+h1PCFymw/s39lFwJgJIrqZTf4SanHYCIpGfTzk2YGVMnTwVg486NKUckcVLCFymwgRsG0g5BEqQqHRGRglDCFxEpCCV8EZGCUMIXESkIJXwRkYJQwhcRKQglfBGRglDCFxEpCCV8EZGCUMIXESkIJXyRCGhOecmDthO+mf1nM3vNzF4xs29WbL/NzHb4r3263eOIZJnmlJc8aCvhm9lvAEuBC5xz5wN3+dvnAtcC5wNXAmvMrKvNWEUySXPKS160W8K/CVjtnDsE4Jzb529fCnzfOXfIOfcWsANY0OaxRDJJc8pLXrSb8M8DFpnZc2a20cwu9refAeyq2O8df5tIx6mcU97MNKe8ZFbD+fDN7EfAaQEvrfLffzJwKXAx8JCZnQNYwP6uxuffCNwIMHv27HBRi2SI5pSXvGiY8J1zV9R6zcxuAtY75xzwvJmVgF68Ev2ZFbvOAt6t8flrgbUAfX19gRcFERFpX7tVOv8L+CSAmZ0HTAGGgceAa82s28zOBs4Fnm/zWCIi0oZ2lzj8e+Dvzexl4DBwnV/af8XMHgJeBY4AK5xzR9s8loiItKGthO+cOwx8rsZr3wC+0c7ni4hIdDTSVkSkIJTwRUQKQglfRKQgzGtjzQYz2w+8h9fTJw96UaxxyVO8ijUeijW8s5xzMxrtlKmED2Bmg865vrTjCEOxxidP8SrWeCjW6KlKR0SkIJTwRUQKIosJf23aATRBscYnT/Eq1ngo1ohlrg5fRETikcUSvoiIxCDRhG9my/ylEEtm1lf1WuCSiGZ2pb9th5ndWrH9bH8e/p+b2YNmNiXm2OeZ2RYzGzKzQTNb4G83M/tLP76XzOyiivdc58f3czO7Ls74AuJtaunJWuc5wXj/2MycmfX6zzN3Xs3sW2b2f/x4/qeZnVTxWibPa6UsxeLHc6aZPW1m2/2/0y/7208xsyf93++TZnayv73m30SCMXeZ2Qtm9gP/eWAe8ieOfNCP9Tkzm5N0rIGcc4n9AL8OfBj4MdBXsX0u8CLQDZwNvAF0+T9vAOfgzcT5IjDXf89DwLX+4/uBm2KO/QngM/7jq4AfVzx+HG8NgEuB5/ztpwBv+v+e7D8+OaHz/BvAj4Bu//mprZ7nhOI9E9gA7AR6M3xePwVM9h//BfAXWT6vVbFnJpaKmE4HLvIf9wCv++fym8Ct/vZbK85z4N9EwjHfAjwA/MB/HpiHgC8C9/uPrwUeTPNcl38SLeE757Y7514LeKnWkogLgB3OuTedN1Hb94GlZmZ40zI/7L//O8Bvxx0+MN1/fCLH5vdfCnzXebYAJ5nZ6cCngSedc//knPu/wJN46/smodmlJwPPc0KxAtwN/CnjF8nJ3Hl1zj3hnDviP92Ct85DOdYsntdKWYoFAOfcbufcNv/xQWA73sp4S/H+T8P4/9u1/iYSYWazgN8C/tZ/Xi8PVX6Hh4El/v6pykodfq0lEWtt/1fAP1f850tiCcWvAN8ys114i7Xf5m9vNvYkNLv0ZGqxmtk1wC+dcy9WvZS5WKt8Hq+0CdmPlYzFMoFf5XEh8Bww0zm3G7yLAnCqv1va3+EevIJJyX9eLw+Nxeq//i/+/qlqdz78CazOkojOuUdrvS1gmyP4guTq7N+WerEDS4CbnXOPmNm/A/4OuKJOLLHEGDLWZpeerHWeI9Eg1q/iVZVMeFuNmFI7r+W/XzNbhbfOw/fKb6sRU6zntUmxnrd2mNkJwCPAV5xzB+oUhFP7DmZ2NbDPObfVzD4RIp5Mnu/IE76rsyRiHfWWRAzaPox3OzfZv3rWXEKxGfViN7PvAl/2n67Dv62rE/s7wCeqtv+43RjLGsTaytKToZakjDJWM/soXp33i/5/8lnANr9BPHPn1Y/5OuBqYIl/fqkTK3W2Jy30sqNJMrPj8JL995xz6/3Ne83sdOfcbr/KplwlmeZ3WAhcY2ZXAVPxqnfvoXYeKsf6jplNxqsG/qeEYq0tjYYDJjbans/4Rq838RqZJvuPz+ZYQ9P5/nvWMb6x5Isxx7wd+IT/eAmw1X/8W4xvSHre334K8BZeSftk//EpCZ3f/wT8N//xeXi3ltbKeU747+IXHGu0zeJ5vRJvFbcZVdszfV79GDMTS0VMBnwXuKdq+7cY32j7zXp/EynE/QmONdoG5iFgBeMbbR9K81yPxZ7wifodvCvfIWAvsKHitVV4vQhew+8N42+/Cq/1/g282+ry9nPw1snd4Z/07phjvxzY6v9HeQ6YX/FHe68f388YfyH7vB/fDuD6BM/zFOAfgJeBbcAnWz3PCf99VCb8LJ7XHXgXzyH/5/48nNcsxuLHczleNcdLFef0KrwfZXCbAAAAW0lEQVS67qeAn/v/ntLobyLhuCsTfmAewrsLWOdvfx44J+3z7ZzTSFsRkaLISi8dERGJmRK+iEhBKOGLiBSEEr6ISEEo4YuIFIQSvohIQSjhi4gUhBK+iEhB/H9rnLc4vJfvUAAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "U1,S1,Vt1 = la.svd(samples1);\n",
    "Us1 = U1[:,0:2];\n",
    "feat1 = np.matmul(np.transpose(Us1),samples1);\n",
    "feat10 = feat1[np.ix_((0,1),np.where(labels==0)[0])];\n",
    "feat11 = feat1[np.ix_((0,1),np.where(labels==1)[0])];\n",
    "feat12 = feat1[np.ix_((0,1),np.where(labels==2)[0])];\n",
    "print(\"There are the two-dimensional features of the centered dataset.\")\n",
    "plt.scatter(feat10[0,:].tolist(),feat10[1,:].tolist(),15,'r','s');\n",
    "plt.scatter(feat11[0,:].tolist(),feat11[1,:].tolist(),15,'g','^');\n",
    "plt.scatter(feat12[0,:].tolist(),feat12[1,:].tolist(),15,'b','o');"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### PCA on Standardized Data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[17.57782045 16.30549754 16.25608694 ... 16.39196609 16.26843959\n",
      "  17.45429396]\n",
      " [ 1.53500398  1.5978404   2.11848502 ...  3.84199826  2.32494755\n",
      "   3.6804189 ]\n",
      " [ 8.88247841  7.82242955  9.75976023 ...  8.26107046  8.66315796\n",
      "  10.0156341 ]\n",
      " ...\n",
      " [ 4.56283226  4.60670565  4.51895888 ...  2.58852984  2.63240323\n",
      "   2.67627662]\n",
      " [ 5.53677584  4.80230557  4.47744372 ...  2.20341079  2.28815736\n",
      "   2.2599085 ]\n",
      " [ 3.39148583  3.34371843  3.77362508 ...  2.65905227  2.67497474\n",
      "   1.78331649]]\n"
     ]
    }
   ],
   "source": [
    "samplesstd = np.repeat(np.reshape(np.asmatrix(samples1.std(1)),(13,1)),178,axis=1);\n",
    "samples2 = samples / samplesstd;\n",
    "print(samples2);"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "There are the two-dimensional features of the centered and normalized dataset.\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXsAAAD8CAYAAACW/ATfAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAHddJREFUeJzt3X+QXWV5B/DvNxuyQdiVTrMNaliDLTBAdJJmG0WaTIdfYVoK9QczVrTOaInMkKmmcVSadhQ7FixUcCpoM63UFq3TtKBWpUBaNeWX7IYsGAwoUkIQEnZ1bICskd19+se9C3dv7j3n3HPec877nvP9zNzZ7O655zz35u5z3vO873lfmhlERKTaFpQdgIiI5E/JXkSkBpTsRURqQMleRKQGlOxFRGpAyV5EpAaU7EVEasBZsifZR3IXyW+42qeIiLjhsmX/AQB7HO5PREQcWehiJySXAfg9AJ8E8Kdx2y9ZssSWL1/u4tAiIrWxc+fOSTMbSvNcJ8kewPUAPgxgIMnGy5cvx9jYmKNDi4jUA8m9aZ+buYxD8gIAz5rZzpjtNpAcIzk2MTGR9bAiItIDFzX7MwFcSPIJAF8BcBbJm9s3MrOtZjZiZiNDQ6muQkREJKXMyd7MrjCzZWa2HMA7APy3mb0rc2QiIuKMxtmLiNSAqw5aAICZfQfAd1zuU0REslPLXkSkBpTs0xocBMiXH4ODZUckItKVkn1azz0X/b2IiEeU7EVEakDJXkSkBpTs0xoYiP5eRMQjTode1srBg2VHICKSmFr2IiI1oGQvIlIDSvZSL7o/QmpKyV7qRfdHSE0p2YuI1ICSvYhIDSjZS73o/gipKY2zl3rR/RFSU2rZi4jUgIsFxxeTvJ/kgyQfJnmli8BERMQdF2WcwwDOMrPnSR4F4C6St5nZfQ72LSIiDmRO9mZmAJ5vfntU82FZ9ysiIu44qdmT7CM5DuBZAHea2fdc7DdYvdylqTs6RaQATpK9mc2Y2UoAywCsIbmifRuSG0iOkRybmJhwcVg/DQ72dpem7ugUkQI4HY1jZj8H8B0A53f43VYzGzGzkaGhIZeH7a6MVrOStaQweWgSq7euxuShybJDkYpyMRpniORxzX8fDeAcAI9k3a8TajVLIK666yqM7x/H1XddXXYoUlEuWvavAvBtkg8BGEWjZv8NB/utjqi7NHVHZ+1NHprE50Y/h1mbxY2jN6p1L7nInOzN7CEzW2VmbzCzFWb2CReBBatT8o66a/PgQcDs5Yfu8Kydq+66CrM2CwCYtVm17iUX1b6DNkmr2XVdX8lberRj7w6QxOKFi0ES39373bJDkgqq9tw4SRJtp7o+2fh3XKs8qfYROr7vV3oyeWgS629ej9vfdTuWvGJJz88fvXQ0h6hE5qt2yz4rVx26rjqK269CiuiA1n0AsdS5KiFQsg9JGaOJNKIpkjpXJRRK9hr9Ihmoc1VCUe2afRKtNe5ONXAXBgby2W+n40gh5ur0szb7UucqAHWuireU7Fvl1bnpar+dThpZ9x3XyVvUiSowc3X6TW/ahF3v31V2OCKxVMbxVaeO0TyGdcbV5Ks8lDRl57Pq9BIiJXtfqWM0fwne405z1qhOLyFSss8ixGGJ7TFLpE7DKnUTlISIjbVHijUyMmJjY2OFH9e5TsnS1fuZ177jEnydbsyKeY8nD01i+LphTE1P4eiFR+PJTU+mumlKxBWSO81sJM1z1bL3lYMJ0hJPm1vVmnycmPdY5RqpEiX7IqQp9yTpGG3fb9v+dWdnjJj3WOUaqRKVcbJIOjdNt9JJ1pJJRElm8oWJziUIzacjEiyVccqSdVji3KRrvXbwzrXoI3QtQVR5KGUBtKKUhErJ3idJh1cm2C5zCcKHkUY+xNBGpTEJVeYyDskTAPwTgOMBzALYamafiXpOZco4SXVahLybqP+PpPtxUZrJc6RRSDG00OgcKVvZZZxpAJvN7FQAbwJwOcnTHOy3OuZKJ1mnGohK9CrN5E6jcyRkLpYlfMbMHmj++zkAewC8Jut+K6m9Xh43vDLpDVCar6YQGp0jIXM6ERrJ5QBWAfheh99tALABAIaHh10eNlxzLfC58sxch+1cGSZJySavsoYPE6D5EEML31aUyrpCltSLsw5akscC+HcAHzSzI+oIZrbVzEbMbGRoaMjVYash7Tw4eSY/H0bt+BCDx9RZLL1wkuxJHoVGov+Smd3iYp9RvBr+VvSIkTySX7fX4OFomDx59bmKoZk3pVeZkz1JAvgHAHvM7NPZQ4pXWIsmSbLLc3ZKB1MmJNLtNdRs5s2QWsrqLJZeuWjZnwng3QDOIjnefPyug/12VGiLpqhk1y2pq4xRmNBaykk6i0O6UpH8uRiNc5eZ0czeYGYrm49vuQiuk0q2aFwn9ZqVX+LEJr3BQVz1liHMTk0BAGZ/MeX952r00lFMbZl66dGp8zikKxXJX3B30Ho3/K2oUktSnW68irsi6fYafHttKcUmveeew47XAgSw+EWAFv5asqFdqUj+NBFalLSThpU52Vi38fgl3nmaSE7vWaK7Xj27U9eFzXdsxg3334DDM4fR39ePjWs24trzri07LMmo7DtoqyttecXHjk3fSzs5vWeVLPsl4N0VsJROyb4O2m9OAvw4ARUgUdLrsVzlQ8dnXAy3XXIbThs6Dfs27eta05d6UbIH/G/19iLQunpeknRk9noF50PHZ1wMSWL04aQlxVGyB3pr9SY5MZTZsdkpUSVpxUe9riJOhoF0BvvQ8RkXQ9IYfThpSXGU7HuV5MSQZknBvK8m4pJp1OsqogQUyD0FPvQBxMWQJEYfTlpSLCX7shRdQw8kmfrOh47PuBiSxOjDSUuK5XTWy2AlmV2xlwVIYkwemsT6DcDtNwNLDjnZ5XydXo/Wnk2lfWZJHzo642JIEmPrCQEI/74CiadkDyRLet0SfYra8lV3XYXx44GrzwSuvbPnp8dLsuh5++uJOuF5NtVwIo5Obq117SqNU/fhpCXFUhknixTlkJdqpQuAG38LmHxF8xdlJ9CoMk+IJSAHZbIkHaFJR7No5IuUTcm+YPNqpUf34+pbN4eTQGsmSUdo0tEscdsmORnohCFZhJnsyxgX72hoYGkdfIEMbfTJjr07AABEowTW+n/Vy2iWJNsmOXFoqKRkEWayL3oki8POzUQ3+eQhxFJMFg5ObqOXjuLyNZeDJDau2Tjv/6qX0Sxx2yY5GWiopGQVZrIvWk2nGgiag5NbVILt5Qqtfdvtj2+fV45JcuLQUEnJSslepIuoBNvLFVr7tme/7ux55ZgkJw4X5T/V/OvNydBLkl8AcAGAZ81shYt9Ruo0sVfr0EKNIRcH8hiL3n618NHf/miiUp6Lcl9Vh5FKMq5a9v8I4HxH+4rXfoneznWZxefOzaSd1b5M9pZzHC5br3n0r5RVjlHNX5wkezPbAeBnLvblJZ87N5P2J/jS75BzHL6PWClrNJZq/lJYzZ7kBpJjJMcmJiaKOqzUSAit19arhX2b9mHWZguJ04c5faRchSV7M9tqZiNmNjI0NOR25z6XWaQwobVei7wK0WImUo3ROD6XWfKW9ESX5oSYR309xxNzSK3Xoq9CfC9vSf40EVrokp7Y0pwA86iv53giDqm12ukqJK8RMp1GAB2x6LpUnpOWPcl/AXAvgFNIPkXyfS72K1JVRV6FhFbeknw4admb2R+62I9IXeR1FdI+/z6gueulQWUcRzr9kQUvxHnsa67TjVMhlbckP9XooPVAJTvAfOn49uWGMM+FMPRUyqNk74D+yHLmyw1hnlNtXqIo2TugPzLxQR6dvpo8rTpUs3dAHWDigzxq85o8rTponSYSy9nIyIiNjY0VflwJlMPFYyS5yUOTGL5uGFPTUzh64dF4ctOT1Rl8ECiSO81sJM1zVcbJizoV3fGlo7hgZZdQVJ6sFiV71+aSvDoVJaOyR3iFNP1ECMo+eSvZu6akXjt5/BH7MMKrtPWSK6rsk7eSvZTHs1JX2qSdxx+xSijV4sPJW8m+KBW4+9R5C9azUleapJ3XH7FKKNXiw8lbyd61TlP4VqRTsezL0Dy1Ju1P3/tpPDr5aKLn5fVHnKWEUnZtWI7kw8lbyd61io4c8eEyNE+tSdtguOSWSxI9L+sfcR6JueyTsk42R/Kh/0PJXhLJpQXr0QpjO/bumPf9A888kChZZf0jdp2Yezkp55WUyz7ZSGfhJHvPOvPqJpfLUI+ugkYvHcXlay5Hf18/AGBR36Lck1XWq6VOyTrqpNy+fR5JuepXgCELJ9l71plXNz5chuat6Lpq1qulTsk66jW0bp9XUvahI1I6czI3DsnzAXwGQB+Avzcz/Q9LEFrXISj6BJZlTqVuSw12ew3t2x968VAuyyJqnih/ZU72JPsA3ADgXABPARgl+XUz+0HWfYvkrcyJvqISc9xCOL2uYdu+/S17bsklKSc9YVZysR/PuSjjrAHwmJk9bma/BPAVABc52O98HnXmSTX4Wl9OUkvvteTUvv0Jrzyh1LKcOnGLl3nWS5JvB3C+mf1x8/t3A3ijmW1s224DgA0AMDw8vHrv3r2ZjiuS1eY7NuOG+2/A4ZnD6O/rx8Y1GxO17vNsldZhpsk6vMa8lD3rJTv87IgziJltNbMRMxsZGhpycFiRbNJ2yLZ3dLocvliHDs46vEYfueigfQrACS3fLwPwtIP9iuSq17tS19+8Hl9+65eP6Oh0WfOvQwdnHV6jj1yUcRYC+CGAswH8BMAogHea2cPdnqPFSyQ0m+/YjOvvux6rjl+F3c/uxuGZw1i0YBFmMYvp2WmVI6QQpZZxzGwawEYAtwPYA+BfoxK9iC+SlmBaO3IfeOYBAMDihYsxYzOYmZ0BoHKE+M/JTVVm9i0zO9nMft3MPulinyJ5SzoipLXGvKhvETau2YipLVNY9apV6F/Yr5kpJQhacFxqqdtNSZ1G2nSrMVfxLmKprnCmSxBJKEl5ptuIkE6t/TpMFSHVp2QvlZP2pqQsN1lpWl/xnZK9VErShN2ptZ5l/LfuCBXfKdlLpXRL2Ela3mlvsvJ12gWRVuqglUrp1pmaZMKztLX4XiclEymDkr1USqeE3W3kjSu6I1RCoGQvQctjOuBeaXSOhEA1ewlaHtMBi1RR5rlx0tDcOOJC1afK1QIf0q7sKY5FSlH1qXI1nFNcUrKXYFW5PJN2OKdu7pJu1EErwapyx2jaTuUy19QVv6llL8GoU6s1zVWLbu6SKGrZSzDq1GpNc9Wim7skilr2EgS1WuNVuQ9DssvUsid5MYCPAzgVwBoz03hKyYVarfGq3Ich2WVt2e8G8FYAOxzEItKVWq0vq1PfhbiTqWVvZnsAgKSbaES6UKv1ZXXquxB3CqvZk9xAcozk2MTERFGHFakU9V1IWrHJnuR2krs7PC7q5UBmttXMRsxsZGhoKH3EIhWUtDRT9buGJT+xyd7MzjGzFR0eXysiQJEQHTgArFsHDA42vh44cOQ2rQk+6dQI6ruQtDTOXiQHF18M3HsvMD3d+HrxxcCOtmEMcwn+Y9/+GG4avyl2vv3JQ5OYtVns27RPE6NJzzLV7Em+heRTAM4A8E2St7sJSyRs4+ONRA80vo6Pz/99a+196wNbMTM7AyC6NKOJ0SSLTMnezG41s2Vm1m9mS81svavAREK2ciWwsHndvHBh4/tWrbX3mdkZzNhMZGlGHbOSle6gleD4Os68Na5t24AzzgAGBhpft22bv21r7b1/YT9WvWoVprZMYWrLVMdhpuqYlaxUs5dEfFpIw9dx5u1xtdfoW/V634DWuZWs1LKXRHypF/tazsg7rtFLR19q+c+1/n29whE/KdlLLJ8SrK/lDBdx9Zq8fTkBSxiU7CWWTwm2rHHmcYnYRVy9JG+fTsASBtXsJZZP9eKy5siJ6ydIE1drPwiAecm721j71ng0C6j0QsleYtV9ErL2VnRcIk6q9QRisJ6St08nYAmDkr14ybfRP65b0e0nkFOWnNJT8q77CVh6p2QvXvJpeGUerej2E8jZJ56NXe/flXm/It0o2Yt38iqbpJVHK1plGCmakr14pw6djyrDSNGU7MU7avWKuKdkL95Rq9d/PnWgSzK6qUpEeqa7d8OjZC8iPdHdu2FSsheRnvg0fYYkl3WlqmtIPkLyIZK3kjzOVWAi4ietgxumrB20dwK4wsymSX4KwBUAPpI9LBHp1YEDjbVux8cbK2Nt2wYsXer+OOpAD1PWZQnvMLPmSpu4D8Cy7CGJhOvAAWDdOmBwsPH1wIHijj23yPlzz728yLnIHJc1+/cCuK3bL0luIDlGcmxiYsLhYUX8UWbCjVvkXOotNtmT3E5yd4fHRS3bbAEwDeBL3fZjZlvNbMTMRoaGhtxEn1GZrTCppqwJN8tnMm6Rc6m32GRvZueY2YoOj68BAMn3ALgAwCVmZnkH7JIue8W1rAk3y2cybpFzqbeso3HOR6ND9kIzO+QmpOLosvdIutrJJmvCzfKZXLoU2LEDOHiw8TWPzlkJV9aa/WcBDAC4k+Q4yc87iKkwuuw9UlzLsi4ng7SvM2vC1WdS8pJ1NM5vmNkJZray+bjMVWBF0GXvkeJalnUpfWV9nWlPFvpMSl5YRpl9ZGTExsbGCj+uxFu3rpHcpqcbLcszzmi0UOcMDjYS4JyBgUYrtmqyvs6491EkDZI7zWwkzXM1XYLME9eyrEuZIevrVH+Q+EbJXuaJqznXpcyQ9XXW5aQo4dB89tKTuZNB1WV9ndu2HTl1gUiZlOxFclCXk6KEQ2UcEY/VZair5E/JXkrnMqFNHprE6q2rc11Qo8gEXJehrpI/JXspTLck6TKh5bVcXmvsJ59cXALWqB5xRcleYrlqyXZL6q4SWp7L5bXGfvBgcQlYo3rEFSV7ieWq5d0tqbtKaHkul9cae6u8E3BdhrpK/pTsJZarlne3pO4qoeW5XF5r7H19jaucIhKwJjcTVzT0UmKtXDn/1v+0LdluY89dDVPMc7m8TrEr8UpI1LKvkbIn52pvpQLlDivs5f1QC1tCp2RfI2lr73klurKHFeZ9fF/HyPsal+RLyb5G0tbe80oOZS7h5+L4cXo9mRSVhMs+yUo5lOxrJO2ol6TJoddkVeYSfi6OH6f9ZHL33dHvSVFJWGP36ynrsoR/SfKh5ipVd5B8tavAxL20tfekyaHXZFXmEn4ujh+n/eQxOxv9nhSVhDV2v6bMLPUDwGDLv/8EwOeTPG/16tUm4Vi71mzhQjOg8XXt2s7bDQw0tpl7LFjQ+NnatWb795cXV1n27zfr65v/ngwMdN++qNezf39j33n+30g+AIxZynyddVnC1rV7jgFQ/LJXkrukLeDWFiPQaMnmWZJoj+vGGxvlo4EB4JWvBI49ttwOyKVLgTe/OXkruqgbqDSyqKbSniXmHgA+CWAfgN0AhiK22wBgDMDY8PBwrmc/ma+ollzrcZK2aF3G1toynnskaSHn+f502rda1pIWMrTskyTz7c1E3v64qG27KwBcmeSgKuNk10vCKKPckfSYLmNrLyMlKZ24jiEJ38tP4q8syd7ZguMkXwvgm2a2Im5bLTieXS8LWpexSPiBA8nuOG2PbcEC4Omn05UWWt+TOUkW+y76/anLou3iXmkLjpM8qeXbCwE8kmV/klwvIzfKGH2RtC7c64gVoPsQz23bgNWrGycMACAb38fVvot+fzQaRsqQdZz91SR3k3wIwHkAPuAgJkmgl4Th88yJ27Y1JhZrdc890WP1uw3xXLoUWLTo5WTf19f4Pu4qoej3x+f/D6kuZ2WcXqiMk13SMkkIx+9UfgG6l2CiyiAqkUiVlVbGkfKUPXzO5d2erS3dBS2fyG7lqairGpVIRDpTspdUXN7t2XriOvPM+GQdVQZRiUSkM81nL6m4muO+Xbc571tFzX/vam58kapRspdUkiTlNJSsRfKhZC+pKCmLhEU1e5HAaPERSUPJXiQwWnxE0lCyFwmMFh+RNJTsRQKjewkkDSV7kcDoXgJJQ8leglBGp6SvHaFl3z0tYVKyFy/EJdYyOiWLOqavJxWpFiV78UJcYi2jU7KoY2p0jRRByV5y0WtrNS6xltEpWdQxNbpGiqBkL7notbUal1jL6JQs6pgaXSNFULKXXPTaWo1LrHOdkj/6UeP7k07Kv75dVEeoRtdIEZzMjUPyQwCuATBkZpMu9ilh63VWzKRz7cxdMUxPv3zFEPocPZpnSIqQuWVP8gQA5wJ4Mns4UhV5tVZV3xZJx0UZ5zoAHwZQ/PqGFVK14Xd5lUBU3xZJJ1OyJ3khgJ+Y2YMJtt1Acozk2MTERJbDVpKG3yWj+rZIOrELjpPcDuD4Dr/aAuDPAJxnZv9H8gkAI0lq9lpw/EhaKFtE4mRZcDy2g9bMzuly0NcDOBHAgyQBYBmAB0iuMbP9aYKps7yW+RMRATKUcczs+2b2a2a23MyWA3gKwG8q0aej8kRD6H0Xoccv1RVbxkm8I5VxxIF16+Zf4ZxxRljDEkOPX/yWaxknqWbrXiST0IdWhh6/VJfuoBWvhD60MvT4pbqU7MUrofddhB6/VJezMo6IC6FPHRB6/FJdatmLiNSAkr2ISA0o2YuI1ICSvYhIDSjZi4jUgJK9iEgNKNmLiNSAkr2ISA04mwitp4OSEwD2OtjVEgChrXmrmIuhmIuhmIsxF/NrzWwozQ5KSfaukBxLOwNcWRRzMRRzMRRzMVzErDKOiEgNKNmLiNRA6Ml+a9kBpKCYi6GYi6GYi5E55qBr9iIikkzoLXsREUkguGRP8hqSj5B8iOStJI9r/nw5ySmS483H58uOtVW3uFt+P0zyeZIfKivGdhHv9ZqW9/lBkm8pO9Y5ETGfS3Inye83v55VdqxzImL+VZLfbn4uPlt2nK2iPs8kryD5GMlHSa4vM85WJC8m+TDJWZIjLT9fRPKm5mfjQZK/U2KY80TEfBTJLzZj3kPyirh9BZfsAdwJYIWZvQHADwG0vsgfm9nK5uOycsLrKipuALgOwG2FRxWtW8y70VhcfiWA8wH8HUlfFsLpFvMkgN83s9cDeA+Afy4pvk66xfwLAH8BwJsGQIuOMZM8DcA7AJyOxmfjRpJ9pUU5324AbwXQvrzMpQDQ/GycC+BvSPqSG7vFfDGA/mbMqwG8n+TyqB358oISM7M7zKy5pDPuA7CszHiSioqb5B8AeBzAw2XE1k23mM3sUMvPFwPwpuMnIuZdZvZ08+cPA1hMsr+MGNtFxPyCmd2FRtL3SsTn+SIAXzGzw2b2vwAeA7CmjBjbmdkeM3u0w69OA/BfzW2eBfBzAF6Mw4+I2QAc02xkHQ3glwAORu0ruGTf5r2Y3xo+keQukt8lubasoBJ4KW6SxwD4CIArS40o3rz3muQbST4M4PsALmv5w/dJ++djztsA7DKzwwXHk0S3mH3WGvNrAOxr+d1TzZ/57EEAF5FcSPJENFrKJ5QcU5x/A/ACgGcAPAngWjP7WdQTfLn0nofkdgDHd/jVFjP7WnObLQCmAXyp+btnAAyb2U9JrgbwVZKnm1nk2c6llHFfCeA6M3ueZDGBtkgZM8zsewBOJ3kqgC+SvM3MCmmBpo25+fPTAXwKwHl5x9l23NQxlyVlzJ0+xIVd+SWJuYMvADgVwBga07jcg8ZrKkTKmNcAmAHwagC/AuB/SG43s8e7HcfLZG9m50T9nuR7AFwA4Gxrjh1tttION/+9k+SPAZyMxn9gIdLEDeCNAN5O8q8BHAdgluQvzKyQDrmUMbc+fw/JFwCsQEHvddqYSS4DcCuAPzKzH+cb5XxZ3+cypIz5KcxvFS8D8HT7c/MSF3OX50wD2DT3Pcl7APzIZVwxx+85ZgDvBPCfZvYigGdJ3o1G6alrsg+ujEPyfDTKHhea2aGWnw/NdQSRfB2AkxDxwovWLW4zW2tmy81sOYDrAfxVUYk+TsR7feJchyzJ1wI4BcATpQTZJiLm4wB8E8AVZnZ3WfF10i1mn0XE/HUA7yDZ3yyJnATg/jJiTIrkK5rlVJA8F8C0mf2g5LDiPAngLDYcA+BNAB6JekJwN1WRfAxAP4CfNn90n5ldRvJtAD6BxuXXDICPmdl/lBTmEbrF3bbNxwE8b2bXFhxeRxHv9bsBfBTAiwBmAXzCzL5aUpjzRMT852iMGGltsZ3X7JArVdRng+QTAAYBLEKj4/A8HxJRTMxb0KjjTwP4oJl50QfBxhDhvwUwhMZ7OW5m65ujWG5H47P8EwDvMzMXs/JmFhHzsQBuQqNzmQBuMrNrIvcVWrIXEZHeBVfGERGR3inZi4jUgJK9iEgNKNmLiNSAkr2ISA0o2YuI1ICSvYhIDSjZi4jUwP8Dvmq+ff3JzaIAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "U2,S2,Vt2 = la.svd(samples2);\n",
    "Us2 = U2[:,0:2];\n",
    "feat2 = np.matmul(np.transpose(Us2),samples2);\n",
    "feat20 = feat2[np.ix_((0,1),np.where(labels==0)[0])];\n",
    "feat21 = feat2[np.ix_((0,1),np.where(labels==1)[0])];\n",
    "feat22 = feat2[np.ix_((0,1),np.where(labels==2)[0])];\n",
    "print(\"There are the two-dimensional features of the centered and normalized dataset.\");\n",
    "plt.scatter(feat20[0,:].tolist(),feat20[1,:].tolist(),15,'r','s');\n",
    "plt.scatter(feat21[0,:].tolist(),feat21[1,:].tolist(),15,'g','^');\n",
    "plt.scatter(feat22[0,:].tolist(),feat22[1,:].tolist(),15,'b','o');"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The second plot will lead to better classification. The data in the first plot does not separate easily, because the blue and green data are interspersed together. However, the data in the second plot is mostly separable because one can easily identify distinct red, blue, and green regions."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Naive Bayes Classification"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "def Gaussian(x, mean, std) :\n",
    "    return 1/np.sqrt(2*np.pi*std**2) * np.exp(-(x-mean)**2/(2*std**2));\n",
    "def NBClassifier(training, trlabels, test) :\n",
    "    # Split training data and calculate means and standard deviations\n",
    "    train0 = training[np.ix_((0,1),np.where(trlabels==0)[1])];\n",
    "    train1 = training[np.ix_((0,1),np.where(trlabels==1)[1])];\n",
    "    train2 = training[np.ix_((0,1),np.where(trlabels==2)[1])];\n",
    "    means = np.concatenate((np.concatenate((train0.mean(1),train1.mean(1)),1),train2.mean(1)),1);\n",
    "    stds = np.concatenate((np.concatenate((train0.std(1),train1.std(1)),1),train2.std(1)),1);\n",
    "    # Classify each test datum based on Gaussian output\n",
    "    result = [0]*np.shape(test)[1];\n",
    "    for i in range(len(result)) :\n",
    "        x1 = test[0,i];\n",
    "        x2 = test[1,i];\n",
    "        prob = [0]*3;\n",
    "        prob[0] = Gaussian(x1,means[0,0],stds[0,0]) * Gaussian(x2,means[1,0],stds[1,0]);\n",
    "        prob[1] = Gaussian(x1,means[0,1],stds[0,1]) * Gaussian(x2,means[1,1],stds[1,1]);\n",
    "        prob[2] = Gaussian(x1,means[0,2],stds[0,2]) * Gaussian(x2,means[1,2],stds[1,2]);\n",
    "        result[i] = np.argmax(prob);\n",
    "    return result;"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Case I: My Classifier on Centered Data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "The accuracies for the five folds are:  [0.7222222222222222, 0.6944444444444444, 0.7777777777777778, 0.8333333333333334, 0.7777777777777778]\n",
      "Mean:  0.7611111111111111\n",
      "Standard deviation:  0.048432210483785286\n"
     ]
    }
   ],
   "source": [
    "labels = np.reshape(labels,(1,178));\n",
    "feat1a = np.concatenate((feat1,labels));\n",
    "np.random.shuffle(feat1a.T);\n",
    "\n",
    "kcross1 = NBClassifier(feat1a[0:2,36:178],feat1a[2,36:178],feat1a[0:2,0:36]);\n",
    "train2 = np.concatenate((feat1a[0:2,0:36],feat1a[0:2,72:178]),1);\n",
    "labels2 = np.concatenate((feat1a[2,0:36],feat1a[2,72:178]),1);\n",
    "kcross2 = NBClassifier(train2,labels2,feat1a[0:2,36:72]);\n",
    "train3 = np.concatenate((feat1a[0:2,0:72],feat1a[0:2,108:178]),1);\n",
    "labels3 = np.concatenate((feat1a[2,0:72],feat1a[2,108:178]),1);\n",
    "kcross3 = NBClassifier(train3,labels3,feat1a[0:2,72:108]);\n",
    "train4 = np.concatenate((feat1a[0:2,0:108],feat1a[0:2,143:178]),1);\n",
    "labels4 = np.concatenate((feat1a[2,0:108],feat1a[2,143:178]),1);\n",
    "kcross4 = NBClassifier(train4,labels4,feat1a[0:2,108:143]);\n",
    "kcross5 = NBClassifier(feat1a[0:2,0:143],feat1a[2,0:143],feat1a[0:2,143:178]);\n",
    "\n",
    "accs = [];\n",
    "accs.append(len(np.where(kcross1-feat1a[2,0:36] == 0)[1])/36);\n",
    "accs.append(len(np.where(kcross2-feat1a[2,36:72] == 0)[1])/36);\n",
    "accs.append(len(np.where(kcross3-feat1a[2,72:108] == 0)[1])/36);\n",
    "accs.append(len(np.where(kcross4-feat1a[2,108:143] == 0)[1])/36);\n",
    "accs.append(len(np.where(kcross5-feat1a[2,143:178] == 0)[1])/36);\n",
    "\n",
    "print(\"The accuracies for the five folds are: \",accs);\n",
    "print(\"Mean: \",np.mean(accs));\n",
    "print(\"Standard deviation: \",np.std(accs));"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Case II: My Classifier on Standardized Data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "The accuracies for the five folds are:  [0.9722222222222222, 0.8888888888888888, 0.9444444444444444, 0.8611111111111112, 0.8888888888888888]\n",
      "Mean:  0.9111111111111111\n",
      "Standard deviation:  0.04082482904638629\n"
     ]
    }
   ],
   "source": [
    "feat2a = np.concatenate((feat2,labels));\n",
    "np.random.shuffle(feat2a.T);\n",
    "\n",
    "kcross1 = NBClassifier(feat2a[0:2,36:178],feat2a[2,36:178],feat2a[0:2,0:36]);\n",
    "train2 = np.concatenate((feat2a[0:2,0:36],feat2a[0:2,72:178]),1);\n",
    "labels2 = np.concatenate((feat2a[2,0:36],feat2a[2,72:178]),1);\n",
    "kcross2 = NBClassifier(train2,labels2,feat2a[0:2,36:72]);\n",
    "train3 = np.concatenate((feat2a[0:2,0:72],feat2a[0:2,108:178]),1);\n",
    "labels3 = np.concatenate((feat2a[2,0:72],feat2a[2,108:178]),1);\n",
    "kcross3 = NBClassifier(train3,labels3,feat2a[0:2,72:108]);\n",
    "train4 = np.concatenate((feat2a[0:2,0:108],feat2a[0:2,143:178]),1);\n",
    "labels4 = np.concatenate((feat2a[2,0:108],feat2a[2,143:178]),1);\n",
    "kcross4 = NBClassifier(train4,labels4,feat2a[0:2,108:143]);\n",
    "kcross5 = NBClassifier(feat2a[0:2,0:143],feat2a[2,0:143],feat2a[0:2,143:178]);\n",
    "\n",
    "accs = [];\n",
    "accs.append(len(np.where(kcross1-feat2a[2,0:36] == 0)[1])/36);\n",
    "accs.append(len(np.where(kcross2-feat2a[2,36:72] == 0)[1])/36);\n",
    "accs.append(len(np.where(kcross3-feat2a[2,72:108] == 0)[1])/36);\n",
    "accs.append(len(np.where(kcross4-feat2a[2,108:143] == 0)[1])/36);\n",
    "accs.append(len(np.where(kcross5-feat2a[2,143:178] == 0)[1])/36);\n",
    "\n",
    "print(\"The accuracies for the five folds are: \",accs);\n",
    "print(\"Mean: \",np.mean(accs));\n",
    "print(\"Standard deviation: \",np.std(accs));"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Case III: sklearn Classifier on Standardized Data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "The accuracies for the five folds are:  [1.0, 0.9166666666666666, 0.9444444444444444, 0.8611111111111112, 0.8888888888888888]\n",
      "Mean:  0.9222222222222222\n",
      "Standard deviation:  0.04779069592801459\n"
     ]
    }
   ],
   "source": [
    "from sklearn.naive_bayes import GaussianNB;\n",
    "\n",
    "gauss = GaussianNB();\n",
    "gauss.fit(feat2a[0:2,36:178].T,np.ravel(feat2a[2,36:178].T));\n",
    "kcross1 = gauss.predict(feat2a[0:2,0:36].T);\n",
    "gauss.fit(train2.T,np.ravel(labels2.T));\n",
    "kcross2 = gauss.predict(feat2a[0:2,36:72].T);\n",
    "gauss.fit(train3.T,np.ravel(labels3.T));\n",
    "kcross3 = gauss.predict(feat2a[0:2,72:108].T);\n",
    "gauss.fit(train4.T,np.ravel(labels4.T));\n",
    "kcross4 = gauss.predict(feat2a[0:2,108:143].T);\n",
    "gauss.fit(feat2a[0:2,0:143].T,np.ravel(feat2a[2,0:143].T));\n",
    "kcross5 = gauss.predict(feat2a[0:2,143:178].T);\n",
    "\n",
    "accs = [];\n",
    "accs.append(len(np.where(kcross1-feat2a[2,0:36] == 0)[1])/36);\n",
    "accs.append(len(np.where(kcross2-feat2a[2,36:72] == 0)[1])/36);\n",
    "accs.append(len(np.where(kcross3-feat2a[2,72:108] == 0)[1])/36);\n",
    "accs.append(len(np.where(kcross4-feat2a[2,108:143] == 0)[1])/36);\n",
    "accs.append(len(np.where(kcross5-feat2a[2,143:178] == 0)[1])/36);\n",
    "\n",
    "print(\"The accuracies for the five folds are: \",accs);\n",
    "print(\"Mean: \",np.mean(accs));\n",
    "print(\"Standard deviation: \",np.std(accs));"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "The mean accuracy is significantly higher for the standardized data set than for the only centered data set.\n",
      "The accuracy of my Bayesian classifier was about the same as the accuracy of the built-in classifier in scikit learn.\n",
      "The standard deviations were about the same for every classification.\n"
     ]
    }
   ],
   "source": [
    "print(\"The mean accuracy is significantly higher for the standardized data set than for the only centered data set.\");\n",
    "print(\"The accuracy of my Bayesian classifier was about the same as the accuracy of the built-in classifier in scikit learn.\");\n",
    "print(\"The standard deviations were about the same for every classification.\");"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "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.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}