{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Problem sheet 3\n",
    "The previous exercises gave an introduction to Python, Numpy and Pandas. Beginning with this exercise, we shift our focus to statistical learning itself. To this end, we will employ the module scikit-learn which offers many functions we will cover over the remaining semester.\n",
    "\n",
    "If not already done, please download the file [Advertising.csv](https://www.tu-chemnitz.de/mathematik/numa/lehre/ds-2018/exercises/Advertising.csv) and move it into a subfolder called `datasets`."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Exercise 1:\n",
    "We start this exercise with the Advertising dataset known from the lecture."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We read the dataset using Pandas:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "import pandas as pd\n",
    "import numpy as np\n",
    "\n",
    "adv = pd.read_csv('./datasets/Advertising.csv', index_col=0)\n",
    "\n",
    "# Print first entries of adv\n",
    "print(adv.head(3))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For convenience, we extract the values from this pandas-DataFrame"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "X = adv.values[:,0:3]\n",
    "tv, radio, newspaper = np.hsplit(X,3)\n",
    "Y = adv.values[:,3]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Part (a)\n",
    "Compute for each of the 3 predictor variables **TV**, **radio** and **newspaper** simple (1-dimensional) linear regressions, e.g.\n",
    "\n",
    "\n",
    "$$ y^{TV}_i \\approx \\beta_0^{TV} + \\beta_1^{TV} \\, x_i^{TV}$$\n",
    "\n",
    "Use the following function:\n",
    "\n",
    "    from sklearn.linear_model import LinearRegression\n",
    "    \n",
    "You can use a command similar to\n",
    "\n",
    "    print('y = %5.4f + %5.4f x TV' % (intercept, lincoef))\n",
    "    \n",
    "to print your results in a nice fashion."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Put your code always into these code blocks"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "You should observe, that the regression coefficients for **TV** and **newspaper** are very similar.\n",
    "As you already know from the lecture, it is not satisfying from a mathematical point of view to restrict our investigation to the absolute values of the coefficients.\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Part (b)\n",
    "\n",
    "In the lecture you learned about different measures for assessing the quality of a linear fit.\n",
    "In the last exercise, we already implemented a function to compute the mean squared error (MSE).\n",
    "\n",
    "This time, we want to compare the $R^2$ scores. You can use the method `score()` of a `LinearRegression` to get the $R^2$ values.\n",
    "Remember that this value is the proportion of variability in $Y$ explained using **TV**, **radio** or **newspaper** as predictor in a 1-dimensional linear regression fit."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Put your code for part (b) here."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Part (c)\n",
    "Now we want to compute the predicted value of sales if we restrict our prediction to one input, i.e. **TV**, **radio** or **newspaper**, resp.\n",
    "Predict the values $\\hat{y}^{TV}$ $\\hat{y}^{radio}$ and $\\hat{y}^{newspaper}$ using the method `predict()`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# and that of part (c) here..."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Part (d)\n",
    "\n",
    "Plot the datapoints as well as the corresponding regression line for each of the inputs  **TV**, **radio** or **newspaper**.\n",
    "\n",
    "You can use the functions `subplots` or `fig.add_subplot` to arrange the plots in one figure."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# We plot our findings using subplots\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "fig = plt.figure()\n",
    "fig.add_subplot(1,3,1)\n",
    "# I guess you know what you have to put in here ...\n",
    "\n",
    "fig.add_subplot(1,3,2)\n",
    "# ... and here ...\n",
    "\n",
    "fig.add_subplot(1,3,3)\n",
    "# and here.\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Part (e)\n",
    "Take a closer look at the correlation matrix.\n",
    "You can use the method `corr()` that is implemented for pandas `DataFrames`.\n",
    "Which features are correlated most strongly?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Now, I think you know how to proceed."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Answer**: "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Part (f)\n",
    "Investigate the statistical significance of the medium **newspaper** in a linear regression involving only this feature. Use a **t-test** for this purpose as described on slide 80 in the lecture notes.\n",
    "\n",
    "You should observe the following values:\n",
    "\n",
    "|Coefficient | Estimate | SE | t-statistic | p-value|\n",
    "|:-----------|----------|----|-------------|--------|\n",
    "| $\\beta_0$  | 12.351 | 0.621 | 19.88 | < 0.0001 |\n",
    "| $\\beta_{newspaper}$ | 0.055 | 0.017 | 3.30 | 0.00115\n",
    "\n",
    "You should use `scipy` to get the $t$-distribution using\n",
    "\n",
    "    from scipy.stats import t\n",
    "    \n",
    "The cumulative distribution function at a point `x` for `n` degrees of freedom can than be called by\n",
    "\n",
    "    t.cdf(x, n)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Part (g)\n",
    "Now construct a linear regression on all three predictor variables, i.e.\n",
    "\n",
    "$$y_i ≈ \\beta_0 + \\beta_{TV} x^{TV}_i + \\beta_{radio} x^{radio}_i + \\beta_{newspaper} x^{newspaper}_i$$ \n",
    "\n",
    "What do you observe? Compare your results with your findings from above."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Answer**: "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Part (h)\n",
    "What portion of the variance is explained by this linear regression fit?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Part (i)\n",
    "Now perform a linear regression that incorporates only the predictors **TV** and\n",
    "**radio**.\n",
    "Compute also the $R^2$-value and compare it to the full multiple linear regression.\n",
    "\n",
    "**Extra task**: Present the datapoints and the regression plane in a 3-dimensional plot.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Answer**: "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Homework:\n",
    "\n",
    "We have already seen, that the **t-test** comes in handy when one has to decide whether a coefficient for a single feature is significant or not.\n",
    "As has been outlined in the lecture, one can also use the t-test in a multiple linear regression fit\n",
    "\n",
    "$$ Y = X \\beta + \\varepsilon $$\n",
    "\n",
    "while the intercept is incorporated into $X$, i.e. a column containing only ones is stacked in front of the original matrix $X$.\n",
    "\n",
    "The formula to compute the test statistic in this generalized setting is\n",
    "\n",
    "$$ t_j = \\frac{\\hat{\\beta}_j}{\\hat{\\sigma} \\sqrt{v_j}} $$\n",
    "\n",
    "while $\\beta_j$ is the $j$-th entry of the coefficient vector\n",
    "\n",
    "$$ \\beta = (X^\\top X)^{-1} X^\\top y, $$\n",
    "\n",
    "$\\hat{\\sigma}$ is the unbiased estimate of $\\sigma$, which is determined by\n",
    "\n",
    "$$ \\hat{\\sigma} = \\sqrt{\\frac{1}{n-p-1} \\, \\sum_{i=1}^n (y_i - \\hat{y}_i)^2} $$\n",
    "\n",
    "and $v_j$ is the $j$-th diagonal element of the matrix $(X^\\top X)^{-1}$.\n",
    "\n",
    "Then $t_j$ is distributed according to a $t$-distribution with $n-p-1$ degrees of freedom (dofs). \n",
    "\n",
    "**Task**: Compute the values in the following statistic and try to print it in a similar way. \n",
    "\n",
    "|  Coefficient     | Estimate  | SE    | t-statistic | p-value |\n",
    "|:-----------------|-----------|-------|-------------|---------|\n",
    "| $\\beta_0$        | 2.939     |0.3119 | 9.42        | < 0.0001|\n",
    "| $\\beta_{TV}$     | 0.046     |0.0014 | 32.81       | < 0.0001|\n",
    "| $\\beta_{radio}$  | 0.189     |0.0086 | 21.89       | < 0.0001|\n",
    "| $\\beta_{news}$   | −0.001    |0.0059 | −0.18       | 0.8599  |\n",
    "\n"
   ]
  },
  {
   "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.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}