Подводные камни линейной регрессии и как их избежать

Подводные камни линейной регрессии и как их избежать

Что делать, если предположения линейной регрессии не верны

Новичка в науке о данных всегда можно отличить по скорости, с которой они подбирают нейронную сеть.

Нейронные сети - это круто, и они могут делать потрясающие вещи, которые для многих из нас (включая меня) являются причиной, по которой мы вообще занялись наукой о данных. Я имею в виду, кто занимается наукой о данных, чтобы поиграть со старыми моделями линейной регрессии?

Тем не менее, ирония ситуации заключается в том, что, если вы не работаете в специальной области, такой как компьютерное зрение или обработка естественного языка, в большинстве случаев простые модели, такие как линейная регрессия, на самом деле обеспечивают лучшее решение вашей проблемы, чем сложные модели черного ящика, такие как нейронные сети и опорные векторные машины.

В конце концов, модели линейной регрессии:

• быстро обучать и запрашивать;
• не склонен к переобучению и эффективно использует данные, поэтому может применяться к относительно небольшим наборам данных; а также
• легко объяснить даже людям, не имеющим технического образования.

Я слышал, как старшие специалисты по данным, имеющие опыт работы с передовым искусственным интеллектом, хвалят линейную регрессию именно по этим причинам.

Однако одним из немногих недостатков линейной регрессии являются лежащие в ее основе строгие допущения. Эти допущения могут сделать модели линейной регрессии непригодными для использования в ряде очень распространенных обстоятельств.

К счастью, линейная регрессия существует так давно (если быть точным, с начала 19 века), что статистики давно нашли способ обойти любые нарушения допущений, когда они действительно происходят, при этом сохраняя при этом многие преимущества, связанные с линейной регрессией. регресс.

Тем не менее, чтобы справиться с ситуацией, когда нарушается одно или несколько предположений линейной регрессии, вам сначала необходимо определить, когда такое нарушение происходит. Это означает понимание предположений, лежащих в основе линейной регрессии.

Допущения линейной регрессии

Линейная регрессия основана на пяти ключевых предположениях, которые необходимо соблюдать, чтобы модель давала надежные прогнозы. Конкретно:

• Линейность: отношения между входными и выходными переменными линейны. То есть мы выражаем ожидаемое значение нашей выходной переменной Y как линейную комбинацию наших входных переменных, X₁,…, Xₖ: E (Y) = b₀ + b₁X₁ + b₂X₂ +… + bₖXₖ, где b₀, b₁, …, Bₖ обозначают параметры модели, которые необходимо подогнать.
• Отсутствие мультиколлинеарности: ни одна из входных переменных, X₁,…, Xₖ, не коррелирует друг с другом в сильной положительной или отрицательной степени;
• Нормальность: предполагается, что наблюдения выходной переменной Y основаны на одном и том же нормальном распределении. То есть Y ~ iid N (m, s²), где m = E (Y), а s² обозначает дисперсию Y;
• Гомоскедастичность: дисперсия s² выходной переменной Y считается постоянной, независимо от значений входных переменных; а также
• Независимость: предполагается, что наблюдения выходной переменной Y не зависят друг от друга. То есть в нашей выходной переменной нет автокорреляции.

Если ваши данные соответствуют всем этим предположениям, тогда хорошо для вас. Пойдите и установите свою модель линейной регрессии, а затем отдохните остаток дня.

Однако, цитируя Rolling Stones, «вы не всегда можете получить то, что хотите», и если вы имеете дело с чем-то более сложным, чем простой «игрушечный» набор данных, скорее всего, по крайней мере одно из этих предположений будет нарушено.

На этом этапе у вас есть два варианта: (а) дуться или (б) найти способ обойти любое нарушенное предположение.

Предполагая, что вы выберете вариант B, вот четыре способа обойти нарушение одного из предположений линейной регрессии.

Решение №1: удалите входные переменные, чтобы иметь дело с мультиколлинеарностью

Одна из самых простых проблем, которую нужно выявить и решить, - это мультиколлинеарность. Как правило, если две входные переменные имеют (абсолютный) коэффициент корреляции больше 0,8, велика вероятность возникновения проблемы мультиколлинеарности.

Это может затруднить интерпретацию коэффициентов модели и определение их статистической значимости, поскольку модель разделяет влияние того, что для всех практических целей является действительно одной переменной под двумя разными именами, на две отдельные входные переменные.

Например, рассмотрим набор данных, такой как этот из Kaggle, в котором указаны рост (в дюймах), пол и вес (в фунтах) 10 000 человек.

Подгонка модели к этому набору данных с использованием пакета Python statsmodels дает следующие подогнанные параметры:

import pandas as pd
import statsmodels.formula.api as smf
# Load data
multi_data = pd.read_csv('weight-height.csv')
# Fit linear regression model to data
multi_model = smf.ols(formula='Weight ~ Height + C(Gender)', data=multi_data).fit()
# Print summary statistics
multi_model.summary()

Предположим, мы создаем новую переменную Height_cm, которая дает рост каждого человека в сантиметрах (что явно на 100% коррелирует с исходной переменной Height, поскольку высота в сантиметрах = высота в дюймах / 2,54), и обновляем нашу модель, чтобы включить эта новая переменная (в дополнение к исходным переменным).

multi_model2 = smf.ols(formula='Weight ~ Height + C(Gender) + Height_cm', data=multi_data).fit()
multi_model2.summary()

Мы обнаружили, что коэффициент высоты из первой модели теперь разделен между Height и Height_cm (можно проверить, что 5,9769 = 5,1748 + 2,0373 / 2,54), что влияет на интерпретируемость коэффициентов обеих переменных.

Самое простое решение проблемы мультиколлинеарности - просто исключить из модели одну из сильно коррелированных входных переменных (не имеет значения, какая именно).

Решение # 2: используйте проектирование функций для устранения нелинейности

Линейная регрессия работает путем подбора (прямой) линии наилучшего соответствия вашим данным. Все это хорошо и хорошо в простых «игрушечных» примерах, которые они обычно дают вам на вводных занятиях по статистике, но что произойдет, если ваши данные будут выглядеть примерно так, как синие точки на этой диаграмме?

import numpy as np
import matplotlib.pyplot as plt
# Create dataset
np.random.seed(1)
x = np.random.uniform(-10, 10, 200)
y = x + x**2 + np.random.normal(0, 5, 200)
nl_data = pd.DataFrame({'X':x, 'Y':y})
# Fit linear regression to data
nl_model = smf.ols(formula='Y ~ X', data=nl_data).fit()
# Plot data and regression line
plt.scatter(x, y, alpha=0.5)
plt.xlabel('x')
plt.ylabel('y')
plt.plot(nl_data['X'], nl_model.predict(nl_data['X']), color = 'red')
plt.show()

Подгонка прямой линии регрессии к этим данным (красная линия) приведет к завышению прогноза выходной переменной (y) для значений входной переменной (x) в середине рассматриваемого диапазона и недооценке значений x на любом из крайних значений. диапазон.

Чтобы зафиксировать истинную структуру этих данных, нам действительно нужно подогнать полиномиальную кривую к нашим данным, но это невозможно сделать в рамках ограничений линейной регрессии, не так ли?

Что ж, на самом деле, разработав новые функции, которые являются функциями наших существующих входных переменных (включая мощности, журналы и произведения пар переменных), можно использовать линейную регрессию для соответствия чему-то другому, кроме прямая линия к нашим данным.

Например, в приведенном выше примере мы могли бы создать новую переменную z = x², а затем подогнать нашу модель линейной регрессии, используя как x, так и z в качестве входных переменных. Наша получившаяся модель будет иметь вид:

E(Y) = b₀+ b­₁x + b₂z = b₀+ b­₁x + b2x²

который по-прежнему является линейной комбинацией наших входных переменных x и z. Однако, поскольку z является просто функцией x, мы можем построить нашу подобранную линию регрессии в двух измерениях и получить что-то вроде этого:

# Create engineered variable
nl_data.loc[:, 'X2'] = nl_data['X'].apply(lambda x: x**2)
# Refit model to data, including the new variable
nl_model2 = smf.ols(formula='Y ~ X + X2', data=nl_data).fit()plt.scatter(x, y, alpha=0.5)
# Plot the data and the new regression line
plt.xlabel('x')
plt.ylabel('y')
plt.plot(nl_data['X'], nl_model2.predict(nl_data[['X', 'X2']]), color = 'red')
plt.show()

Кроме того, этот пример иллюстрирует важность построения графиков ваших данных перед попыткой подогнать под них модель, поскольку, визуализируя наши данные, вы можете получить представление о тех функциях, которые было бы полезно разработать.

Решение № 3: используйте преобразование переменных или обобщенную линейную модель, чтобы справиться с ненормальностью и гетероскедастичностью

Линейная регрессия предполагает, что наша выходная переменная получена из нормального распределения. То есть он симметричен, непрерывен и определен по всей числовой прямой.

На практике нарушения двух последних характеристик не представляют большого значения. Например, в приведенном выше примере роста и веса, несмотря на то, что человеческий вес обычно попадает только в относительно узкий диапазон и не может быть отрицательным, мы все же можем без особого беспокойства применить линейную регрессию к данным.

