Прогноз проката велосипедов с использованием регрессии

В этом посте я сравниваю модели регрессии на наборе данных совместного использования велосипедов из Репозитория машинного обучения UCI. Цель задачи регрессии состоит в том, чтобы предсказать количество велосипедов, арендованных в каждый час. Проект был написан на Python с использованием Jupyter Notebooks и доступен по ссылке в конце этого поста.

Содержание

Набор данных
Исследование и предварительная обработка данных
Показатели оценки
Линейная регрессия
Полиномиальная регрессия
Деревья регрессии
Сравнение результатов

Набор данных

Данные состоят из характеристик окружающей среды, таких как температура,
влажность и количество осадков, регистрируемых каждый час в течение годичного периода с декабря 2017 г. по ноябрь 2018 г. Существует 8760 наблюдений и 13 характеристик, за исключением целевого признака, который Rented Bike Count. есть. в наборе данных нет нулевых значений. Большинство признаков являются непрерывными, но есть также три категориальных признака. Первые 5 строк необработанного набора данных можно увидеть ниже.

Исследование данных и предварительная обработка

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

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

Сразу видно, что Dew point temperature и Temperature сильно коррелированы (это может показаться очевидным, но рекомендуется подтверждать эти предположения количественным методом!) Учитывая коллинеарность между этими признаками, я оставляю только Temperature, потому что он имеет более сильную корреляцию с целевым признаком.

Тепловая карта сравнивает только непрерывные функции, как мы можем оценить важность нечисловых функций? Seasons, Holiday и Functioning Day — это категориальные признаки, не имеющие количественного значения. Нам нужно преобразовать их в числовые переменные. При преобразовании типов данных нужно быть осторожным! Заманчиво применить простое преобразование, при котором каждая категория заменяется числом. Времена года Весна, Лето, Осень и Зима, например, могут быть закодированы как 1, 2, 3 и 4. Однако числа подразумевают соотношение количества и/или порядка. Это кодирование ошибочно подразумевает, что Весна «меньше» категории «Зима» или предшествует ей*. Мы не хотим сообщать нашей модели, что есть отношения там, где их нет, поэтому вместо этого мы можем использовать метод преобразования, называемый двоичным или горячим кодированием.

*это может показаться нелогичным, поскольку времена года можно интерпретировать как порядковые. Однако без дополнительной информации нет неявного порядка, т. е. весна может быть как до, так и после зимы (весна 2017 года, зима 2017 года и весна 2018 года).

Двоичное и горячее кодирование

Оба эти метода представляют собой способы преобразования нечисловых переменных в последовательность нулей и единиц; разница заключается в количестве новых функций (столбцов), которые вводятся в набор данных. В качестве примера возьмем категориальную переменную Seasons, для которой есть четыре возможных значения: Весна, Лето, Осень и Зима.

Горячее кодирование

Этот метод заменяет номинальное значение (например, «Зима») бинарным вектором n-элементов, где n — количество возможных значений в категории. Ниже показан пример, в котором столбец Seasons был заменен n=4 фиктивными столбцами.

Двоичная кодировка

Этот метод больше подходит, когда в категории имеется большое количество возможных значений, поскольку он уменьшает количество фиктивных столбцов. Категориальный признак сначала преобразуется в числовые значения с помощью порядкового кодировщика, а затем в двоичное число. В нашем примере Зима, Осень, Весна и Лето сначала преобразуются в 1, 2, 3 и 4, а затем кодируются как 001, 010, 011 и 100 соответственно.

Для этой задачи у нас есть три категориальные переменные для преобразования: Seasons, Holiday и Functioning Day. Я использую метод горячего кодирования, потому что два из них уже являются бинарными функциями, а Seasons содержит только четыре возможных значения.

