Руководство по PyTorch для новичков: создаём модель множественной регрессии с нуля

Практический PyTorch: строим трёхслойную нейросеть для множественной регрессии

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

Про машинное обучение обычно говорили в контексте scikit-learn, XGBoost и им подобных, тогда как PyTorch и TensorFlow доминировали на сцене, когда речь заходила о глубоком обучении.

После взрывного роста интереса к ИИ я всё чаще вижу, что PyTorch заметно опережает TensorFlow по популярности. Оба фреймворка очень мощные и позволяют дата-сайентистам решать самые разные задачи, включая обработку естественного языка, что вновь подогрело интерес к глубокому обучению.

Однако в этой статье я не собираюсь говорить про NLP. Вместо этого мы разберём задачу множественной линейной регрессии с двумя целями:

• Изучить, как создать модель с использованием PyTorch.

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

Приступим.

Подготовка данных

Давайте я избавлю вас от ещё одной вычурной формальной дефиниции линейной регрессии. Скорее всего, вы уже видели её десятки раз в бесконечных туториалах по всему Интернету. Достаточно сказать так: когда у вас есть переменная Y, которую вы хотите предсказывать, и другая переменная X, которая может объяснить изменение Y с помощью прямой линии, — в сущности, это и есть линейная регрессия.

Набор данных

В этом упражнении мы будем использовать набор данных Abalone [1].

Nash, W., Sellers, T., Talbot, S., Cawthorn, A., & Ford, W. (1994). Abalone [Dataset]. UCI Machine Learning Repository. https://doi.org/10.24432/C55C7W

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

Итак, давайте загрузим данные. Дополнительно мы применим one-hot-кодирование к переменной Sex, так как это единственная категориальная переменная в наборе.

# Загрузка данных from ucimlrepo import fetch_ucirepo import pandas as pd import numpy as np import matplotlib.pyplot as plt import seaborn as sns sns.set_style('darkgrid') from feature_engine.encoding import OneHotEncoder  # Загружаем набор данных abalone = fetch_ucirepo(id=1)  # данные (в виде pandas DataFrame) X = abalone.data.features y = abalone.data.targets  # One-Hot кодирование признака Sex ohe = OneHotEncoder(variables=['Sex']) X = ohe.fit_transform(X)  # Объединяем признаки и целевую переменную в один DataFrame df = pd.concat([X, y], axis=1)

Вот набор данных.

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

Исследуем данные

Первые шаги, которые я обычно выполняю при разведочном анализе данных, такие:

1. Проверить распределение целевой переменной.

# Анализ распределения целевой переменной plt.hist(y) plt.title('Rings [Target Variable] Distribution')

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

2. Посмотреть статистическое описание.

Статистика позволяет заметить важные характеристики — среднее значение, стандартное отклонение, а также быстро увидеть аномалии в минимальных или максимальных значениях. Объясняющие переменные выглядят вполне корректно: находятся в небольшом диапазоне и на одной шкале. Целевая переменная (Rings) — на другой шкале.

# Статистическое описание данных df.describe()

Теперь проверим корреляции.

# Анализ корреляций (df  .drop(['Sex_M', 'Sex_I', 'Sex_F'], axis=1)  .corr()  .style  .background_gradient(cmap='coolwarm') )

Объясняющие переменные имеют среднюю и сильную корреляцию с Rings. Также видно, что есть мультиколлинеарность между Whole_weight и Shucked_weight, Viscera_weight и Shell_weight. Length и Diameter тоже коллинеарны. Позже можно проверить, улучшится ли модель, если их убрать.

sns.pairplot(df);

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

• Нарушено предположение гомоскедастичности. Это означает, что дисперсия ошибки непостоянна.

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

• Переменная Height имеет минимум два явных выброса при Height > 0.3.

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

Есть ещё один простой способ изучить данные — построить графики связи переменных, сгруппированные по Sex.

Переменная Diameter выглядит наиболее линейной при Sex = I, но на этом всё.

# Создаём FacetGrid с диаграммами рассеяния sns.lmplot(x="Diameter", y="Rings", hue="Sex", col="Sex", order=2, data=df)

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

# Создаём FacetGrid с диаграммами рассеяния sns.lmplot(x="Shell_weight", y="Rings", hue="Sex", col="Sex", data=df)

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

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

Моделирование: используем Scikit-Learn

Давайте запустим модель из sklearn и оценим её по RMSE (корню из среднеквадратичной ошибки).

from sklearn.linear_model import LinearRegression from sklearn.metrics import root_mean_squared_error  df2 = df.query('Height < 0.3 and Rings > 2 ').copy() X = df2.drop(['Rings'], axis=1) y = np.log(df2['Rings'])  lr = LinearRegression() lr.fit(X, y)  predictions = lr.predict(X)  df2['Predictions'] = np.exp(predictions) print(root_mean_squared_error(df2['Rings'], df2['Predictions']))

2.2383762717104916

Если посмотреть на первые строки таблицы (head), можно убедиться, что модели сложно даются оценки для больших значений целевой переменной (например, строки 0, 6, 7 и 9).

Шаг назад: пробуем другие трансформации

Итак, что мы можем сделать дальше?

Скорее всего, стоит удалить ещё часть выбросов и попробовать снова. Попробуем применить несупервизируемый алгоритм, чтобы найти дополнительные выбросы. Мы воспользуемся методом Local Outlier Factor и отбросим 5% наблюдений, признанных выбросами.

Заодно уберём мультиколлинеарность, удалив признаки Whole_weight и Length.

from sklearn.neighbors import LocalOutlierFactor from sklearn.preprocessing import StandardScaler from sklearn.pipeline import Pipeline  # загрузим набор данных abalone = fetch_ucirepo(id=1)  # данные (как pandas DataFrame) X = abalone.data.features y = abalone.data.targets  # One Hot Encode для признака Sex ohe = OneHotEncoder(variables=['Sex']) X = ohe.fit_transform(X)  # удалим Whole_weight и Length (мультиколлинеарность) X.drop(['Whole_weight', 'Length'], axis=1, inplace=True)  # посмотрим на датафрейм df = pd.concat([X,y], axis=1)  # создадим Pipeline для масштабирования данных и поиска выбросов методом Local Outlier Factor (на основе k ближайших соседей) steps = [ ('scale', StandardScaler()), ('LOF', LocalOutlierFactor(contamination=0.05)) ] # обучим модель и получим предсказания (метки выбросов) outliers = Pipeline(steps).fit_predict(X)  # добавим столбец с метками выбросов df['outliers'] = outliers  # моделирование df2 = df.query('Height < 0.3 and Rings > 2 and outliers != -1').copy() X = df2.drop(['Rings', 'outliers'], axis=1) y = np.log(df2['Rings'])  lr = LinearRegression() lr.fit(X, y)  predictions = lr.predict(X)  df2['Predictions'] = np.exp(predictions) print(root_mean_squared_error(df2['Rings'], df2['Predictions']))

2.238174395913869

Результат тот же. Хм…

Ладно, мы можем продолжать экспериментировать с переменными и заниматься генерацией признаков (feature engineering) — и увидим небольшие улучшения. Например, если добавить квадраты Height, Diameter и Shell_weight. В сочетании с обработкой выбросов это снизит RMSE до 2.196.

# признаки второй степени X['Diameter_2'] = X['Diameter'] ** 2 X['Height_2'] = X['Height'] ** 2 X['Shell_2'] = X['Shell_weight'] ** 2

Важно отметить: каждый новый признак, добавляемый в модель линейной регрессии, влияет на R? и иногда искусственно «раздувает» его, создавая ложное ощущение, что модель стала лучше, хотя это не так. В данном случае модель действительно улучшается, потому что мы добавляем в неё нелинейные компоненты за счёт признаков второй степени. Это можно подтвердить, посчитав скорректированный R?: он вырос с 0.495 до 0.517.

# Скорректированный R? from sklearn.metrics import r2_score   r2 = r2_score(df2['Rings'], df2['Predictions']) n= df2.shape[0] p = df2.shape[1] - 1 adj_r2 = 1 - (1 - r2) * (n - 1) / (n - p - 1) print(f'R?: {r2}') print(f'Adjusted R?: {adj_r2}')

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

Моделирование: использование PyTorch

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

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

import torch import torch.nn as nn import torch.optim as optim from torch.utils.data import DataLoader, TensorDataset

Нужно подготовить данные для моделирования в PyTorch. Здесь нам потребуется преобразовать данные в формат, который PyTorch «понимает», так как он не работает напрямую с pandas DataFrame.

• Возьмём тот же датафрейм, что и для базовой модели.

• Разделим X и Y.

• Преобразуем целевую переменную Y с помощью логарифма.

• Преобразуем обе части в массивы numpy, так как PyTorch не принимает DataFrame.

df2 = df.query('Height < 0.3 and Rings > 2 and outliers != -1').copy() X = df2.drop(['Rings', 'outliers'], axis=1) y = np.log(df2[['Rings']])   # Преобразуем X и Y в Numpy X = X.to_numpy() y = y.to_numpy()

Далее, используя TensorDataset, мы превращаем X и Y в тензоры и выводим пример.

# Подготовка с помощью TensorDataset # TensorDataset помогает преобразовать набор данных в объект Tensor dataset = TensorDataset(torch.tensor(X).float(), torch.tensor(y).float())   input_sample, label_sample = dataset[0] print(f'** Input sample: {input_sample}, 
** Label sample: {label_sample}')

** Input sample: tensor([0.3650, 0.0950, 0.2245, 0.1010, 0.1500, 1.0000,  0.0000, 0.0000, 0.1332, 0.0090, 0.0225]),  ** Label sample: tensor([2.7081])

Затем, с помощью функции DataLoader, мы можем разбить данные на батчи. Это означает, что нейросеть будет обрабатывать данные порциями размера batch_size.

# Теперь воспользуемся DataLoader batch_size = 500 dataloader = DataLoader(dataset, batch_size=batch_size, shuffle=True)

Модели в PyTorch лучше всего определять в виде классов.

• Класс наследуется от nn.Module — базового класса PyTorch для нейросетей.

• В методе init мы задаём слои модели. 

  • Вызов super().__init__() гарантирует, что класс будет вести себя как torch-объект.

• Метод forward описывает, что происходит с входными данными при прохождении через модель.

Здесь мы пропускаем их через линейные слои, определённые в init, и используем функции активации ReLU, чтобы добавить нелинейность при прямом проходе.

# 2. Создаём класс class AbaloneModel(nn.Module):   def __init__(self):     super().__init__()     self.linear1 = nn.Linear(in_features=X.shape[1], out_features=128)     self.linear2 = nn.Linear(128, 64)     self.linear3 = nn.Linear(64, 32)     self.linear4 = nn.Linear(32, 1)     def forward(self, x):     x = self.linear1(x)     x = nn.functional.relu(x)     x = self.linear2(x)     x = nn.functional.relu(x)     x = self.linear3(x)     x = nn.functional.relu(x)     x = self.linear4(x)     return x   # Создаём экземпляр модели model = AbaloneModel()

Далее попробуем модель в деле, используя скрипт, который имитирует Random Search по гиперпараметрам.

• Создадим функцию ошибки для оценки модели.

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

• Зададим диапазон для шага обучения: будем использовать степени от ?2 до ?4 (то есть от 0.01 до 0.0001).

• Зададим диапазон для параметра momentum от 0.9 до 0.99.

• Получим данные.

• Обнулим градиенты, чтобы очистить результаты вычислений с предыдущих шагов.

• Обучим модель.

• Посчитаем функцию потерь и запишем показатели для лучшей модели.

• Посчитаем градиенты весов и смещений на обратном проходе.

• Повторим эти шаги N раз и в конце выведем лучшую модель.

# Среднеквадратичная ошибка (MSE) — стандартная функция потерь для задач регрессии criterion = nn.MSELoss()   # Random Search values = [] best_loss = 999 for idx in range(1000):   # Случайным образом выбираем показатель степени для шага обучения в диапазоне от 2 до 4   factor = np.random.uniform(2,5)   lr = 10 ** -factor     # Случайным образом выбираем momentum в диапазоне от 0.90 до 0.99   momentum = np.random.uniform(0.90, 0.99)     # 1. Получаем данные   feature, target = dataset[:]   # 2. Обнуляем градиенты, чтобы очистить результаты прошлых обратных проходов   optimizer = optim.SGD(model.parameters(), lr=lr, momentum=momentum)   optimizer.zero_grad()   # 3. Прямой проход: считаем предсказание   y_pred = model(feature)   # 4. Считаем функцию потерь   loss = criterion(y_pred, target)   # 4.1. Регистрируем лучшую (минимальную) ошибку   if loss < best_loss:     best_loss = loss     best_lr = lr     best_momentum = momentum     best_idx = idx     # 5. Обратный проход: считаем градиенты функции потерь по W и b   loss.backward()   # 6. Обновляем параметры: корректируем W и b с учётом вычисленных градиентов   optimizer.step()   values.append([idx, lr, momentum, loss])   print(f'n: {idx},lr: {lr}, momentum: {momentum}, loss: {loss}')

n: 999,lr: 0.004782946959508322, momentum: 0.9801209929050066, loss: 0.06135804206132889

Когда мы нашли лучший шаг обучения и momentum, можно двигаться дальше.

# --- 3. Функция потерь и оптимизатор ---   # Среднеквадратичная ошибка (MSE) — стандартная функция потерь для задач регрессии criterion = nn.MSELoss()   # Стохастический градиентный спуск (SGD) с небольшим шагом обучения (lr) optimizer = optim.SGD(model.parameters(), lr=0.004, momentum=0.98)

Дальше мы переобучим модель, используя те же шаги, что и раньше, но теперь фиксируя шаг обучения и momentum.

Обучение модели в PyTorch требует более развёрнутого кода, чем простой вызов метода fit() в Scikit-Learn. Но это не проблема: структура почти всегда будет такой:

1. перевести модель в режим обучения model.train();

2. создать цикл по числу итераций (эпох), которое вы хотите выполнить;

3. в каждой итерации обнулять градиенты через optimizer.zero_grad();

4. получать батчи данных из dataloader;

5. считать предсказания y_pred = model(X);

6. вычислить функцию потерь loss = criterion(y_pred, target);

7. выполнить обратный проход loss.backward(), чтобы посчитать градиенты весов и смещений;

8. обновить веса и смещения с помощью optimizer.step().

Мы будем обучать эту модель в течение 1000 эпох (итераций). Дополнительно добавим шаг по сохранению лучшей версии модели (с минимальной ошибкой), чтобы в конце использовать именно её.

# 4. Обучение torch.manual_seed(42) NUM_EPOCHS = 1001 loss_history = [] best_loss = 999   # Переводим модель в режим обучения model.train()   for epoch in range(NUM_EPOCHS):   for data in dataloader:       # 1. Получаем данные     feature, target = data       # 2. Обнуляем градиенты перед обратным проходом     optimizer.zero_grad()       # 3. Прямой проход: считаем предсказание     y_pred = model(feature)       # 4. Считаем функцию потерь     loss = criterion(y_pred, target)     loss_history.append(loss)       # Сохраняем лучшую модель     if loss < best_loss:       best_loss = loss       best_model_state = model.state_dict()  # сохраняем состояние лучшей модели       # 5. Обратный проход: считаем градиенты функции потерь по W и b     loss.backward()       # 6. Обновляем параметры: корректируем W и b с учётом вычисленных градиентов     optimizer.step()       # Загружаем лучшее состояние модели перед тем, как использовать её для предсказаний     model.load_state_dict(best_model_state)     # Печатаем прогресс каждые 200 эпох   if epoch % 200 == 0:     print(epoch, loss.item())     print(f'Best Loss: {best_loss}')

0 0.061786893755197525 Best Loss: 0.06033024191856384 200 0.036817338317632675 Best Loss: 0.03243456035852432 400 0.03307393565773964 Best Loss: 0.03077109158039093 600 0.032522525638341904 Best Loss: 0.030613820999860764 800 0.03488151729106903 Best Loss: 0.029514113441109657 1000 0.0369877889752388 Best Loss: 0.029514113441109657

Отлично. Модель обучена. Теперь пора её оценить.

Оценка качества

Давайте проверим, стала ли эта модель лучше по сравнению с обычной линейной регрессией. Для этого я переведу модель в режим оценки с помощью model.eval(), чтобы PyTorch понял, что нужно изменить поведение: из режима обучения перейти в режим инференса. Например, в этом режиме отключаются слои Dropout и batch normalization переводится в режим оценки.

# Получаем признаки и целевую переменную features, targets = dataset[:]   # Получаем предсказания model.eval() with torch.no_grad():   predictions = model(features)   # Добавляем предсказания в датафрейм df2['Predictions'] = np.exp(predictions.detach().numpy())   # RMSE print(root_mean_squared_error(df2['Rings'], df2['Predictions'])) 2.1108551025390625

2.1108551025390625

Улучшение получилось скромным — примерно 4%.

Давайте посмотрим на некоторые предсказания обеих моделей.

Обе модели дают очень похожие результаты. Им становится сложнее справляться с предсказаниями по мере роста значения Rings. Причина — в «конусной» форме распределения целевой переменной.

Если немного об этом подумать:

• По мере увеличения числа колец (Rings) растёт и разброс, вносимый объясняющими переменными.

• Морское ушко с 15 кольцами может иметь гораздо более широкий диапазон значений признаков, чем экземпляр с 4 кольцами.

• Это сбивает модель с толку, потому что ей приходится проводить одну-единственную линию посередине данных, которые по сути не особенно линейны.

Заключение

В этом проекте мы разобрали много всего:

• Как исследовать данные.

• Как понять, подходит ли линейная модель для конкретной задачи.

• Как создать модель на PyTorch для задачи многомерной линейной регрессии.

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

Мы попытались использовать глубокое обучение, чтобы улучшить результат, но всей этой мощности оказалось недостаточно, чтобы заметно снизить ошибку. Я бы, вероятно, выбрал модель на Scikit-Learn: она проще и лучше интерпретируется.

Другой путь улучшения качества — собрать кастомный ансамбль, например Random Forest + линейная регрессия. Но эту часть я оставляю вам, если захотите поэкспериментировать.

Научиться работать с важнейшими моделями машинного обучения, NLP, DL, рекомендательными системами на практике с реальными данными можно на курсе "Machine Learning. Professional". Пройдите вступительный тест, чтобы узнать, подойдет ли вам программа курса.

А если не хватает «живых» кейсов и хочется увидеть, как ML ведёт себя в продакшене, приходите на демо-уроки:

• 11 декабря 20:00 — Ускоряем работу с видеопотоками при помощи TorchCodec. Записаться

• 17 декабря 18:00 — Оптимизируем построение модели через Pipeline. Записаться

• 22 декабря 18:00 — ML как основа современного AI. Записаться