Подбор модели, корреляция, p-значение, t-статистика, доверительные интервалы и визуализация в Python
Линейная регрессия — это наиболее часто используемая регрессионная модель. Причина в том, что он прост в использовании, он может давать полезную информацию и его легко понять. В этой статье мы обсудим подгонку модели линейной регрессии к данным, вывод из нее и некоторую полезную визуализацию.
Используемые инструменты:
- Язык программирования Python и несколько его популярных библиотек. Если вы не знаете всех этих библиотек, вы все равно сможете прочитать эту статью и понять концепцию.
- Среда Jupyter Notebook.
Выбор функции
По моему опыту, лучший способ научиться — это использовать пример. Я буду использовать набор данных и продолжу объяснять процесс подгонки модели к данным и выводить информацию. Я использую набор данных обследования под названием NHANES dataset. Это очень хороший набор данных для практики. Пожалуйста, не стесняйтесь скачать набор данных по этой ссылке:
Давайте импортируем пакеты и набор данных:
%matplotlib inline import matplotlib.pyplot as plt import seaborn as sns import pandas as pd import statsmodels.api as sm import numpy as npdf = pd.read_csv('nhanes_2015_2016.csv') Этот набор данных довольно большой, и нам не нужны все функции для этой статьи. Мы сосредоточимся на регрессионной модели, в которой переменной результата будет размер талии. Если вы хотите работать с любой другой переменной, убедитесь, что вы выбрали количественную, а не категориальную переменную. Поскольку линейная регрессия предназначена для прогнозирования только количественной переменной, а не категориальной переменной. Но вы можете использовать категориальные переменные в качестве независимых переменных. Проверьте названия всех столбцов, чтобы иметь представление о наборе данных.
df.columns#Output: Index(['SEQN', 'ALQ101', 'ALQ110', 'ALQ130', 'SMQ020', 'RIAGENDR', 'RIDAGEYR', 'RIDRETH1', 'DMDCITZN', 'DMDEDUC2', 'DMDMARTL', 'DMDHHSIZ', 'WTINT2YR', 'SDMVPSU', 'SDMVSTRA', 'INDFMPIR', 'BPXSY1', 'BPXDI1', 'BPXSY2', 'BPXDI2', 'BMXWT', 'BMXHT', 'BMXBMI', 'BMXLEG', 'BMXARML', 'BMXARMC', 'BMXWAIST', 'HIQ210'], dtype='object')
Итак, просто сохраните размер талии и несколько других переменных, которые, похоже, связаны с размером талии, и создайте новый DataFrame меньшего размера. Инстинктивно любой догадается, что вес может иметь сильную корреляцию с размером талии. Пол, ИМТ, рост тоже могут сыграть важную роль. Итак, создайте новый DataFrame меньшего размера только с этими столбцами и для удобства удалите все нулевые значения.
keep = ['BMXWT', 'BMXHT', 'BMXBMI', 'BMXWAIST', 'RIAGENDR'] db = df[keep].dropna() db.head()
Линейная регрессия и интерпретация
У нас есть набор данных, содержащий пять столбцов: вес, рост, индекс массы тела (ИМТ), размер талии и пол. Как упоминалось ранее, размер талии — это выходная переменная, которую мы попытаемся предсказать, используя другие переменные. Первоначально используйте только одну переменную или одну ковариату, чтобы предсказать размер талии. Вес (BMXWT) может быть хорошей ковариатой для начала, потому что при более высоком весе ожидается, что объем талии будет выше. Хотя есть и другие факторы, такие как рост или пол. Но о них мы подумаем позже. Подойдет модель, в которой размер талии будет выражаться как функция веса.
model = sm.OLS.from_formula("BMXWAIST ~ BMXWT", data=db) result = model.fit() result.summary() 
Эта таблица выше может показаться вам пугающей. Но большая часть информации для нас не так важна. Нам понадобится только эта часть таблицы:
В первом столбце у нас есть коэффициенты. Помните формулу линейной регрессии:
В таблице выше 42,7189 — это B, а 0,6991 — это A. И мы знаем, что A — это наклон. Итак, наш наклон равен 0,6991. Это означает, что если человек имеет одну единицу лишнего веса, его талия будет на 0,6991 единицы больше, и это основано на p-значении, упомянутом в P ›| t | столбец. Затем стандартная ошибка составляет 0,005, что указывает на расстояние этого оцененного наклона от истинного наклона. t-статистика говорит, что предполагаемый наклон 0,6991 составляет стандартную ошибку 144,292 выше нуля. Последние два столбца — это уровни достоверности. По умолчанию это уровень достоверности 95%. Доверительный интервал составляет 0,69 и 0,709, что является очень узким диапазоном. Позже мы нарисуем полосу доверительного интервала.
Стандартное отклонение составляет 16,85, что кажется намного выше, чем наклон регрессии 0,6991. Но наклон регрессии — это среднее изменение размера талии за сдвиг веса на одну единицу. Это означает, что если у человека на 10 единиц больше веса, чем у другого человека, его талия будет на 0,6991 * 10 или на 6,99 единиц больше.
Корреляция
Помимо этих значений в небольшой подтаблице, важен еще один параметр из сводки результатов. Это значение R-квадрата в верхней строке сводки результатов. Здесь значение R-квадрата составляет 0,795, что означает, что 79,5% объема талии можно объяснить весом. Теперь проверьте этот коэффициент регрессии с квадратом коэффициента Пирсона.
cc = db[["BMXWAIST", "BMXWT"]].corr() print(cc)BMXWAIST BMXWT BMXWAIST 1.000000 0.891828 BMXWT 0.891828 1.000000
Чтобы найти значение R-квадрата:
Это снова возвращает значение R в квадрате 0,795. Самая важная часть, а именно прогнозируемые значения размера талии в зависимости от веса, может быть найдена следующим образом:
Это лишь часть результата. Исходный результат намного больше.
Давайте добавим вторую ковариату и посмотрим, как она влияет на эффективность регрессии. Я выбираю пол в качестве второй ковариаты. Я хочу использовать столбец с пометкой «пол»:
db["RIAGENDRx"] = db.RIAGENDR.replace()
Вот модель и краткое описание модели:
model = sm.OLS.from_formula("BMXWAIST ~ BMXWT + RIAGENDRx", data=db) result = model.fit() result.summary() 
В приведенном выше коде BMXWT + RIAGENDRx не означает, что эти два столбца объединены или математически сложены. Это просто указывает на то, что они оба включены в модель. В этой новой модели размер талии выражается как функция веса и пола. Из приведенного выше результата мы можем найти, что коэффициент веса (BMXWT) составляет 0,7272, что немного выше, чем раньше. На этот раз коэффициент подразумевает, что если вес двух людей одного пола отличается на одну единицу, их талия будет отличаться на 0,7272 единицы. С другой стороны, коэффициент пола (RIAGENDRx) -5,0832 означает, что если мы сравним мужчину и женщину с одинаковым весом, у мужчины будет размер талии на 5,0832 единицы меньше, чем у женщины.
Все коэффициенты выражены как средние. Если мы сравним мужчину веса 70 и женщину веса 50, талия мужчины будет примерно в -5,0832 + (70–50) * 0,7272 раза больше, чем талия женщины.
Коэффициент регрессии веса (BMXWT) изменился немного, не слишком сильно после добавления пола в модель. Добавление второй ковариаты изменяет первую ковариату, если они в некоторой степени коррелированы. Давайте проверим корреляцию между двумя ковариатами:
Как видите, корреляция составляет -0,2375. Так что корреляция есть, но не слишком сильная. Вы можете найти прогнозируемые значения размера талии, используя метод подобранной стоимости, как я показал ранее. Я не буду этого снова демонстрировать.
Добавим третью ковариату. В качестве третьей ковариаты я выбрал BMXBMI. Вы можете попробовать и другие переменные.
model = sm.OLS.from_formula("BMXWAIST ~ BMXWT + RIAGENDRx + BMXBMI", data=db) result = model.fit() result.summary() 
Обратите внимание, что после добавления BMXBMI коэффициент для гендерной переменной значительно изменился. Можно сказать, что ИМТ работает как маскирующая часть связи между размером талии и гендерной переменной. С другой стороны, коэффициент веса также существенно изменился. ИМТ также работал как маскирующая часть во взаимосвязи между талией и весом.
Вы можете добавить в модель дополнительные ковариаты и увидеть эффект каждой ковариаты.
Визуализация модели
В этом разделе мы визуализируем результат регрессионной модели. Мы построим линию регрессии, которая представляет собой подогнанные значения или предсказанные значения с доверительным интервалом. Если вам нужно освежить в памяти концепции доверительного интервала, ознакомьтесь с этой статьей:
Для этого графика мы зафиксируем пол как женский, а ИМТ как 25. Кроме того, нам нужно сохранить одну независимую переменную в качестве переменной фокуса. Мы сохраним его как вес (BMXWT). Итак, график покажет прогнозируемый размер талии женщины с ИМТ 25 всех возрастов.
from statsmodels.sandbox.predict_functional import predict_functional values = pr, cb, fv = predict_functional(result, "BMXWT", values=values, ci_method="simultaneous")#Here, pr is the predicted values(pr), cb is the confidence band and #the fv is the function valuesax = sns.lineplot(fv, pr, lw=4) ax.fill_between(fv, cb[:, 0], cb[:, 1], color='grey', alpha=0.4) ax.set_xlabel("BMXWT") _ = ax.set_ylabel("BMXWAIST") 
Серая область на картинке — это полосы уверенности. Это означает, что истинный размер талии будет в этой области. Ширина серой области варьируется вдоль линии регрессии. Итак, доверительный интервал с возрастом меняется.
Вы можете зафиксировать вес и увидеть результат для определенного веса. Давайте зафиксируем вес на 65 и построим график зависимости ИМТ от талии для женского населения. Для этого нам нужно изменить параметр «значения». Потому что мы установили здесь значение ИМТ на 25. Теперь мы хотим исправить вес. Итак, нам нужно удалить значение BMI из параметра values и добавить в него вес.
del values["BMXBMI"] # Delete this as it is now the focus variable #del values['BMXWT'] values["BMXWT"] = 65 pr, cb, fv = predict_functional(result, "BMXBMI", values=values, ci_method="simultaneous")ax = sns.lineplot(fv, pr, lw=4) ax.fill_between(fv, cb[:, 0], cb[:, 1], color='grey', alpha=0.4) ax.set_xlabel("BMI") _ = ax.set_ylabel("BMXWAIST") 
На графиках выше мы построили только средние значения. Модель среднего размера талии для данного веса, пола или ИМТ. Используя ту же модель регрессии, можно также оценить структуру дисперсии, которая покажет, насколько наблюдаемые значения отклоняются от среднего. Для этого мы можем построить график остатков против подобранного значения. Помните, что подобранные значения — это предсказанные значения или наблюдаемые средние значения, а остатки — это разница между наблюдаемыми средними и истинными значениями.
pp = sns.scatterplot(result.fittedvalues, result.resid) pp.set_xlabel("Fitted values") _ = pp.set_ylabel("Residuals") 
Похоже, что когда наблюдаемые значения ниже, дисперсия немного выше.
Также возможно визуализировать эффект только одной ковариаты. Мы можем увидеть график компонент плюс остаток или график частичного остатка, используя только одну ковариату, чтобы увидеть, сохраняем ли мы другие коварианты фиксированными, как изменяются зависимые переменные:
from statsmodels.graphics.regressionplots import plot_ccprax = plt.axes() plot_ccpr(result, "BMXWT", ax) ax.lines[0].set_alpha(0.2) # Reduce overplotting with transparency _ = ax.lines[1].set_color('orange') 
Теперь мы можем увидеть влияние ИМТ таким же образом, когда весовая переменная фиксирована.
ax = plt.axes() plot_ccpr(result, "BMXBMI", ax) ax.lines[0].set_alpha(0.2) ax.lines[1].set_color("orange") 
Влияние ИМТ намного сильнее, чем вес.
В этой статье вы узнали, как подобрать модель линейной регрессии, различные статистические параметры, связанные с линейной регрессией, и некоторые хорошие методы визуализации. Методы визуализации включали построение доверительной полосы линии регрессии, построение остатков и построение графика эффекта одной ковариаты.
Дополнительные рекомендации для чтения: