Регрессионный анализ данных на языке Python

Однако мы все еще имеет широкий доверительный интервал и в целом видно, что исходные данные не могут бить идеально описаны функцией только одоной переменной LSTAT.Рисунок 2.8. Визуализация модели полиномиальной регрессии.Множественная регрессииВ предыдущем разделе мы выяснили, что для построения эффективной модели недостаточно использовать только один признак. Теперь мы построим более сложную модель, в которую войдут уже все доступные признаки.Реализация модели множественной регрессии аналогична. В отчете модели на рисунке 2.9 мы видим, что коэффициент детерминации увеличился до 0.741, что уже является хорошим результатом.Эта модель позволяет нам изучить влияние каждого признака на зависимую переменную. Например, коэффициент при CRIMравен -0.108, то есть увеличения уровня преступности на душу населения на 1 медианна цена на жилье уменьшится на 108$. Также интересно дать интерпретацию индикаторной переменной CHAS: наличие рядом реки увеличивает медианную стоимость жилья на в среднем на 2686$. При этом доверительный интервал говорит о том, что надбавка за наличие реки варьируется от 994$ до 4380$.Кроме того, в отчете мы можем видеть, что переменные AGEи INDUSявляются незначимыми (по колонке P>|t| видно, что p-valueсоответствующей гипотезы больше 0.05). Рисунок 2.9. Результат модели множественной регрессии.Лучшая модельПосле нескольких экспериментов было выявлено, что наилучший результат обеспечивает комбинация ранее описанных моделей множественной и полиномиальной регрессии (степени 2). То есть к уже имеющимся признаком добавляются все возможные полиномы степени 1 и 2. Подготовка данных и обучение модели представлены на рисунке 2.10.Рисунок 2.10. Реализация модели.Визуализировать зависимость от сразу 13 признаков практически невозможно, поэтому мы построим наглядный график зависимости от одной переменной, но вместо линии регрессии для всех значений признаков будем вычислять только предсказания модели в уже известных точках.На рисунке 2.11 представлен результат работы описанной модели.Рисунок 2.11. Результат работы лучшей модели.Как мы видим, использование дополнительных признаков позволило построить более сложную структуру линии регрессии и аппроксимировать облако точек значительно лучше. Значение коэффициента детерминации достигло 0.929, что значительно лучше всех предыдущих моделей. Такой показатель говорит о том, что модель можно использовать для предсказания медианных цен недвижимости в новых районах. Однако для этого также необходимо убедится в отсутствии переобучения модели.Проблема переобученияДля проверки отсутствия переобучения был использован метод отложенной выборки. Все данные были разделены на train/test, после чего модель была обучения только на части train, реализация представлена на рисунке 2.12.Рисунок 2.12. Обучение модели на части выборки.После этого были получены предсказания на признаках изtrainи test. Были получены коэффициенты детерминации 0.929и 0.888соответственно. Это говори о том, что если переобучение и присутствует, то оно минимально. Визуализация предсказаний представлена на рисунке 2.13.Как мы видим, качество предсказаний не тестовой выборке в большинстве точек не отличается от обучающей выборки.Рисунок 2.13. Сравнение качества предсказаний на trainи test.В попытке еще улучшить модель было решено увеличить максимальную степень полиномов до 3. Это привело к явному переобучению модели: на обучающей выборке коэффициент детерминации составил 0.999, а на тестовой 0.004. Визуализация этого переобучения представлена на рисунке 2.14. Мы видим, что на обучающей выборке предсказания идеально совпадают с исходными данными, в то время как на тестовой выборке модель ведет себя неадекватно.Рисунок 2.14. Сравнение качества предсказаний на trainи test с помощью переобученной модели.ЗаключениеВ данной работе была подробно изучена модель линейной регрессии для решения задач машинного обучения. Мы рассмотрели метод наименьших квадратов и метод максимального правдоподобия, дали вероятностную постановку задачи и описали идея получения доверительных интервалов и проверки гипотез о значимости признаков. Также мы рассмотрели проблему переобучения и изложили основные способы борьбы с ней.В практической части по разработали несколько моделей для регрессии на наборе данных Bostondataset. После изучения и визуализации данных мы построили базовую модель и улучшили ее с помощью полиномиальных признаков. После этого мы построили модель множественной регрессии и дали подробную интерпретацию результатов.Наилучший результат показала модель множественной регрессии с добавлением полиномиальных признаков, которая показала значение коэффициента детерминации 0.92. Мы проверили ее на отсутствие переобучения, а также привели пример переобученной модели, которая была получена при увеличении степени полиномов.Список литературыМуравьев Ю.Н., Панов А.П. Математическая статистика: Учебное пособие для вузов. - Москва: Юрайт, 2019.Маккинни У. Python и анализ данных: примеры и задачи. - Санкт-Петербург: Питер, 2018.Крис Элбон. Практическая статистика для специалистов Data Science: на примерах на языке Python. - Москва: ДМК Пресс, 2020.Хасте Т., Тибширани Р., Фридман Дж. Элементы статистического обучения: В 2-х ч. - Москва: Мир, 2019.James G., Witten D., Hastie T., Tibshirani R. An Introduction to Statistical Learning: with Applications in R. - New York: Springer, 2017.Murphy K.P. Machine Learning: A Probabilistic Perspective. - Cambridge, MA: MIT Press, 2012.Bishop C.M. Pattern Recognition and Machine Learning. - New York: Springer, 2006.Friedman J., Hastie T., Tibshirani R. The Elements of Statistical Learning: Data Mining, Inference, and Prediction. - New York: Springer, 2009.Wes McKinney. Python для анализа данных. - Москва: ДМК Пресс, 2017.Kelleher J.D., Tierney B. Data Science An Introduction. - Boca Raton, FL: CRC Press, 2018.Приложениеimport numpy as npimport pandas as pdimport matplotlib.pyplot as pltimport seaborn as snssns.set(rc={'figure.figsize':(10,10)})from sklearn.datasets import load_bostonboston_dataset = load_boston()boston = pd.DataFrame(boston_dataset.data, columns=boston_dataset.feature_names)boston['MEDV'] = boston_dataset.targetboston.head()# Визуализацияданныхsns.distplot(boston['MEDV'], bins=30)plt.show()correlation_matrix = boston.corr().round(2)sns.heatmap(data=correlation_matrix, annot=True)plt.figure(figsize=(10, 5))features = ['LSTAT', 'RM']target = boston['MEDV']for i, col in enumerate(features):plt.subplot(1, len(features) , i+1) x = boston[col] y = targetplt.scatter(x, y, marker='o')plt.title(col)plt.xlabel(col)plt.ylabel('MEDV')# Модель простой регрессии от LSTATboston = boston.sort_values('LSTAT')X = boston['LSTAT']Y = boston['MEDV']import statsmodels.api as smX = sm.add_constant(X)model = sm.OLS(Y, X)results = model.fit()print(results.summary())print("Parameters: ", results.params)print("R2: ", results.rsquared)pred_ols = results.get_prediction()iv_l = pred_ols.summary_frame()["obs_ci_lower"]iv_u = pred_ols.summary_frame()["obs_ci_upper"]fig, ax = plt.subplots(figsize=(8, 6))ax.plot(boston["LSTAT"], Y, "o", label="data")ax.plot(boston["LSTAT"], results.fittedvalues, "g-", linewidth=4, label="OLS")ax.plot(boston["LSTAT"], iv_u, "r--", label="Confidence interval")ax.plot(boston["LSTAT"], iv_l, "r--")ax.legend(loc="best")ax.set_xlabel("LSTAT")ax.set_ylabel("MEDV")# Модельполиномиальнойрегрессииот LSTATfrom sklearn.preprocessing import PolynomialFeaturesX = PolynomialFeatures(3).fit_transform(np.array(boston['LSTAT']).reshape(-1,1))Xmodel = sm.OLS(Y, X)results = model.fit()print(results.summary())print("Parameters: ", results.params)print("R2: ", results.rsquared)pred_ols = results.get_prediction()iv_l = pred_ols.summary_frame()["obs_ci_lower"]iv_u = pred_ols.summary_frame()["obs_ci_upper"]fig, ax = plt.subplots(figsize=(8, 6))ax.plot(boston['LSTAT'], Y, "o", label="data")ax.plot(boston['LSTAT'], results.fittedvalues, "g-", linewidth=4, label="OLS")ax.plot(boston['LSTAT'], iv_u, "r--", label="Confidence interval")ax.plot(boston['LSTAT'], iv_l, "r--")ax.legend(loc="best")ax.set_xlabel("LSTAT")ax.set_ylabel("MEDV")# Множестваеннаярегрессияотвсехпеременныхboston = boston.sort_values('LSTAT')X = boston.drop('MEDV', axis=1)Y = boston['MEDV']X = sm.add_constant(X)model = sm.OLS(Y, X)results = model.fit()print(results.summary())# Лучшаямодельboston = boston.sort_values('LSTAT')X = boston.drop('MEDV', axis=1)Y = boston['MEDV']X = PolynomialFeatures(2).fit_transform(X)model = sm.OLS(Y, X)results = model.fit()results.rsquaredpred_ols = results.predict(X)fig, ax = plt.subplots(figsize=(8, 6))ax.plot(boston['LSTAT'], Y, "o", label="data")ax.plot(boston['LSTAT'], pred_ols, "gx", linewidth=4, label="OLS")ax.legend(loc="best")ax.set_xlabel("LSTAT")ax.set_ylabel("MEDV")from sklearn.model_selection import train_test_splitX_train, X_test, Y_train, Y_test = train_test_split(X, Y, test_size=0.2, random_state=42)model = sm.OLS(Y_train, X_train)results = model.fit()from sklearn.metrics import r2_scoreprint("R2: ", r2_score(results.predict(X_train), Y_train))print("R2: ", r2_score(results.predict(X_test), Y_test))pred_train = results.predict(X_train)pred_test = results.predict(X_test)plt.figure(figsize=(12, 6))plt.subplot(1,2,1)plt.plot(X_train[:,13], Y_train, "o", label="data")plt.plot(X_train[:,13], pred_train, "gx", linewidth=4, label="OLS")plt.legend(loc="best")plt.xlabel("LSTAT")plt.ylabel("MEDV")plt.title("Train")plt.subplot(1,2,2)plt.plot(X_test[:,13], Y_test, "o", label="data")plt.plot(X_test[:,13], pred_test, "gx", linewidth=4, label="OLS")plt.legend(loc="best")plt.xlabel("LSTAT")plt.ylabel("MEDV")plt.title("Test")# Примерпереобученияboston = boston.sort_values('LSTAT')X = boston.drop('MEDV', axis=1)Y = boston['MEDV']X = PolynomialFeatures(3).fit_transform(X)X_train, X_test, Y_train, Y_test = train_test_split(X, Y, test_size=0.2, random_state=42)model = sm.OLS(Y_train, X_train)results = model.fit()print("R2: ", results.rsquared)print("R2: ", r2_score(results.predict(X_test), Y_test))pred_train = results.predict(X_train)pred_test = results.predict(X_test)plt.figure(figsize=(12, 6))plt.subplot(1,2,1)plt.plot(X_train[:,13], Y_train, "o", label="data")plt.plot(X_train[:,13], pred_train, "gx", linewidth=4, label="OLS")plt.legend(loc="best")plt.xlabel("LSTAT")plt.ylabel("MEDV")plt.title("Train")plt.subplot(1,2,2)plt.plot(X_test[:,13], Y_test, "o", label="data")plt.plot(X_test[:,13], pred_test, "gx", linewidth=4, label="OLS")plt.legend(loc="best")plt.xlabel("LSTAT")plt.ylabel("MEDV")plt.title("Test")