Построение модели линейной регрессии
Сорокина Валерия
2 марта 2019 г
Введение
В файле checkout.xls содержатся данные о времени обслуживания на кассе 18 покупателей:
- customer - номер покупателя;
- time - время обслуживания;
- amount - сумма покупки; - items - количество товаров в чеке.
library(readxl) # считывание данных из файлов Excel
library(tidyverse) # манипулирование данными
library(GGally) # построение матрицы диаграмм рассеяния
library(ggfortify) # визуализация диагностических графиков
library(modelr) # вспомогательные функции для работы с моделями
library(broom) # преобразование результатов моделирования в табличный вид1. Загрузка данных.
checkout <- read_excel('checkout.xls')2. Корреляционная матрица
Постройте корреляционную матрицу времени обслуживания, суммы покупки и количества товаров.
Коэффициент корреляции Пирсона r позволяет количественно оценить степень линейной связи между двумя непрерывными переменными.
options(digits = 3) # число знаков после запятой в выводе
checkout %>%
select(-customer) %>%
cor()## time amount items
## time 1.000 0.959 0.876
## amount 0.959 1.000 0.923
## items 0.876 0.923 1.000Наблюдается сильная зависимость переменных между собой.
3. Линейная регрессия
Постройте две модели простой линейной регрессии для зависимостей времени обслуживания от суммы покупки и от количества товаров. Запишите для каждой модели формулу для прогнозирования времени обслуживания (на основе вывода lm). Проверьте значимость модели с помощью F-критерия. Проверьте значимость коэффициентов модели с помощью t-критериев.
Модели простой линейной регрессии для зависимости времени обслуживания от суммы покупки
m_amount <- lm(time ~ amount, data=checkout)
# Выделение коэффициентов
coef_amount <- coef(m_amount)
coef_amount## (Intercept) amount
## 0.2633 0.0821\(time = 0.2633 + 0.0821*amount\)
Таким образом, при увеличении стоимости покупки на 100, время обслуживания, в среднем, увеличивается на 8,21.
checkout %>%
add_predictions(m_amount) %>%
ggplot() +
geom_point(aes(amount, time)) +
geom_line(aes(amount, pred), color = 'red') +
labs(title='Регрессионная прямая для оценки времени обслуживания\nв зависимости от суммы покупки')
Проверим значимость модели:
anova(m_amount)## Analysis of Variance Table
##
## Response: time
## Df Sum Sq Mean Sq F value Pr(>F)
## amount 1 128.9 128.9 186 3.2e-10 ***
## Residuals 16 11.1 0.7
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Модель на основе переменной amount значима, т.к. p-value < .05.
Проверим значимость коэффициентов модели:
confint(m_amount)## 2.5 % 97.5 %
## (Intercept) -0.4761 1.0027
## amount 0.0693 0.0948Первый интервал включает 0, поэтому гипотезу о равенстве нулю соответствующего коэффициента для совокупности нельзя отвергнуть на 5% уровне значимости.
Второй коэффициент оказался значимым, поскольку соответствующие p-значения меньше порога 5%, а доверительный интервал для него не включает 0.
Модели простой линейной регрессии для зависимости времени обслуживания от количества товаров
m_items <- lm(time ~ items, data=checkout)
# Выделение коэффициентов
coef_items <- coef(m_items)
coef_items## (Intercept) items
## -1.277 0.526\(time = -1.277 + 0.526*amount\)
Таким образом, при увеличении количества товаров на 10, время обслуживания, в среднем, увеличивается на 5,26.
checkout %>%
add_predictions(m_items) %>%
ggplot() +
geom_point(aes(items, time)) +
geom_line(aes(items, pred), color = 'red') +
labs(title='Регрессионная прямая для оценки времени обслуживания\nв зависимости от количества товаров')
Проверим значимость модели:
anova(m_items)## Analysis of Variance Table
##
## Response: time
## Df Sum Sq Mean Sq F value Pr(>F)
## items 1 107.5 108 52.8 1.9e-06 ***
## Residuals 16 32.5 2
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Модель на основе переменной items значима, т.к. p-value < .05.
Проверим значимость коэффициентов модели:
confint(m_items)## 2.5 % 97.5 %
## (Intercept) -3.023 0.469
## items 0.373 0.680Первый интервал включает 0, поэтому гипотезу о равенстве нулю соответствующего коэффициента для совокупности нельзя отвергнуть на 5% уровне значимости.
Второй коэффициент оказался значимым, поскольку соответствующие p-значения меньше порога 5%, а доверительный интервал для него не включает 0.
4. Доверительные интервалы
Постройте доверительные интервалы для среднего значения и отдельных наблюдений. Какая из двух моделей наиболее точная?
Доверительные интервалы для среднего значения
checkout %>%
add_predictions(m_amount) %>%
ggplot(aes(amount, time)) +
geom_point() +
geom_smooth(method = 'lm',
color = 'red') +
geom_smooth(method = 'lm',
level = .8,
color = 'red') +
labs(title = '80% и 90% доверительные интервалы для регрессионной прямой\nзависимости оценки времени обслуживания от суммы покупки')
checkout %>%
add_predictions(m_items) %>%
ggplot(aes(amount, items)) +
geom_point() +
geom_smooth(method = 'lm',
color = 'red') +
geom_smooth(method = 'lm',
level = .8,
color = 'red') +
labs(title = '80% и 90% доверительные интервалы для регрессионной прямой\nзависимости оценки времени обслуживания от количества товаров')
Вывод: в данном случае наблюдается большая зависимость времени обслуживания от суммы покупки нежели от количества товара.
Доверительные интервалы для среднего значения и отдельных наблюдений
# 95% доверительный интервал для среднего
ci_amount <- predict(m_amount,
newdata = checkout,
interval = 'confidence',
level = 0.95) %>%
as.data.frame() %>%
rename(lci = lwr, uci = upr)
# 95% доверительный интервал для отдельных наблюдений
pi_amount <- predict(m_amount,
newdata = checkout,
interval = 'prediction',
level = 0.95) %>%
as.data.frame() %>%
select(-fit) %>% # точечный прогноз уже есть в наборе, удаляем этот столбец
rename(lpi = lwr, upi = upr)
# Склеиваем все в одну таблицу
pred_amount <- bind_cols(checkout, ci_amount, pi_amount)
ggplot(data = pred_amount) +
geom_ribbon(aes(amount, ymin = lpi, ymax = upi),
fill = 'lightskyblue', alpha = 0.5) +
geom_ribbon(aes(amount, ymin = lci, ymax = uci),
fill = 'darkgray', alpha = 0.5) +
geom_line(aes(amount, fit), color = 'red') +
geom_point(aes(amount, time)) +
labs(title = paste('Интервальный прогноз (95% интервалы) оценки времени обслуживания',
'в зависимости от суммы покупки',
sep = '\n'),
y = 'Оценка времени обслуживания' )
# 95% доверительный интервал для среднего
ci_items <- predict(m_items,
newdata = checkout,
interval = 'confidence',
level = 0.95) %>%
as.data.frame() %>%
rename(lci = lwr, uci = upr)
# 95% доверительный интервал для отдельных наблюдений
pi_items <- predict(m_items,
newdata = checkout,
interval = 'prediction',
level = 0.95) %>%
as.data.frame() %>%
select(-fit) %>% # точечный прогноз уже есть в наборе, удаляем этот столбец
rename(lpi = lwr, upi = upr)
# Склеиваем все в одну таблицу
pred_items <- bind_cols(checkout, ci_items, pi_items)
ggplot(data = pred_items) +
geom_ribbon(aes(items, ymin = lpi, ymax = upi),
fill = 'lightskyblue', alpha = 0.5) +
geom_ribbon(aes(items, ymin = lci, ymax = uci),
fill = 'darkgray', alpha = 0.5) +
geom_line(aes(items, fit), color = 'red') +
geom_point(aes(items, time)) +
labs(title = paste('Интервальный прогноз (95% интервалы) оценки времени обслуживания',
'в зависимости от количества товаров',
sep = '\n'),
y = 'Оценка времени обслуживания' )
Интерпретация построенных графиков следующая:
Можно с 95% уверенностью утверждать, что регрессионная прямая совокупности находится внутри доверительного интервала для среднего (серая область на графике)
С 95% вероятностью фактические значения будут находиться внутри доверительного интервала для отдельных наблюдений (синяя область на графике)
Ширина доверительного интервала для отдельных наблюдений существенно больше, чем для среднего.
Источники: