Линейная регрессия (LR) - это основной метод статистического моделирования, при котором числовые результаты прогнозируются с помощью линейной комбинации предикторов с коэффициентами βᵢ. Модель линейной регрессии была чрезвычайно популярной и популярной среди практиков машинного обучения из-за ее простоты и хорошо работает в большинстве задач, если выполняется предположение о линейности. Регуляризация ограничивает или сужает оценки коэффициентов до нуля. Оценка коэффициента сокращения может уменьшить дисперсию.
В конце этого поста вы узнаете, как:
1) Разработайте модель линейной регрессии и модели регуляризованной регрессии (например, гребень, лассо и эластичная сеть) в R.
2) Выполните повторную тысячную перекрестную проверку (CV), чтобы оценить производительность регрессионных моделей.
Сбор и визуализация данных
Сначала загрузите необходимые пакеты и импортируйте данные (данные о белом вине находятся в свободном доступе в репозитории машинного обучения UCI и Kaggle) в глобальную среду Rstudio. Затем используйте str() и summary() для отображения типов атрибутов и статистической сводки каждой переменной. Отмечено, что качество вина (переменная ответа) находится в целочисленной форме, в которой высокая оценка указывает на высокое качество и наоборот.
library(ggplot2) library(rsample) # stratified sampling library(reshape2) # data wrangling library(knitr) library(car) # variance inflation factor library(glmnet) # regularized linear regression library(glmnetUtils) library(gridExtra) library(dplyr) # data wrangling data=read.csv(‘winequality-white.csv’,sep = “;”) str(data) # statistical summary of each attribute summary(data)
Выполните стратифицированное разбиение, чтобы разделить данные на обучающий и тестовый наборы данных.
set.seed(2) split=initial_split(data,prop=0.8,strata = “quality”) train=training(split) test=testing(split) ggplot(data=train,aes(x=quality)) + geom_bar(fill=”steelblue”) + theme_classic() + ggtitle(‘distribution of response variable in training dataset)
Большинство белых вин имеют рейтинг качества от 5 до 7, при этом небольшое количество образцов находится на низком и высоком уровне.
Всегда полезно вычислить матрицу корреляции предикторов и визуализировать ее на тепловой карте, чтобы узнать коллинеарность (степень корреляции) между предикторами.
vars=setdiff(colnames(train),’quality’) corr=cor(train[,vars]) library(reshape2) corr[upper.tri(corr)]=NA melted_cormat=melt(corr,na.rm = TRUE) # heatmap ggplot(data=melted_cormat,aes(Var1,Var2,fill=value)) + geom_tile(color=’white’) + scale_fill_gradient2(low = “blue”, high=”red”,mid=”white”,midpoint = 0, limit=c(-1,1),space = “Lab”,name=”Correlation”) + theme_minimal()+ theme(axis.text.x = element_text(angle = 45,vjust = 1,size = 12,hjust = 1)) + coord_fixed() + geom_text(aes(Var1,Var2,label=round(value,2)),color=”black”,size=3) + xlab(“”) + ylab(“”)
Существует сильная положительная линейная корреляция между переменными «плотность» и «остаточный сахар». Переменные «плотность» и «алкоголь» имеют отрицательную корреляцию.
Линейная регрессия
Наличие коллинеарности может создать проблемы в контексте регрессии, поскольку может быть трудно выделить индивидуальное влияние коллинеарных переменных на отклик.
Гарет Джеймс и др. др. Введение в статистическое обучение: с приложениями в Р. Нью-Йорк: Springer, 2013 г. (стр. 99)
Коллинеарность снижает точность оценок коэффициентов регрессии. Помимо матрицы корреляции, более надежным способом является расчет коэффициента инфляции дисперсии (VIF). Как правило, VIF более 10 указывает на высокую корреляцию.
Затем я подбираю модель линейной регрессии, используя все предикторы, и вычисляю VIF для каждого предиктора.
# Train linear regression model model_lm=lm(quality~.,data = train) summary(model_lm) kable(data.frame(vif=vif(model_lm)))
VIF переменных «плотность» и «остаточный сахар» больше 10. Таким образом, одним из способов смягчения коллинеарности является удаление одного из предикторов, как показано ниже.
drops=’density’ model_lm=lm(quality~.,data=train[,!(colnames(train) %in% drops)]) summary(model_lm) kable(data.frame(vif=vif(model_lm)))
Удаление переменной «плотность» приводит к изменению оценок коэффициентов, коэффициента детерминации (статистика R-квадрата) и остаточной стандартной ошибки модели линейной регрессии (LR). Из сводки модели мы знаем, что «летучие кислоты» отрицательно коррелируют с качеством белого вина, тогда как «сульфаты» и «алкоголь» положительно коррелируют с качеством белого вина.
Регулярная линейная регрессия
Регуляризованная линейная регрессия может быть реализована в R с помощью пакета glmnet. glmnet решает следующую проблему:
по сетке значений λ, охватывающей весь диапазон возможных решений.
Регрессия хребта
Модель гребневой регрессии можно обучить, установив входной аргумент в функции 'cv.glmnet', альфа равной 0. Стандартные оценки коэффициентов линейной регрессии не зависят от масштаба, но оценки коэффициентов гребневой регрессии могут существенно измениться, когда умножение заданного предиктора на константу. Таким образом, перед обучением модели гребневой регрессии мы должны стандартизировать предикторы.
# First standardize the predictors X_stand=scale(train[,vars],center = FALSE, scale = TRUE) train_stand=data.frame(X_stand,quality=train$quality) X1_stand=scale(test[,vars],center = FALSE, scale = attributes(X_stand)$’scaled:scale’) test_stand=data.frame(X1_stand,quality=test$quality) model_ridge=cv.glmnet(quality~.,data=train_stand,alpha = 0) # visualize coefficients coef_ridge=coef(model_ridge) df_coef_ridge=data.frame(coef_name=rownames(coef_ridge)[-1], coef=coef_ridge[-1,1]) ggplot(data = df_coef_ridge,aes(x=coef_name,y=coef)) + geom_bar(stat = ‘identity’, width=0.5) + coord_flip() # visualized cross-validation error with repect to different lambda roundoff=function(x) sprintf(“%.3f”,x) df_cvm_ridge=data.frame(lambda=model_ridge$lambda,cvm=model_ridge$cvm, cvm_up=model_ridge$cvup,cvm_down=model_ridge$cvlo) ggplot(data=df_cvm_ridge,aes(x=lambda,y=cvm,ymin=cvm_down,ymax=cvm_up)) + geom_line() + geom_point() + geom_errorbar() + scale_x_continuous(trans = “log2”,labels = roundoff) + geom_vline(xintercept = model_ridge$lambda.min, linetype=2, color=”red”) + geom_vline(xintercept = model_ridge$lambda.1se, linetype=3, color=”blue”) + ggtitle(‘ridge regression’)
Функция cv.glmnet() выполняет перекрестную проверку, чтобы выбрать лямбду, дающую минимальную ошибку перекрестной проверки. При вызове такой функции, как предсказание(), объект cv.glmnet по умолчанию использует модель lambda.1se, которая представляет собой наибольшее значение лямбда, при котором возникает ошибка. в пределах 1 стандартной ошибки от минимума, lamda.min.
Лассо
Модель регрессии Лассо можно обучить, установив входной аргумент в функции ‘cv.glmnet’, альфа равной 1, как показано ниже.
model_lasso=cv.glmnet(quality~.,data=train_stand,alpha = 1) # visualize coefficients coef_lasso=coef(model_lasso) df_coef_lasso=data.frame(coef_name=rownames(coef_lasso)[-1], coef=coef_lasso[-1,1]) ggplot(data = df_coef_lasso,aes(x=coef_name,y=coef)) + geom_bar(stat = ‘identity’, width=0.5) + coord_flip() # visualized cross-validation error with repect to different lambda # plus the number of nonzero coefficients df_cvm_lasso=data.frame(lambda=model_lasso$lambda,cvm=model_lasso$cvm, cvm_up=model_lasso$cvup,cvm_down=model_lasso$cvlo) p1=ggplot(data=df_cvm_lasso,aes(x=lambda,y=cvm,ymin=cvm_down,ymax=cvm_up)) + geom_line() + geom_point() + geom_errorbar() + scale_x_continuous(trans = “log2”,labels = roundoff) + geom_vline(xintercept = model_lasso$lambda.min, linetype=2, color=”red”) + geom_vline(xintercept = model_lasso$lambda.1se, linetype=3, color=”blue”) + ggtitle(‘Lasso’) df_nzero_lasso=data.frame(lambda=model_lasso$lambda,nzero=model_lasso$nzero) p2=ggplot(data=df_nzero_lasso,aes(x=lambda,y=nzero)) + geom_line(color=”steelblue”) + scale_x_continuous(trans = “log2”,labels = roundoff) + ylab(‘Number of non-zero coefficients’) grid.arrange(p1,p2,ncol=2)
Эластичная сетка
Эластичная сеть — это компромисс между гребнем и лассо, в котором альфа находится между 0 и 1. Я использую функцию cva.glmnet() в glmnetUtils, чтобы определить наилучшую альфа-версию, затем используйте cv.glmnet(), установив альфу как лучшую альфа-версию.
elastic_net=cva.glmnet(quality~.,train_stand)
# function to get cvm of the model$lambda.1se in each alpha
get_cvm=function(model) {
index <- match(model$lambda.1se, model$lambda)
model$cvm[index]
}
# data frame that contains alpha and its corresponding cross-validation error
enet_performance=data.frame(alpha=elastic_net$alpha)
models=elastic_net$modlist
enet_performance$cvm=vapply(models,get_cvm,numeric(1))
# get the best alpha and train the model
best_alpha=enet_performance[which.min(enet_performance$cvm),’alpha’]
model_enet=cv.glmnet(quality~.,train_stand,alpha = best_alpha)
Как показано на гистограмме выше, оценки коэффициентов для разных методов регуляризации различались. Важно отметить, что регуляризация влияет на интерпретируемость модели, поскольку у нас нет стандартной ошибки и значения p для оценок коэффициентов, позволяющих оценить их статистическую значимость. Лассо и эластичная сетевая регрессия могут рассматриваться как своего рода переменная селекция. Чем больше лямбда, тем больше оценок коэффициентов будет стремиться к нулю.
Коды и соответствующие результаты для оценки модели с использованием набора тестов можно найти в RPubs.
k-кратная перекрестная проверка
По сути, k-кратная перекрестная проверка (CV) — это метод повторной выборки, который включает в себя разделение набора данных на k отдельных групп. Всего для обучения модели используется (k-1) групп экземпляров данных, а оставшийся один блок данных используется для измерения производительности модели. Этот шаг повторяется k раз. Таким образом, мы получим k показателей производительности, что позволит нам рассчитать стандартное отклонение и построить доверительный интервал. Повторение k-кратного CV просто означает, что процедура k-кратного CV повторяется l раз, где l — произвольный параметр. Коды можно найти на Github и RPubs.
Резюме
Среднеквадратическая ошибка (RMSE), средняя абсолютная ошибка (MAE) и статистика R-квадрата — вот некоторые из наиболее часто используемых показателей производительности для регрессии. Чем ниже RMSE и MAE и чем ближе R-квадрат приближается к 1, тем лучше производительность. Однако R2 может вводить в заблуждение, если дисперсия ошибки велика. Показатели каждой модели, оцененные при стратифицированной выборке и повторном k-кратном CV, представлены в таблице, как показано ниже:
В этом случае линейная регрессия работает немного лучше, чем модели регуляризованной регрессии. Низкий R² означает, что более 70% вариаций данных остаются необъяснимыми линейными моделями. Нелинейные методы, такие как нейронные сети и метод опорных векторов, являются хорошими кандидатами для повышения точности прогнозирования.
Ссылка
П. Кортес, А. Сердейра, Ф. Алмейда, Т. Матос и Дж. Рейс. Моделирование винных предпочтений путем извлечения данных из физико-химических свойств. В системах поддержки принятия решений, Elsevier, 47 (4): 547–553. ISSN: 0167–9236.