Однако, если наши данные искажены, это может привести к другим нарушениям наших предположений о линейной регрессии, если мы не исправим их.

Например, рассмотрим набор данных (gamma_data) с положительно искаженными значениями выходной переменной Y, как показано ниже:

# Create gamma dataset
np.random.seed(1)
x1 = np.random.uniform(-1, 1, 200)
x2 = np.random.uniform(-1, 1, 200)
mu = np.exp(1 + x1 + 2*x2 + np.random.randn())
y = np.random.gamma(shape = 2, scale = mu/2, size = 200)
gamma_data = pd.DataFrame({'X1':x1, 'X2':x2, 'Y':y})
# Plot data
plt.hist(gamma_data['Y'], bins=30)
plt.ylabel('Count')
plt.xlabel('Y')
plt.show()

Подгонка модели линейной регрессии к этому набору данных без предварительного преобразования нашей выходной переменной, а затем построения графика остатков против подобранных значений нашей выходной переменной дает нам следующий график остатков:

# Fit linear regression
non_norm_model = smf.ols(formula='Y ~ X1 + X2', data=gamma_data).fit()
# Calculate residuals
resid = gamma_data['Y'] - non_norm_model.predict(gamma_data[['X1', 'X2']])
# Plot residuals
plt.scatter(non_norm_model.predict(gamma_data[['X1', 'X2']]), resid, alpha=0.5)
plt.xlabel('fitted')
plt.ylabel('residual')
plt.show()

Разброс точек на этой диаграмме рассеяния увеличивается слева направо, что является явным признаком гетероскедастичности (то есть нарушения предположения о гомоскедастичности).

Чтобы справиться с этой проблемой, мы можем преобразовать нашу выходную переменную, чтобы сделать ее симметричной перед подгонкой нашей модели, например, взяв ее логарифм или квадратный корень, как показано ниже:

# Transform Y by taking the log of it
gamma_data.loc[:, 'log_y'] = gamma_data['Y'].apply(lambda x: np.log(x))
# Refit linear regression to transformed data
non_norm_model2 = smf.ols(formula='log_y ~ X1 + X2', data=gamma_data).fit()
# Calculate residuals
resid = gamma_data['log_y'] - non_norm_model2.predict(gamma_data[['X1', 'X2']])
# Plot residuals
plt.scatter(non_norm_model2.predict(gamma_data[['X1', 'X2']]), resid, alpha=0.5)
plt.xlabel('fitted')
plt.ylabel('residual')

Остатки этой модели теперь демонстрируют постоянный разброс, что указывает на гомоскедастичность.

В качестве альтернативы мы можем подобрать модель, специально разработанную для нестандартных данных, такую ​​как обобщенная линейная модель (GLM). Я подробно обсуждаю GLM в своей статье За пределами линейной регрессии: Введение в GLM.

В случае положительно искаженных данных гамма-модель GLM обычно является лучшим выбором.

# Fit GLM to data
gamma_model = sm.GLM(gamma_data['Y'], sm.add_constant(gamma_data[['X1', 'X2']]), family=sm.families.Gamma(sm.genmod.families.links.log)).fit()
# Calculate residuals
resid = gamma_model.resid_deviance
# Plot residuals
plt.scatter(gamma_model.predict(sm.add_constant(gamma_data[['X1', 'X2']]),), resid, alpha=0.5)
plt.xlabel('fitted')
plt.ylabel('residual')
plt.show()

Подбор гамма-модели GLM с лог-ссылкой на наши данные также решает проблему гетероскедастичности, выявленную ранее.

Решение # 4: используйте модель временных рядов для работы с автокорреляцией

Если вы имеете дело с данными, собираемыми через регулярные промежутки времени, то есть данными временных рядов (например, погодными или финансовыми данными), то, скорее всего, ваши данные автокоррелированы. Автокорреляция означает, что знание значения одного наблюдения вашей выходной переменной может дать вам представление о значении предыдущего или следующего наблюдения.

Например, в контексте данных о погоде, если бы вы знали максимальную температуру сегодня для вашего города, то вы имели бы приблизительное представление о максимальной температуре завтра в том же месте (поскольку температура обычно не слишком сильно меняется изо дня в день. ).

Самый простой способ определить, присутствует ли автокорреляция в ваших данных, - это построить график зависимости выходной переменной от времени.

Рассмотрим набор данных S&P 500, полученный от Kaggle. Если мы сосредоточимся на одной акции (с кодом фондовой биржи AAL) и построим график цены закрытия в зависимости от времени на протяжении всего набора данных, мы получим:

# Read in dataset
sandp_data = pd.read_csv('all_stocks_5yr.csv')
# Only keep data for AAL
sandp_data = sandp_data[sandp_data['Name'] == 'AAL']
# Transform date to datetime format
sandp_data['date'] = pd.to_datetime(sandp_data['date'])
# Sort data by date
sandp_data.sort_values(by = ['date'], inplace = True)
# Plot data
plt.plot(sandp_data['date'], sandp_data['close'])
plt.show()

Из этого графика мы видим, что в этих данных явно присутствует автокорреляция, с которой стандартная модель линейной регрессии не может справиться.

Если бы мы попытались подогнать к этим данным модель линейной регрессии, используя год и месяц в качестве входных переменных, мы бы получили показанную ниже красную линию, которая не идеально соответствовала нашим данным:

# Create year and month variables
sandp_data['Year'] = sandp_data['date'].map(lambda x: x.year)
sandp_data['Month'] = sandp_data['date'].map(lambda x: x.month)
# Fit linear regression to data
sandp_model = smf.ols(formula='close ~ Year + C(Month)', data=sandp_data).fit()
# Plot data with regression line overlay
plt.plot(sandp_data['date'], sandp_data['close'])
plt.plot(sandp_data['date'], sandp_model.predict(sandp_data[['Year', 'Month']]), color = 'red')
plt.show()

Чтобы обеспечить автокорреляцию в этом наборе данных, нам необходимо использовать модель, специально разработанную для работы с временными рядами, такую ​​как модель ARIMA.

Учитывая выходную переменную в момент времени t, Y (t), модель ARIMA (p, d, q) сначала берет разницу между последовательными наблюдениями нашей выходной переменной d несколько раз подряд, чтобы создать преобразованную выходную переменную y (t). Это делается для того, чтобы удалить из данных любые долгосрочные тенденции, такие как тенденция к увеличению на графике, показанном выше.

Например, если d = 1, y (t) = Y (t) - Y (t-1), если d = 2, y (t) = z (t) - z (t - 1), где z ( t) = Y (t) - Y (t - 1) и т. д.

После того, как мы взяли разность d наших данных, мы затем моделируем полученную преобразованную выходную переменную как линейную комбинацию p непосредственно предшествующих наблюдений y (t) и q непосредственно предшествующих остаточных значений модели (то есть разницы между фактические и прогнозируемые значения y (t)).

Существует много теории о том, как установить соответствующие значения для параметров p, d и q, что выходит за рамки этой статьи.

Для этого примера предположим, что d = 1, p = 5 и q = 0. То есть после однократного дифференцирования, чтобы учесть общую тенденцию к увеличению наших данных, разностные выходные значения могут быть смоделированы как линейная комбинация. из пяти непосредственно предшествующих (разностных) выходных значений.

Эту модель можно подогнать с помощью функции tsa.ARIMA в пакете Python statsmodels и получить следующую подобранную модель:

# Index dataset by date
sandp_data2 = sandp_data.set_index('date')
# Fit ARIMA model to the data
ts_model = sm.tsa.ARIMA(sandp_data2['close'], (5, 1, 0)).fit()
# Plot actual vs fitted values
ts_model.plot_predict(dynamic=False)

Эта модель очень хорошо отслеживает наши данные, при этом прогноз и фактические данные отслеживаются настолько точно, что различить их практически невозможно.

Кроме того, следует отметить, что такая высокая производительность частично связана с тем, что мы делаем прогнозы в пределах выборки, которую мы использовали для соответствия нашим данным, и такая производительность вряд ли будет продолжаться вне выборки.

Принцип бритвы Оккама гласит, что если существует два объяснения происшествия, то обычно правильным является то, которое требует наименьшего количества предположений.

В контексте науки о данных это означает, что чем сложнее модель, тем меньше вероятность того, что она объяснит ваши данные, или с практической точки зрения, если вы пытаетесь смоделировать набор данных, то сначала начните с простейшей модели и только переходите к рассмотрению более сложных моделей, если более простые модели оказываются непригодными.

Чаще всего для задачи регрессии самой простой моделью, с которой вы можете начать, является модель линейной регрессии. Однако во многих случаях нарушение одного или нескольких предположений о строгой линейной регрессии может сделать использование этой модели неприемлемым.

В этой статье мы предложили несколько обходных путей для предположений линейной регрессии, которые могут позволить вам продолжать использовать эту универсальную и легко понимаемую модель (или связанные модели, такие как GLM или модели временных рядов), прежде чем перейти к ресурсоемкому черному. коробки, такие как нейронные сети.

Возможно, модели линейной регрессии и не круты, но у них есть проверенный послужной список, и для специалиста по данным это все, что действительно должно иметь значение.