04-LDA-lsg.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Quadratische Diskriminantenanalyse\n",
    "\n",
    "In dieser Übung werden Sie selbst eine quadratische Diskriminantenanalyse (QDA) implementieren. Zur Erinnerung: Die QDA berechnet $p(x|y)=\\frac{p(y|x)*p(x)}{p(y)}$. Die Likelihood $p(y|x)$ wird als normalverteilt angenommen."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Aufgabe 1\n",
    "Eine Fischerin benötigt Ihre Hilfe bei der Klassifikation von Fischen. Vor kurzem hat sie folgende Fische gefangen: \n",
    "\n",
    "| Länge (m)    | Art          | \n",
    "| ------------- |-------------  |\n",
    "| 1.3           | Barsch       |\n",
    "| 0.7           | Lachs       |\n",
    "| 0.62           | Lachs      |\n",
    "| 0.9           | Lachs       |\n",
    "| 0.91          | Barsch       |\n",
    "| 0.31          | Hering       |\n",
    "| 0.26           | Hering       |\n",
    "\n",
    "* Berechnen Sie die Priors $p(\\omega)$ für jede Fischart\n",
    "* Berechnen Sie die Parameter $\\mu$ und $\\sigma^2$ für die Lkelihoods $p(x|\\omega)$. \n",
    "* Die Fischerin fängt einen neuen Fisch mit Länge $x = 0.82 m$. Berechen Sie die Posterior-Wahrscheinlichkeit $p(\\omega|x)$ für jede Klasse. Wie wird der Fisch klassifiziert?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "p_barsch = 2/7\n",
    "p_lachs = 3/7\n",
    "p_hering = 2/7\n",
    "\n",
    "# Das sind die Schulvarianten von Mittelwert und Varianz. Nicht die Bias-bereinigten\n",
    "mean_barsch = (1.3 + 0.91) / 2\n",
    "mean_lachs = (0.7 + 0.62 + 0.9) / 3\n",
    "mean_hering = (0.31 + 0.26) / 2\n",
    "\n",
    "var_barsch = 1/2 * (1.3 - mean_barsch) ** 2 + (0.91 - mean_barsch) ** 2\n",
    "var_lachs = 1/3 * (0.7 - mean_lachs) ** 2 + (0.62 - mean_lachs) ** 2 + (0.9 - mean_lachs) ** 2\n",
    "var_hering = 1/2 * (0.31 - mean_hering) ** 2 + (0.26 - mean_hering) ** 2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "def gaussian(x, mean, var):\n",
    "    return 1./(np.sqrt(2.*np.pi*var))*np.exp(-np.power((x - mean), 2.)/(2*var))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.23416954059335204\n",
      "0.7847669667187298\n",
      "1.8808582763433165e-66\n",
      "1.0189365073120817\n",
      "0.2298175979689676\n",
      "0.7701824020310325\n",
      "1.8459033147266004e-66\n",
      "1.0\n"
     ]
    }
   ],
   "source": [
    "x = 0.82\n",
    "\n",
    "t_barsch = gaussian(x, mean_barsch, var_barsch) * p_barsch\n",
    "t_lachs = gaussian(x, mean_lachs, var_lachs) * p_lachs\n",
    "t_hering = gaussian(x, mean_hering, var_hering) * p_hering\n",
    "\n",
    "p_x = sum([t_barsch, t_lachs, t_hering])\n",
    "\n",
    "w_barsch = t_barsch / p_x\n",
    "w_lachs = t_lachs / p_x\n",
    "w_hering = t_hering / p_x\n",
    "\n",
    "print(t_barsch, t_lachs, t_hering, p_x, sep=\"\\n\")\n",
    "\n",
    "print(w_barsch, w_lachs, w_hering, sum([w_barsch, w_lachs, w_hering]), sep=\"\\n\")\n",
    "#Der Fisch wird also als Lachs klassifiziert"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Aufgabe 2\n",
    "Implementieren Sie eine Funktion `priors(classes)`, die für einen Vektor von Klassen-Labels den Prior $p(x)$ für jede Klasse ausgibt.\n",
    "Die Eingabe soll ein Array von Klassen sein (z.b. `np.array([\"stand\",\"sit\",\"sit\",\"stand\"])`). Die Ausgabe soll ein Data Frame mit den Spalten `class` und `prior` sein."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "   class prior\n",
      "0    sit   0.6\n",
      "1  stand   0.4\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "array(['sit', 'stand'], dtype=object)"
      ]
     },
     "execution_count": 2,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "import numpy as np\n",
    "import pandas as pd\n",
    "\n",
    "def priors(classes):\n",
    "    unique, counts = np.unique(classes, return_counts=True)\n",
    "    counts = counts/counts.sum()\n",
    "    df = pd.DataFrame(np.array([unique,counts]).transpose(),columns=[\"class\",\"prior\"])\n",
    "    return(df)\n",
    "    \n",
    "pp = priors(np.array([\"stand\",\"sit\",\"sit\",\"sit\",\"stand\"]))\n",
    "print(pp)\n",
    "np.array(pp[\"class\"])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Aufgabe 3\n",
    "Implementieren Sie eine Funktion `likelihood(data)`, die für ein Data Frame, bestehend aus einer Spalte $y$ und einer Spalte $x$, die Likelihood $p(y|x)$ für jede Klasse $x$ mit einer Normalverteilung approximiert, d.h. es soll für jede Klasse ein Mittelwert und eine Varianz ausgegeben werden.\n",
    "Die Ausgabe soll also die Spalten `class`, `mean` und `variance` besitzen.\n",
    "\n",
    "Plotten Sie die Likelihood für jede Klasse."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "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>mean</th>\n",
       "      <th>variance</th>\n",
       "      <th>class</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <td>0</td>\n",
       "      <td>523.545587</td>\n",
       "      <td>15.843081</td>\n",
       "      <td>b'bike'</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <td>1</td>\n",
       "      <td>517.316002</td>\n",
       "      <td>18.523919</td>\n",
       "      <td>b'downstairs'</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <td>2</td>\n",
       "      <td>467.870175</td>\n",
       "      <td>56.697558</td>\n",
       "      <td>b'lie'</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <td>3</td>\n",
       "      <td>527.009294</td>\n",
       "      <td>96.483175</td>\n",
       "      <td>b'run'</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <td>4</td>\n",
       "      <td>486.721946</td>\n",
       "      <td>71.790221</td>\n",
       "      <td>b'sit'</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <td>5</td>\n",
       "      <td>528.594727</td>\n",
       "      <td>1.057665</td>\n",
       "      <td>b'stand'</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <td>6</td>\n",
       "      <td>545.413484</td>\n",
       "      <td>62.520043</td>\n",
       "      <td>b'upstairs'</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <td>7</td>\n",
       "      <td>525.170962</td>\n",
       "      <td>1.221023</td>\n",
       "      <td>b'walk'</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "         mean   variance          class\n",
       "0  523.545587  15.843081        b'bike'\n",
       "1  517.316002  18.523919  b'downstairs'\n",
       "2  467.870175  56.697558         b'lie'\n",
       "3  527.009294  96.483175         b'run'\n",
       "4  486.721946  71.790221         b'sit'\n",
       "5  528.594727   1.057665       b'stand'\n",
       "6  545.413484  62.520043    b'upstairs'\n",
       "7  525.170962   1.221023        b'walk'"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "def likelihood(data):\n",
    "    uc = np.unique(data[\"class\"])\n",
    "    def getMV(c):\n",
    "        dsub = data.loc[data[\"class\"]==c,\"x\"]\n",
    "        return([dsub.mean(),dsub.var()])\n",
    "    mvs = list(map(getMV, uc))\n",
    "    r = pd.DataFrame(mvs,columns = [\"mean\",\"variance\"])\n",
    "    r[\"class\"] = uc\n",
    "    return r\n",
    "    \n",
    "data = arff.loadarff('features1.arff')\n",
    "df = pd.DataFrame(data[0])\n",
    "\n",
    "dat = df.loc[:, [\"AccX_mean\",\"class\"]]\n",
    "dat.columns = [\"x\",\"class\"]\n",
    "lik = likelihood(dat)\n",
    "lik"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Aufgabe 4\n",
    "Implementieren Sie eine Funktion mylda(newdat,lik,priors), die für eine neue Beobachtung `newdat` die wahrscheinlichste Klasse zurückgibt.\n",
    "\n",
    "Testen Sie Ihre Implementierung auf dem Datensatz `features1.arff`. „Trainieren“ Sie die QDA (d.h. berechnen Sie likelihood und prior), und führen Sie dann für die gleichen Daten eine Klassifikation durch. Wie gut ist die Klassifikation?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.6638418079096046\n"
     ]
    }
   ],
   "source": [
    "from scipy.io import arff\n",
    "import scipy.stats\n",
    "\n",
    "def mylda(newdat,lik,prior):\n",
    "    #Berechne für jeden Datenpunkt: Wahrscheinlichkeit pro Klasse als prior*likelihood, waehle wahrscheinlichste\n",
    "    def getClass(d):\n",
    "        def getProb(c,d):\n",
    "            pr = prior.loc[prior[\"class\"]==c,\"prior\"].values[0]\n",
    "            m = lik.loc[lik[\"class\"] == c,\"mean\"].values[0]\n",
    "            v = lik.loc[lik[\"class\"] == c,\"variance\"].values[0]\n",
    "            li = scipy.stats.norm(m, np.sqrt(v)).pdf(d)\n",
    "            return(li*pr)\n",
    "        probs = np.array([getProb(c,d) for c in prior[\"class\"].values])\n",
    "        return prior.loc[np.argmax(probs),\"class\"]\n",
    "    newclasses = np.array([getClass(dat) for dat in newdat])\n",
    "    return newclasses\n",
    "\n",
    "\n",
    "data = arff.loadarff('features1.arff')\n",
    "df = pd.DataFrame(data[0])\n",
    "\n",
    "dat = df.loc[:, [\"AccX_mean\",\"class\"]]\n",
    "dat.columns = [\"x\",\"class\"]\n",
    "\n",
    "lik = likelihood(dat)\n",
    "prior = priors(dat[\"class\"])\n",
    "\n",
    "nc = mylda(dat[\"x\"],lik,prior)\n",
    "print(sum(nc == dat[\"class\"])/len(dat))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [],
   "source": [
    "scipy.stats.norm?"
   ]
  }
 ],
 "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.7.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}