Функция Date содержит интервальные данные, т. е. существует измеримое расстояние (во времени) между двумя значениями, но нет фиксированного нуля. Мы можем разделить даты на три непрерывных элемента: Day, Month и Year. Без дальнейшей предварительной обработки эти признаки показывают слабую корреляцию с целевым признаком. Day, в частности, практически не имеет корреляции (коэффициент корреляции Пирсона 0,022). В этом нет ничего удивительного, если задуматься о том, что это означает на самом деле: почему может существовать связь между количеством арендованных велосипедов и днем месяца? Лучше задать вопрос: относятся ли дни недели (будние или выходные) к количеству арендованных велосипедов? Я проверяю корреляции, используя функцию datetime в pandas, чтобы создать дополнительную функцию is_Weekend из столбца даты, которая принимает значения True или False. Корреляция с целевым объектом очень слабая, поэтому я удаляю и Day, и is_Weekend. Наконец, я думаю, не будет ли излишним сохранять обе функции Month и Season. Другими словами, возможно, больше велосипедов берут напрокат в более теплые месяцы. Это было бы похоже на вопрос «Больше велосипедов берут напрокат весной/летом?». Давайте еще раз проверим, верно ли это предположение.

Мы можем видеть на изображении выше, сравнивая эти функции с Rented Bike Count, что аналогичная тенденция передается как Month, так и Seasons. В принципе, мы могли бы выбрать любую функцию для использования в нашей модели, однако есть преимущества и недостатки использования более или менее детальных данных.

Группировка месяцев по сезонам создает группированные данные, что снижает степень детализации. Биннинг помогает уменьшить влияние выбросов и может улучшить способность модели обобщать невидимые данные (а это то, что нам нужно!). С другой стороны, мы можем потерять потенциально полезную информацию из нашего набора данных, игнорируя столбец Month. Кроме того, данные, которые группируются в соответствии с заданными пользователем пороговыми значениями, могут создавать проблемы из-за произвольного группирования и систематической ошибки. Например, деление осадков на легкие, средние и сильные категории довольно субъективно. В нашем случае месяцы и времена года — это предопределенные категории, которые всегда сгруппированы одинаково, поэтому, к счастью, мы вообще избегаем этой проблемы. Для моей модели я выбираю меньшую степень детализации и отказываюсь от функции Month в пользу Seasons.

Последним шагом в предварительной обработке является отделение целевого объекта от остальных данных и сохранение значений в массиве. Это дает нам массив Xразмеров n наблюдений по d объектам и одномерный массив y значений n целевых ярлыков.

Метрики оценки

Важно разделить данные на обучающий набор и тестовый набор, чтобы каждая модель оценивалась по ее производительности на невидимых данных. Это важно, потому что эффективность модели заключается в ее способности обобщать невидимые данные, а не изучать закономерности, характерные для обучающих данных. Для моделей, в которых необходимо выбрать гиперпараметр, также важно разделить данные на третий набор, известный как набор проверки. Гиперпараметры — это параметры, которые не извлекаются из данных, а скорее параметры, которые необходимо настроить при построении модели. Использование проверочного набора позволяет оптимизировать гиперпараметры, при этом гарантируя, что некоторые данные будут скрыты, чтобы окончательную производительность модели можно было оценить на невидимых данных. Я буду использовать метод оценки, называемый k-кратной перекрестной проверкой, для оптимизации гиперпараметров.

k-кратная перекрестная проверка

При k-кратной перекрестной проверке обучающая выборка далее разбивается на k кратностей без замены. k - 1 кратность используется для обучения модели с гиперпараметром λ, а оставшаяся кратность используется в качестве проверочного набора. Затем эта процедура повторяется k раз, чтобы получить k оценок производительности. Среднее значение этих k оценок производительности вычисляется для получения окончательной оценки производительности для гиперпараметра λ. Этот процесс повторяется для каждого значения λ, и значение с наилучшей средней оценкой производительности выбирается для использования в окончательной оценке модели на тестовом наборе. Преимущество использования k-кратной перекрестной проверки заключается в том, что она помогает гарантировать, что модель может обобщаться, уменьшая вероятность переобучения при начальном разделении обучения/проверки. Недостатком этого подхода являются вычислительные затраты на повторение обучения и тестирования для каждой складки.

Показатели эффективности

Я буду измерять производительность каждой модели по сравнению с базовым уровнем, используя модуль sklearn DummyRegressor, который в качестве своей стратегии по умолчанию всегда предсказывает среднее значение целей обучения. Базовый уровень действует как значимая точка отсчета, которая может помочь нам интерпретировать производительность моделей. Например, если базовая линия делает прогнозы с точностью 90 %, любое значительное улучшение производительности может не оправдывать компромисс между вычислительными затратами и точностью более сложной модели. Чтобы количественно оценить производительность каждой модели, я использую два измерения: среднеквадратичную ошибку и оценку R².

Среднеквадратическая ошибка (MSE) — это мера разницы между фактическими значениями y и значениями, предсказанными моделью. Эта разница также известна как остаток. Мы хотим минимизировать MSE. Вычисление MSE включает в себя возведение остатков в квадрат, чтобы положительные и отрицательные значения не компенсировали друг друга. Однако важно признать, что это означает, что это не является надежной мерой для выбросов в наборе данных из-за влияния больших остатков.
Оценка R² является мерой доли дисперсии, которая объясняется моделью. Он измеряет, насколько хорошо модель соответствует данным, по сравнению с простым получением среднего значения целевых значений. Числитель оценки R² — это MSE (см. здесь). Если MSE=0, то это будет означать, что наша модель всегда предсказывает фактические значения y, и в этом случае R² = 1. Это означает, что мы можем интерпретировать значения R² ближе к 1 как признак повышенной производительности и наоборот. Знаменатель оценки R² выступает в качестве базовой линии, относительно которой мы можем измерить качество подгонки модели. Полезно использовать метрику R² в дополнение к MSE из-за ее интерпретируемости, поскольку она не зависит от масштаба данных.

Линейная регрессия

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

где aᵢ — набор коэффициентов, xᵢ — d-мерные входные векторы и b является термином смещения. Мы можем сформулировать это как линейную систему уравнений и найти решение, используя метод наименьших квадратов, чтобы получить параметры модели aᵢ и b. Используя библиотеку sklearn, мы можем реализовать эту модель и выполнить перекрестную проверку k-fold с помощью нескольких простых строк кода:

from sklearn.model_selection import train_test_split
from sklearn.linear_model import LinearRegression
from sklearn.model_selection import KFold, cross_val_score

# split data into train and test sets
Xtrain, Xtest, Ytrain, Ytest = train_test_split(X, target, test_size=0.2, random_state=0)

# specify number of validation folds (here k=5)
kf = KFold(n_splits=5, shuffle=True, random_state=0) 

# instantiate the model
lin_reg = LinearRegression() 

# compute RMSE and R2 score
lin_reg_mse = -np.mean(cross_val_score(lin_reg, X, target, cv=kf, scoring='neg_mean_squared_error'))
lin_reg_r2 = np.mean(cross_val_score(lin_reg, X, target, cv=kf, scoring='r2'))

Примечание. Метрика регрессии neg_mean_squared_error, совместимая с функцией sklearn cross_val_score, возвращает отрицательное среднеквадратичное отклонение.

Модель возвращает MSE = 187930 (!) и показатель R² 0,549. Эти результаты легче интерпретировать с помощью визуализации — если бы прогнозы были точными на 100%, все значения располагались бы вдоль пунктирной линии.

Полиномиальная регрессия

Полиномиальная регрессия — это адаптация метода линейной регрессии, которая включает добавление дополнительных функций во входное пространство. Эти дополнительные функции представляют собой члены более высокого порядка, построенные из исходных входных функций. Это позволяет нам моделировать нелинейную связь между входными функциями и целевыми значениями. Другими словами, мы моделируем связь между входным вектором x и выходной меткой y как полиномn-й степени. Например, в одномерном пространстве признаков мы могли бы написать эту систему уравнений с добавлением членов более высокого порядка, как показано ниже:

