{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Линейная регрессия\n",
"Шестаков А.В. Майнор по анализу данных - 01/03/2016"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Сегодня мы рассмотрим следующие темы:\n",
"\n",
"1. Постановка задачи линейной регрессии\n",
"2. Преобразование переменных и интерпретация\n",
"3. Переобучение\\недообучение, мультиколлинеарность и регуляризация"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## 1. Постановка задачи"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Дано описание $n$ объектов по $m$ признакам. Обычно оно выражается в виде матрицы размера $n \\times m$: $X = [x^{(i)}_j]^{i=1\\dots n}_{j=1\\dots m} $.
($x^{(i)}_j$ означает $j$-ый признак $i$-го объекта)
\n",
"Дана **вещественная** зависимая переменная, которая тоже имеет отношение к этим объекам: $y \\in \\mathbb{R}^n$.\n",
"\n",
"Наша задача, выявить **линейную** зависимость между признаками в $X$ и значениями в $y$:\n",
"$$\\hat{y} = X\\beta \\quad \\Leftrightarrow \\quad \\hat{y}^{(i)} = \\beta_0 + \\beta_1x^{(i)}_1 + \\dots$$\n",
"То есть необходимо оценить коэффициенты $\\beta_i$.\n",
"\n",
"В случае линейной регрессии коэффициенты $\\beta_i$ рассчитываются так, чтобы минимизировать сумму квадратов ошибок по всем наблюдениям:\n",
"$$ L(\\beta) = \\frac{1}{2n}(\\hat{y} - y)^{\\top}(\\hat{y} - y) = \\frac{1}{2n}(X\\beta - y)^{\\top}(X\\beta - y) \\rightarrow \\min$$ $$ \\Updownarrow $$ $$ L(\\beta_0,\\beta_1,\\dots) = \\frac{1}{2n}\\sum^{n}_{i=1}(\\hat{y}^{(i)} - y^{(i)})^2 = \\frac{1}{2n}\\sum^{n}_{i=1}(\\beta_0 + \\beta_1x^{(i)}_1 + \\dots - y^{(i)})^2 \\rightarrow \\min $$\n",
"\n",
"Мы уже знаем несколько способов решения этой задачи (Семинар 4):\n",
"* Градиентный спуск \n",
"* Normal Equations (Проекционные матрицы)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Библиотеки для расчетов"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Самыми распространенными библиотеками в `Python` для работы с линейными методами (регрессии и классификации) являются [`scikit-learn`](http://scikit-learn.org/stable/) и [`statmodels`](http://statsmodels.sourceforge.net/). Так как в дальнейшем мы скорее всего будем работать со `scikit-learn`, то и примеры по большей части будут демонстрироваться именно в нем."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Пример: Стоимость автомобиля"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Загрузите [данные](http://bit.ly/1gIQs6C) по характеристикам автомобилей Honda Accord. Названия столбцов говорят сами за себя."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"import matplotlib.pyplot as plt\n",
"import numpy as np\n",
"import pandas as pd\n",
"\n",
"plt.style.use('ggplot')\n",
"\n",
"%matplotlib inline"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"df = pd.read_csv('http://bit.ly/1gIQs6C')"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"df.head()"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Выберем одну переменную mileage в качестве предиктора, а переменную price в качестве зависимой переменной"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"y = df.price.values\n",
"X = df.mileage.values.reshape(-1,1)\n",
"# В последних версиях sklearn начинает ругаться на одномерные данные (когда array.shape = (m,))\n",
"# Поэтому с помошщью .reshape(-1,1) мы искусственно добавляем единичную размерность к Х"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"from sklearn.linear_model import LinearRegression"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"model = LinearRegression(fit_intercept=True)\n",
"model.fit(X, y)"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"print 'Модель:\\nprice = %.2f + (%.2f)*mileage' % (model.intercept_, model.coef_[0])"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"x = np.linspace(0, max(X), 100)\n",
"y_line = model.intercept_ + model.coef_[0]*x\n",
"\n",
"fig, ax = plt.subplots(1, 1, figsize=(10,5))\n",
"ax.scatter(X, y)\n",
"\n",
"ax.plot(x, y_line, c='red')\n",
"\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Остатки и меры качества"
]
},
{
"cell_type": "markdown",
"metadata": {
"collapsed": true
},
"source": [
"Давайте взглянем на ошибки (остатки)"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"fig, ax = plt.subplots(1, 2, figsize=(14, 5))\n",
"\n",
"y_hat = model.predict(X)\n",
"res = y - y_hat\n",
"ax[0].hist(res)\n",
"ax[0].set_xlabel('residuals')\n",
"ax[0].set_ylabel('counts')\n",
"\n",
"ax[1].scatter(X, res)\n",
"ax[1].set_xlabel('mileage')\n",
"ax[1].set_ylabel('residuals')\n"
]
},
{
"cell_type": "markdown",
"metadata": {
"collapsed": true
},
"source": [
"Важно смотреть на остатки.
\n",
"Во-первых, они должны быть нормально распределены (Теорема Гаусса-Маркова).
\n",
"Во-вторых, не должно быть ярких зависимостей между значениями признака и остатками.\n",
"\n",
"Посмотрим на меры качества"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"from sklearn.metrics import mean_absolute_error\n",
"from sklearn.metrics import median_absolute_error\n",
"from sklearn.metrics import mean_squared_error\n",
"from sklearn.metrics import r2_score"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Можно посчитать простые варианты агрегирования остатков, например:\n",
"\n",
"* $\\frac{1}{n} \\sum_i |\\hat{y}^{(i)}-y^{(i)}|$ - средняя абсолютная ошибка\n",
"* $\\frac{1}{n} \\sum_i (\\hat{y}^{(i)}-y^{(i)})^2$ - средняя квадратичная ошибка"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"print 'Средняя абсолютная ошибка %.2f' % mean_absolute_error(y, y_hat)\n",
"print 'Средняя квадратичная ошибка %.2f' % mean_squared_error(y, y_hat)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Можно рассмотреть более сложную меру: коэффициент детерминации $R^2$:\n",
"\n",
"* $TSS = \\sum_i (y^{(i)}-\\bar{y})^2$ - общая сумма квадратов (total sum of squares)\n",
"* $RSS = \\sum_i (\\hat{y}^{(i)}-y^{(i)})^2$ - сумма квадратов остатков (residual sum of squares)\n",
"* $ESS = \\sum_i (\\hat{y}^{(i)}-\\bar{y})^2$ - объясненная сумма квадратов (explained sum of squares)\n",
"\n",
"Для простоты будем считать, что\n",
"$$TSS = ESS + RSS$$\n",
"\n",
"Тогда Коэффициент детерминации $R^2=1-\\frac{RSS}{TSS}$\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"print 'R^2 %.2f:' % r2_score(y, y_hat)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Тоже самое через `statmodels`"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"import statsmodels.api as sm\n",
"model = sm.OLS(y, sm.add_constant(X))\n",
"results = model.fit()\n",
"print(results.summary())"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## 2. Преобразование переменных"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Нормализация"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Как мы с вами убедились из ДЗ2, может решить всё)\n",
"Переход к близким или единым шкалам улучшает сходимость градиентного спуска, уменьшает риск переполнения разрядности чисел, однако приходится жертвовать прямой интерпретируемостью..\n",
"\n",
"Нормализацию обычно проделывают для вещественных признаков.\n",
"\n",
"Нормализация z-score:\n",
"1. Вычитаем среднее: $x - \\bar{x}$\n",
"2. Делим на стандартное отклонение: $\\frac{x - \\bar{x}}{std(x)}$\n",
"\n",
"Можно проделать вручную, можно с помошью метода ниже"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"from sklearn.preprocessing import StandardScaler\n",
"\n",
"scaler = StandardScaler()\n",
"X = scaler.fit_transform(X)\n",
"\n",
"model = LinearRegression(fit_intercept=True)\n",
"model.fit(X, y)"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"print 'Модель:\\nprice = %.2f + (%.2f)*mileage`' % (model.intercept_, model.coef_[0])"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Номинальная шкала"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"df.head()"
]
},
{
"cell_type": "markdown",
"metadata": {
"collapsed": true
},
"source": [
"Обратим внимание, что в нашем DataFrame есть и другие переменные - год, особенности внешности, количество цилиндров, тип коробки передач. Как учесть их в модели? Иными словами, как их перекодировать, для использование в модели регресии?\n",
"\n",
"Часть из них можно представить в виде так называемых dummy (фиктивных) переменных"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"from sklearn.feature_extraction import DictVectorizer\n",
"cols = df.columns[1:]\n",
"\n",
"dv = DictVectorizer()\n",
"dv.fit(df[cols].T.to_dict().values())\n",
"\n",
"X = dv.transform(df[cols].T.to_dict().values())"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Давайте поймем, как устроен X"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": []
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"#### Пример: Стоимость автомобиля - больше переменных!"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Обратите внимание, что в случае, когда мы добавляем все категории в матрицу `X`, рассчитывать вес для свободного члена не нужно.\n",
"\n",
"Причиной тому - так называемая Dummy Variable Trap. Матрица $X^\\top X$ перестаёт быть обратимой из-за линейной зависимости столбцов (своего рода мультиколлинеарность).
\n",
"Хотя, даже если поставить `fit_intercept=True`, sklearn вас даже не остановит..."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"model = LinearRegression(fit_intercept=False)\n",
"model.fit(X, y)"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"# Выведите значение коэффициента для каждой переменной и проинтерпретируйте их значение\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"y_hat = model.predict(X)\n",
"\n",
"print 'Средняя абсолютная ошибка %.2f' % mean_absolute_error(y, y_hat)\n",
"print 'Средняя квадратичная ошибка %.2f' % mean_squared_error(y, y_hat)\n",
"print 'R^2 %.2f:' % r2_score(y, y_hat)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Природа зависимости"
]
},
{
"cell_type": "markdown",
"metadata": {
"collapsed": true
},
"source": [
"Далеко не всегда переменные зависят друг от друга именно в том виде, в котором они даны. Никто не запрещает зависимость вида\n",
"$$\\log(y) = \\beta_0 + \\beta_1\\log(x_1)$$\n",
"или\n",
"$$y = \\beta_0 + \\beta_1\\frac{1}{x_1}$$\n",
"или\n",
"$$y = \\beta_0 + \\beta_1\\log(x_1)$$\n",
"или\n",
"$$y = \\beta_0 + \\beta_1 x_1^2 + \\beta_2 x_2^2 + \\beta_3 x_1x_2 $$\n",
"и т.д.\n",
"\n",
"Не смотря на то, что могут возникать какие-то нелинейные функции - всё это сводится к **линейной** регрессии (например, о втором пункте, произведите замену $z_1 = \\frac{1}{x_1}$)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"#### Пример: Вес тела - мозгов"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Загрузите [данные](https://www.dropbox.com/s/8srfeh34lnj2cb3/weights.csv?dl=0) и информацией о весах мозга и тел различных биологических видов. Вес тела задан в килограммах, вес могза в граммах."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"df = pd.read_csv('weights.csv', sep=';', index_col=0)\n",
"df.head()"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"df.plot(x = 'body_w', y='brain_w', kind='scatter')\n",
"for k, v in df.iterrows():\n",
" plt.annotate(k, v[:2])\n",
"# Должно получится что-то несуразное.."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Теперь давайте возьмем логарифм от обеих переменных и сонова нарисуем их на графике"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"# Your Code Here\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Постройде линейную регрессию над логарифмами значений. Найдите коэффициенты и проинтерпретируйте их"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"# Your Code Here\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## 3. Переобучение\\недообучение, мультиколлинеарность и регуляризация"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Одна из важнейших характеристик моделей, будь то линейная регрессия, наивные Байес и др. - их **обобщающая способность**.\n",
"Наша задача не построить \"идеальную\" модель, на имеющихся у нас наблюдениях, которая идеально их будет предсказывать, но и применять эту модель для новых данных.\n",
"\n",
"Ниже приводятся примеры 3х моделей."
]
},
{
"cell_type": "markdown",
"metadata": {
"collapsed": true
},
"source": [
"\n",
"[Andrew's Ng Machine Learning Class - Stanford]"
]
},
{
"cell_type": "markdown",
"metadata": {
"collapsed": true
},
"source": [
"Второй момент, который важен для линейных моделей - **мультиколлинеарность**. Этот эффект возникает, когда пара предикторов близка к взаимной линейной зависимости (коэффициент корреляции по модулю близок к 1). Из-за этого:\n",
"\n",
"* Матрица $X^{\\top} X$ становится плохо обусловленной или необратимой\n",
"* Зависимость $y = \\beta_0 + \\beta_1x_1 + \\beta_2x_2$ перестаёт быть одназначной\n",
"\n",
"С этим эффектом можно бороться несколькими способами\n",
"\n",
"* Последовательно добавлять переменные в модель\n",
"* Исключать коррелируемые предикторы"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Регуляризация"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"В обоих случаях может помочь **регуляризация** - добавление штрафного слагаемого за сложность модели в функцию потерь. В случае линейной регрессии было:\n",
"$$ L(\\beta_0,\\beta_1,\\dots) = \\frac{1}{2n}\\sum^{n}_{i=1}(\\hat{y}^{(i)} - y^{(i)})^2 $$\n",
"Стало (Ridge Regularization)\n",
"$$ L(\\beta_0,\\beta_1,\\dots) = \\frac{1}{2n}[ \\sum^{n}_{i=1}(\\hat{y}^{(i)} - y^{(i)})^2 + \\lambda\\sum_{j=1}^{m}\\beta_j^2]$$\n",
"или (Lasso Regularization)\n",
"$$ L(\\beta_0,\\beta_1,\\dots) = \\frac{1}{2n}[ \\sum^{n}_{i=1}(\\hat{y}^{(i)} - y^{(i)})^2 + \\lambda\\sum_{j=1}^{m}|\\beta_j|]$$"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"# В sklearn эти методы называются так\n",
"from sklearn.linear_model import Lasso\n",
"from sklearn.linear_model import Ridge"
]
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 2",
"language": "python",
"name": "python2"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 2
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython2",
"version": "2.7.11"
}
},
"nbformat": 4,
"nbformat_minor": 0
}