Хотя полиномиальная регрессия означает, что мы больше не подгоняем прямую линию к данным, она по-прежнему считается линейной моделью, поскольку коэффициенты полиномиальных членов имеют линейную зависимость. При построении дополнительных функций важно, чтобы мы не подгоняли модель к обучающим данным. Например, полином более высокого порядка может идеально соответствовать обучающим данным, но плохо обобщать невидимые данные. Чтобы найти оптимальное количество полиномиальных признаков, мы можем использовать модуль sklearn PolynomialFeatures. Степень полинома может быть оптимизирована с помощью метода перекрестной проверки k-кратного размера перед подгонкой модели к обучающим данным с использованием того же модуля линейной регрессии, который мы использовали ранее.

Графики выше сравнивают производительность модели на обучающих данных и на проверочном наборе. Это помогает нам настроить гиперпараметр (количество дополнительных членов для добавления в модель), потому что мы можем легко увидеть, когда модели начинают страдать от проблемы переобучения. Показатели MSE и R² для проверочного набора резко расходятся с данными обучения для n≥2. Поэтому я тестирую окончательную модель, используя n=2 дополнительных термина, т. е. мы добавляем термин x² в пространство признаков. Модель дает более низкий, хотя и не впечатляющий, MSE 127311 и улучшенный показатель R² 0,696.

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

Деревья регрессии

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

Например, мы разделяем данные на области M R₁, R₂ , … Rₘ и пусть cₘ равно среднему значению целевых переменных y в области Rₘ. Затем для каждой функции j и возможных точек разделения s мы определяем пару областей разделения:

чтобы найти функцию j и точку разделения s, которая решает:

Этот процесс повторяется, так что пара (j, s), которая минимизирует RSS для каждой результирующей области Rᵢ, выбирается в качестве следующего узла в дереве. Оптимизация размера дерева — важный шаг, позволяющий избежать переобучения. Распространенный метод нахождения оптимального размера дерева называется сокращением стоимости и сложности. Вместо восходящего построения дерева, в котором вы ограничиваете гиперпараметры, такие как глубина дерева, отсечение по стоимости и сложности включает в себя построение полного дерева, а затем работу в обратном направлении для удаления листьев и узлов, которые предоставляют наименьшую информацию. Этот метод помогает устранить избыточную сложность и превосходит усеченные модели. Стоимость-сложность дерева T определяется как:

где |T| — количество листовых узлов, Qₘ(T) — сумма RSS по всем листовым узлам дерева, а α — член регуляризации, который необходимо изучить. Используя функцию sklearn cost_complexity_pruning_path, мы можем найти значение α и, таким образом, оптимальное поддерево Tₐ, которое минимизирует RSS для проверочных данных. Мы будем использовать тот же метод перекрестной проверки, что и выше. Результаты использования перекрестной проверки для выбора оптимального значения α показаны ниже.

Когда α= 0, мы получаем полное дерево, используя параметры по умолчанию модуля дерева регрессии sklearn. Производительность на обучающих данных почти идеальна, когда α = 0, но значительно хуже на проверочных данных. Поскольку α стремится к бесконечности, мы получаем наименьшее возможное поддерево, дерево с одним узлом, и производительность как для обучающих, так и для проверочных данных низкая, учитывая, что
в дереве с одним узлом нет принимаются решения о разделении данных. Мы используем значение α, которое максимизирует оценку проверки R², чтобы обучить окончательную модель, которая дает нам дерево глубины 13 с 72 листовыми узлами.

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

Сравнение результатов

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

Код Python, написанный в Jupyter Notebooks для этого проекта, доступен здесь: СКОРО