Введение

Для анализа были взяты предварительно обработанные данные из трех источников: СберИндекс, ДомКлик и Росстат, представленные в открытом доступе и отображающие распределение социоэкономических показателей по регионам России. Итоговый объединенный датафрейм, взятый для анализа, включает 15 показателей (Итоговая таблица.xlsx):

Регион
Месяц
Квартал
Индекс потребительской активности
Доля безналичных платежей
Количество внутренних туристов
Общее количество одобренных заявок на кредит
Доля онлайн-заявок на кредит
Доля заявок на кредит в офисе банка
Общее количество ипотечных сделок
Доля сделок, первичное жилье
Доля сделок, вторичное жилье
Среднемесячная заработная плата
Число активных абонентов беспроводного интернета
Уровень безработицы населения

Задачи исследования:

Проанализировать особенности распределения указанных выше показателей
Выявить регионы России, лидирующие по ряду социоэкономических показателей
Выявить ключевые факторы, определяющие спрос населения на банковские услуги

На основании сравнительного анализа социально-экономических показателей выбрать наиболее перспективные регионы для открытия филиалов банка TrustUs

Разведывательный анализ:

# Подгружаем библиотеки:
library(tidyverse)

## Warning in Sys.timezone(): unable to identify current timezone '!0@0B>2 ':
## please set environment variable 'TZ'

library(readxl)
library(psych)
library(lmtest)
library(sandwich)
library(pROC)
library(stargazer)

# Загружаем данные:
dat <- read_xlsx("Итоговая таблица.xlsx")
dat$`Число активных абонентов беспроводного доступа к сети` <- as.numeric(dat$`Число активных абонентов беспроводного доступа к сети`)
dat$Месяц <- factor(dat$Месяц, ordered = TRUE,
                        levels = c("Январь", "Февраль", "Март", 
                                                "Апрель", "Май", "Июнь", "Июль",
                                                "Август", "Сентябрь", "Октябрь",
                                                "Ноябрь", "Декабрь"))
dat$Квартал <- factor(dat$Квартал,
                         levels = c("Квартал 1", "Квартал 2", "Квартал 3", "Квартал 4"))

dat$`Всего одобренных заявок` <- factor(dat$`Всего одобренных заявок`,
                                               ordered = TRUE,
                                               levels = c("100 - 500",
                                                          "500 - 1 000",
                                                          "1 000 - 5 000",
                                                          "5 000 - 10 000",
                                                          "> 10 000"))
dat$`Всего ипотечных сделок` <- factor(dat$`Всего ипотечных сделок`,
                                          ordered = TRUE,
                                          levels = c("10 - 50",
                                                     "50 - 100",
                                                     "100 - 500",
                                                     "500 - 1 000",
                                                     "> 1 000"))

# Выведем описательные статистики для показателей СберИндекс:
dat %>% 
  select(`Индекс потребительской активности` : `Количество внутренних туристов`) %>% 
  describe()

##                                   vars   n  mean    sd median trimmed   mad
## Индекс потребительской активности    1 912 65.51  8.80  67.00   66.20  7.41
## Доля безналичных платежей            2 912 52.33  7.67  53.25   53.22  5.86
## Количество внутренних туристов       3 912 -7.46 28.25  -4.54   -6.96 28.50
##                                      min   max  range  skew kurtosis   se
## Индекс потребительской активности  34.00 88.00  54.00 -0.73     0.61 0.29
## Доля безналичных платежей          16.00 66.10  50.10 -1.85     5.19 0.25
## Количество внутренних туристов    -85.47 86.14 171.61 -0.09    -0.18 0.94

Можно отметить, что три исследуемых показателя (Индекс потребительской активности, доля безналичных платежей, количество внутренних туристов) имеют высокий размах, т.е. значительную разницу между минимальными и максимальными значениями показателей. Вместе с этим, имеются небольшие различия между медианными и средними значениями, что может говорить о наличии нетипичных значений. Значения skew во всех трех случаях < 0, что говорит о скошенности распределения влево.

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

## Построим график, отображающий динамику средних значений данных показателей по месяцам:
dat_mean <- dat %>%
  drop_na %>% 
  group_by(Месяц) %>% 
  summarise(`Индекс потребительской активности` = mean(`Индекс потребительской активности`),
            `Доля безналичных платежей` = mean(`Доля безналичных платежей`),
            `Количество внутренних туристов` = mean(`Количество внутренних туристов`))
dat_mean <- pivot_longer(data = dat_mean, cols = 2:4, names_to = "Индекс", values_to = "Значение")

dat_mean %>% 
  ggplot(aes(x = Месяц, y = Значение, group = Индекс, color = Индекс)) +
  geom_line(size = 0.8) +
  geom_point(size = 1.2) +
  labs(title = "Динамика средних показателей СберИндекс в 2020 году",
       x = "Месяц",
       y = "Среднее значение по России") +
  scale_colour_manual(values = c("tomato", "steelblue", "palegreen3")) +
  theme(axis.text.x = element_text(angle = 90, vjust = 0.5, hjust=1),
        panel.background = element_rect(fill = "white", colour = "grey50"),
        panel.grid.major = element_line(color = "grey95"),
        panel.grid.minor = element_line(color = "grey98"))

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

dat_long <- dat %>%
  drop_na %>% 
  select(1:6) %>% 
  pivot_longer(cols = 4:6,names_to = "Индекс", values_to = "Значение")

dat_long %>% 
  ggplot(aes(y = Значение, group = Индекс, fill = Индекс)) +
  geom_boxplot(color = "gray30",
               outlier.fill = "coral",
               outlier.shape = 21,
               outlier.size = 1.8) +
  scale_fill_manual(values = c("lavender", "lightpink", "steelblue")) +
  facet_wrap(facets = ~Месяц) +
  scale_y_continuous(breaks = c(-50, -25, 0, 25, 50, 75)) +
  coord_cartesian(ylim = c(-80, 85)) +
  labs(title = "Разброс значений показателей СберИндекс по месяцам") +
  theme(axis.ticks.x = element_blank(), 
        axis.text.x = element_blank(),
        panel.background = element_rect(fill = "grey99"),
        panel.grid.major = element_line(color = "grey95"),
        panel.grid.minor = element_line(color = "grey97"),
        strip.background = element_rect(color = "grey70", fill = "grey95"))

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

# Поиск нетипичных значений:
dat_clean <- dat %>% 
  drop_na()

get_outliers <- function(x){
  iqr <- IQR(x)
  q1 <- quantile(x, 0.25)
  q3 <- quantile(x, 0.75)
  lower_b <- q1 - 1.5 * iqr
  upper_b <- q3 + 1.5 * iqr
  extra <- x[x < lower_b | x > upper_b]
  return(extra)
}

consume_out <- dat[dat$`Индекс потребительской активности` %in% get_outliers(dat_clean$`Индекс потребительской активности`), 1:4]
cashless_out <- dat[dat$`Доля безналичных платежей` %in% get_outliers(dat_clean$`Доля безналичных платежей`), c(1:3, 5)]
voyage_out <- dat[dat$`Количество внутренних туристов` %in% get_outliers(dat_clean$`Количество внутренних туристов`), c(1:3, 6)]

dat_out <- full_join(consume_out, cashless_out)
dat_out <- full_join(dat_out, voyage_out)

dat_out <- dat_out %>% 
  pivot_longer(cols = `Индекс потребительской активности` : `Количество внутренних туристов`, names_to = "Индекс", values_to = "Значение")
dat_out <- dat_out %>% 
  drop_na()

Рассмотрим, каким месяцам соответствует наибольшее количество нетипичных значений:

# Визуализируем нетипичные значения:
dat_out %>% 
  ggplot(aes(x = Месяц, fill = Индекс)) +
  geom_bar(color = "grey20") +
  labs(title = "Распределение выбросов по месяцам",
       x = "Месяц",
       y = "Количество") +
  scale_fill_manual(values = c("wheat", "lightpink2", "slategray3")) +
  theme(axis.text.x = element_text(angle = 90),
        panel.background = element_rect(fill = "grey99"),
        panel.grid.major = element_line(color = "grey95"),
        panel.grid.minor = element_line(color = "grey97"))

Нетипичные значения в показателе “Доля безналичных платежей” равномерно распределены и встречаются в каждом месяце. Соответственно, это сможет быть связано не с периодом, а с определенными регионами. Нетипичные значения количества внутренних туристов встречаются в феврале, марте, августе и октябре. Вероятно, внутренний туризм в период с апреля по август был затруднен вследствие ограничений. Выбросы в показателе “Индекс потребительской активности” встречаются в четырех месяцах и с максимумом в апреле, что тоже укладывается в предположение о влиянии ограничений, связанных с пандемией COVID-19.

# По регионам:
dat_out %>% 
  ggplot(aes(x = Регион, fill = Индекс)) +
  geom_bar(color = "grey20") +
  labs(title = "Распределение выбросов по регионам",
       x = "Регион",
       y = "Количество") +
  scale_fill_manual(values = c("wheat", "lightpink2", "slategray3")) +
  theme(axis.text.x = element_text(angle = 90),
        panel.background = element_rect(fill = "grey99"),
        panel.grid.major = element_line(color = "grey95"),
        panel.grid.minor = element_line(color = "grey97"))

Если посмотреть на распределение нетипичных значений по регионам, можно заметить, что выбросы среди значений индекса потребительской активности распределены равномерно, за исключением трех регионов. Так как нетипичные значения в апреле по этому показателю ниже граничного значения (q1 - 1.5 * IQR), то, вероятно, потребительская активность в этих регионах пострадала наиболее сильно. Нетипичные значения доли безналичных платежей соответствуют четырем регионам: Кабардино-Балкарской Республике, Республике Дагестан, республике Карачаево-Черкессии, Республике Северной Осетии-Алании. Так как нетипично низкие значения этого показателя распределены равномерно по месяцам, можно предположить, что в этих регионах традиционно более распространены платежи наличными.

# Находим топ-30 регионов по индексу потребительской активности:
top_consume <- dat %>%
  drop_na %>% 
  group_by(Регион) %>% 
  summarise(mean_consume = mean(`Индекс потребительской активности`)) %>%
  arrange(desc(mean_consume)) %>%
  slice_head(n = 30)

# Находим топ-30 регионов по доле безналичных платежей:
top_cashless <- dat %>%
  drop_na %>% 
  group_by(Регион) %>% 
  summarise(mean_cashless = mean(`Доля безналичных платежей`)) %>%
  arrange(desc(mean_cashless)) %>%
  slice_head(n = 30)

# Находим регионы, общие для обоих датафреймов:
top_both <- inner_join(top_consume, top_cashless, by = "Регион")
top_both

## # A tibble: 15 × 3
##    Регион                  mean_consume mean_cashless
##    <chr>                          <dbl>         <dbl>
##  1 Калининградская область         74.8          58.9
##  2 Воронежская область             73.4          54.5
##  3 Тюменская область               73.2          60.0
##  4 Хабаровский край                71            57.7
##  5 Кировская область               70.3          57.7
##  6 Иркутская область               68.5          56.6
##  7 Свердловская область            68.1          57.3
##  8 Самарская область               67.7          54.9
##  9 Москва                          67.7          57.1
## 10 Нижегородская область           67.6          55.7
## 11 Республика Башкортостан         67.4          55.1
## 12 Камчатский край                 67.3          60.4
## 13 Кемеровская область             66.5          55.3
## 14 Томская область                 66.5          59.0
## 15 Удмуртская Республика           66.4          57.4

# Выберем в основном датафрейме только регионы, попавшие в топ:
top_both <- dat[dat$Регион %in% top_both$Регион, ]

# Среди полученных регионов оставим только те, в которых индекс внутреннего туризма принимал значения выше 0 не менее чем в течение 5 месяцев:
top_both <- top_both %>% 
  mutate(tourindex = ifelse(`Количество внутренних туристов` > 0, 1, 0), .after = `Количество внутренних туристов`)

top_tour <- top_both %>% 
  group_by(Регион) %>% 
  summarise(tour_sum = sum(tourindex)) %>% 
  filter(tour_sum >= 5)

Выберем регионы, лидирующие по всем трем показателям:

# Выберем из основного датафрейма регионы, лидирующие по всем трем показателями:
top <- dat[dat$Регион %in% top_tour$Регион, ]
unique(top$Регион)

## [1] "Воронежская область"     "Калининградская область"
## [3] "Камчатский край"         "Кировская область"      
## [5] "Нижегородская область"   "Самарская область"

Рассмотрим, в какие месяцы наблюдаются максимальные значения исследуемых нами показателей:

## Построим график, отображающий динамику значений исследуемых показателей по месяцам :
top_long <- top %>%
  select(1:6)

top_long <- pivot_longer(data = top_long, cols = 4:6, names_to = "Индекс", values_to = "Значение")

top_long %>% 
  ggplot(aes(x = Месяц, y = Значение, group = Индекс, color = Индекс)) +
  geom_line(size = 0.8) +
  geom_point(size = 1.2) +
  labs(title = "Динамика показателей СберИндекс наиболее перспективных регионов \n в 2020 году",
       x = "Месяц",
       y = "Значение") +
  scale_colour_manual(values = c("tomato", "steelblue", "palegreen3")) +
  facet_wrap(facets = ~`Регион`) +
  theme(axis.text.x = element_text(angle = 90, vjust = 0.5, hjust=1),
        panel.background = element_rect(fill = "white", colour = "grey50"),
        panel.grid.major = element_line(color = "grey95"),
        panel.grid.minor = element_line(color = "grey98"))

Таким образом, наиболее перспективными для открытия филиалов банка являются следующие регионы (для каждого региона указан период, в котором наблюдались максимальные значения исследуемых показателей):

Воронежская область
- Доля безналичных платежей: ноябрь
- Индекс потребительской активности: сентябрь
- Количество внутренних туристов: январь - март, август - сентябрь
Калининградская область
- Доля безналичных платежей: май, ноябрь
- Индекс потребительской активности: июнь - октябрь
- Количество внутренних туристов: январь - февраль, июль - октябрь
Камчатский край
- Доля безналичных платежей: ноябрь
- Индекс потребительской активности: февраль, июль - октябрь
- Количество внутренних туристов: март, август - сентябрь, декабрь
Кировская область
- Доля безналичных платежей: ноябрь
- Индекс потребительской активности: февраль, июнь - октябрь
- Количество внутренних туристов: январь - март, август - сентябрь
Нижегородская область
- Доля безналичных платежей: май, ноябрь
- Индекс потребительской активности: март, июль - октябрь
- Количество внутренних туристов: январь - февраль, август - сентябрь
Самарская область
- Доля безналичных платежей: ноябрь
- Индекс потребительской активности: март, июнь - ноябрь
- Количество внутренних туристов: январь - март, август

# График обобщенных показателей по всем выбранным регионам:
top %>% 
  group_by(Месяц) %>% 
  summarise(`Индекс потребительской активности` = mean(`Индекс потребительской активности`),
            `Доля безналичных платежей` = mean(`Доля безналичных платежей`),
            `Количество внутренних туристов` = mean(`Количество внутренних туристов`)) %>%
  pivot_longer(cols = 2:4, names_to = "Индекс", values_to = "Значение") %>% 
  ggplot(aes(x = Месяц, y = Значение, group = Индекс, color = Индекс)) +
  geom_line(size = 0.8) +
  geom_point(size = 1.2) +
  labs(title = "Динамика средних значений показателей СберИндекс \n выбранных регионов в 2020 году",
       x = "Месяц",
       y = "Значение") +
  scale_colour_manual(values = c("tomato", "steelblue", "palegreen3")) +
  theme(axis.text.x = element_text(angle = 90, vjust = 0.5, hjust=1),
        panel.background = element_rect(fill = "white", colour = "grey50"),
        panel.grid.major = element_line(color = "grey95"),
        panel.grid.minor = element_line(color = "grey98"))

Обобщая для всех регионов, месяцы, в которых показатели были максимальными:

Доля безналичных платежей: май, ноябрь
Индекс потребительской активности: февраль - март, июнь - октябрь
Количество внутренних туристов: январь - март, июль - октябрь

Поиск взаимосвязей:

Проверим наличие связи между двумя показателями: Доля безналичных платежей и Количество внутренних туристов. Для этого построим график и вычислим коэффициент корреляции Пирсона.

top %>% 
  ggplot(aes(x = `Доля безналичных платежей`, y = `Количество внутренних туристов`)) +
  geom_point(size = 3, shape = 21, color = "grey30", fill = "steelblue", alpha = 0.6) +
  labs(title = "Характер связи между двумя показателями") +
  theme_minimal()

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

# Коэффициент корреляции Пирсона:
cor.test(top$`Доля безналичных платежей`, top$`Количество внутренних туристов`)

## 
##  Pearson's product-moment correlation
## 
## data:  top$`Доля безналичных платежей` and top$`Количество внутренних туристов`
## t = 0.045248, df = 70, p-value = 0.964
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  -0.2265443  0.2367800
## sample estimates:
##         cor 
## 0.005408091

Значение p-value значительно превышает 0.05, поэтому на 5%-ном уровне значимости мы не можем отвергнуть нулевую гипотезу о независимости двух показателей. Другими словами, мы не можем говорить о наличии статистически значимой связи между долей безналичных платежей и количеством внутренних туристов.

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

Проверим наличие связи между общим количеством одобренных заявок на кредит и числом ипотечных сделок:

cor.test(as.numeric(top$`Всего одобренных заявок`), as.numeric(top$`Всего ипотечных сделок`), method = "kendall")

## 
##  Kendall's rank correlation tau
## 
## data:  as.numeric(top$`Всего одобренных заявок`) and as.numeric(top$`Всего ипотечных сделок`)
## z = 6.6059, p-value = 3.952e-11
## alternative hypothesis: true tau is not equal to 0
## sample estimates:
##       tau 
## 0.7130607

Мы получили значение p-value меньшее 0.05 - можно отвергнуть нулевую гипотезу о равенстве нулю коэффициента корреляции между двумя исследуемыми показателями. Другими словами, на 5%-ном уровне значимости можно говорить о наличии сильной положительной связи между общим количеством одобренных заявок на кредит и числом ипотечных сделок.

# Вычислим коэффициент корреляции Пирсона для доли безналичных платежей и  количеством заявок на кредиты онлайн:

cor.test(top$`Доля безналичных платежей`, top$`Доля онлайн-заявок`, method = "pearson")

## 
##  Pearson's product-moment correlation
## 
## data:  top$`Доля безналичных платежей` and top$`Доля онлайн-заявок`
## t = 2.0461, df = 70, p-value = 0.04451
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  0.006224355 0.444743383
## sample estimates:
##       cor 
## 0.2375504

При расчете коэффициента корреляции Пирсона выявлена статистически значимая положительная слабая связь на уровне значимости 5% (p-value < 0.05).

Таким образом, анализ взаимосвязей позволил сделать следующие выводы:

Нет достаточных оснований говорить о наличии связи между количеством внутренних туристов в регионе и долей безналичных платежей.
Общее количество одобренных заявок на кредит сильно положительно коррелирует с общим числом ипотечных сделок.
Между долей безналичных платежей и долей онлайн-заявок на кредит существует слабая статистически значимая связь.

# Введем новый столбец Сезон:
top <- top %>% 
  mutate(Сезон = as.factor(case_when(top$Месяц %in% c("Декабрь", "Январь", "Февраль") ~ "Зима",
                              top$Месяц %in% c("Март", "Апрель", "Май") ~ "Весна",
                              top$Месяц %in% c("Июнь", "Июль", "Август") ~
"Лето",
                              top$Месяц %in% c("Сентябрь", "Октябрь", "Ноябрь") ~ "Осень")), .after = `Квартал`)

# Посторим график распределени ипотечных сделок по сезонам:
top %>% 
  ggplot(aes(x = `Сезон`, fill = `Всего ипотечных сделок`)) +
  geom_bar() +
  labs(title = "Распределение количества ипотечных сделок \n по сезонам",
       x = "Сезон",
       y = "Количество") +
  scale_fill_manual(values = c("palegreen3", "wheat", "lightpink", "steelblue")) +
  theme_minimal()

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

# Боксплот доля онлайн-заявок на кредиты ~ сезон:
top %>% 
  ggplot(aes(x = `Сезон`, y = `Доля онлайн-заявок`, fill = Сезон)) +
  geom_boxplot() +
  labs(title = "Количество онлайн-заявок в разное время года",
       x = "Сезон",
       y = "Количество") +
  scale_fill_manual(values = c("lightpink", "steelblue", "palegreen3", "coral")) +
  stat_summary(fun.y = "mean", size = 0.3, shape = 21, color = "black", stroke = 0.5, alpha = 0.8) +
  theme_minimal()

# Из боксплотов не очень понятно, особенно если отдельно проставить средние значения, влияет ли существенно время год на количество онлайн-заявок. Нужно проверить, используя статистический критерий. Так как у нас есть факторная переменная Сезон, разбивающая наблюдения на четрые выборки, для сравнения средних получившихся выборок будем использовать критерий Краскела-Уоллиса:
kruskal.test(top$`Доля онлайн-заявок` ~ top$Сезон)

## 
##  Kruskal-Wallis rank sum test
## 
## data:  top$`Доля онлайн-заявок` by top$Сезон
## Kruskal-Wallis chi-squared = 4.7869, df = 3, p-value = 0.1881

p-value > 0.05, поэтому мы не можем отвергнуть нулевую гипотезу о равенстве средних в каждой выборке и принять альтернативную гипотезу о том, что по крайней мере в одной выборке результаты отличаются. Таким образом, на 5%-ном уровне значимости у нас нет оснований считать, что доля онлайн-сделок на вторичное жильё зависит от времени года.

Проверим, различается ли доля заявок на вторичное жилья в разное время года:

# Боксплот доля ипотечных сделок на вторичное жилье ~ сезон:
top %>% 
  ggplot(aes(x = `Сезон`, y = `Доля сделок, вторичка`, fill = Сезон)) +
  geom_boxplot(outlier.fill = "coral",
               outlier.shape = 21,
               outlier.size = 1.8) +
  labs(title = "Количество ипотечных сделок на вторичное жилье \n в разное время года",
       x = "Сезон",
       y = "Доля сделок - вторичка") +
  scale_fill_manual(values = c("lightpink", "steelblue", "palegreen3", "coral")) +
  stat_summary(fun.y = "mean", size = 0.3, shape = 21, color = "black", stroke = 0.5, alpha = 0.8) +
  theme_minimal()

График показывает, что доля сделок на вторичное жилье, совершенных весной, летом и осенью, примерно одинаково. В то же время летом оно определенно ниже. Попробуем проверить формально с помощью статистического критерия.

# Так как критерий Краскела-Уоллиса применяется для тестирования гипотезы о равенстве средних, предварительно удалим нетипичные значения: 
sec_out <- top %>% 
  filter(Сезон == "Весна") %>% 
  select(`Доля сделок, вторичка`) %>% 
  pull()
iqr <- IQR(sec_out) 
q1 <- quantile(sec_out, 0.25)
q3 <- quantile(sec_out, 0.75)
lower_b <- q1 - 1.5 * iqr
upper_b <- q3 + 1.5 * iqr
sec_out[sec_out < lower_b | sec_out > upper_b]

## [1] 0.46

# Есть одно нетипично низкое значение. Попробуем исключить его перед тестированием гипотезы
top_clean <- top %>% 
  filter(`Доля сделок, вторичка` != 0.46)

# Оказалось, что удаленное значение совсем незначительно влияет на резльтат, В любом случае, p-value на порядок меньше 0.05:
kruskal.test(top$`Доля сделок, вторичка` ~ top$Сезон)

## 
##  Kruskal-Wallis rank sum test
## 
## data:  top$`Доля сделок, вторичка` by top$Сезон
## Kruskal-Wallis chi-squared = 14.122, df = 3, p-value = 0.002743

kruskal.test(top_clean$`Доля сделок, вторичка` ~ top_clean$Сезон)

## 
##  Kruskal-Wallis rank sum test
## 
## data:  top_clean$`Доля сделок, вторичка` by top_clean$Сезон
## Kruskal-Wallis chi-squared = 16.254, df = 3, p-value = 0.001006

Да, действительно, в случае ипотеки на вторичное жилье влияние времен года существенно: значение p-value < 0.05, поэтому мы можем отвергнуть нулевую гипотезу о равенстве средних значений в выборках и утверждать, что на уровне значимости 5% доля ипотечных сделок, совершенных по крайней мере в одном из сезонов, статистически отличается от остальных.

Обобщая вышесказанное, можно отметить, что количество ипотечных сделок неодинаково в разное время года. Больше всего сделок заключается в осенние месяцы. При этом не обнаружено статистически значимых различий между средними значениями количества онлайн-заявок на ипотеку в разное время года. В то же время имеет место влияние времени года на долю заявок на вторичное жилье Полученные данные могут быть полезны, например, при планировании рекламы и акций ипотечных продуктов в разные месяцы.

Построение регрессионных моделей:

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

Множественная линейная регрессия:

Построим линейную модель, описывающую зависимость потребительской активности от средней заработной платы в регионе, числа активных абонентов беспроводного наземного фиксированного доступа к сети Интернет, уровня безработицы населения и времени года:

lmod <- lm(data = top, formula = `Индекс потребительской активности` ~ `Среднемесячная заработная плата` + `Число активных абонентов беспроводного доступа к сети` + `Уровень безработицы населения` + `Сезон`)

summary(lmod)

## 
## Call:
## lm(formula = `Индекс потребительской активности` ~ 
##     `Среднемесячная заработная плата` + 
##         `Число активных абонентов беспроводного доступа к сети` + 
##         `Уровень безработицы населения` + 
##         Сезон, data = top)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -16.3064  -2.7598  -0.3264   2.8794  15.8226 
## 
## Coefficients:
##                                                           Estimate Std. Error
## (Intercept)                                              5.606e+01  6.905e+00
## `Среднемесячная заработная плата`                       -1.847e-05  5.712e-05
## `Число активных абонентов беспроводного доступа к сети`  4.891e-04  4.892e-04
## `Уровень безработицы населения`                          8.516e-01  1.008e+00
## СезонЗима                                                9.749e+00  2.262e+00
## СезонЛето                                                1.480e+01  2.293e+00
## СезонОсень                                               1.457e+01  2.307e+00
##                                                         t value Pr(>|t|)    
## (Intercept)                                               8.118 1.80e-11 ***
## `Среднемесячная заработная плата`                        -0.323    0.747    
## `Число активных абонентов беспроводного доступа к сети`   1.000    0.321    
## `Уровень безработицы населения`                           0.845    0.401    
## СезонЗима                                                 4.309 5.66e-05 ***
## СезонЛето                                                 6.454 1.58e-08 ***
## СезонОсень                                                6.316 2.76e-08 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 6.74 on 65 degrees of freedom
## Multiple R-squared:  0.4969, Adjusted R-squared:  0.4505 
## F-statistic:  10.7 on 6 and 65 DF,  p-value: 3.045e-08

Для этой модели статистически значимыми на любом разумном уровне значимости являются только коэффициенты при показателях, связанных с временами года. Иными словами, в разное время года потребительская активность достоверно различается. Базовой категорией в нашем случая является весна, поэтому коэффициенты напротив других сезонов показывают различие относительно весеннего периода. Интерпретировать результат можно следующим образом (при предположении, что другие условия не изменяются):

Индекс потребительской активности зимой в среднем на 9.75 выше, по сравнению с весной.
Летом индекс потребительской активности выше в среднем на 14.8, чем весной.
Осенью потребительская активность в среднем выше на 14.6 по сравнению с весенним периодом.

Следует обратить внимание на R-квадрат полученной модели. Значение R-squared = 0.4969 говорит о том, что наша модель объясняет примерно 50% дисперсии зависимой переменной, а значит, имеет среднюю предсказательную силу.

Построим модель, описывающую влияние среднемесячной заработной платы, числа активных пользователей беспроводного интернета, уровня безработицы и времени года на долю безналичных платежей:

lmod2 <- lm(data = top, formula = `Доля безналичных платежей` ~ `Среднемесячная заработная плата` + `Число активных абонентов беспроводного доступа к сети` + `Уровень безработицы населения` + `Сезон`)
summary(lmod2)

## 
## Call:
## lm(formula = `Доля безналичных платежей` ~ 
##     `Среднемесячная заработная плата` + 
##         `Число активных абонентов беспроводного доступа к сети` + 
##         `Уровень безработицы населения` + 
##         Сезон, data = top)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -3.6203 -1.2485 -0.2735  1.1891  3.5084 
## 
## Coefficients:
##                                                           Estimate Std. Error
## (Intercept)                                              4.856e+01  1.774e+00
## `Среднемесячная заработная плата`                        8.194e-05  1.467e-05
## `Число активных абонентов беспроводного доступа к сети` -3.714e-04  1.257e-04
## `Уровень безработицы населения`                          1.134e+00  2.590e-01
## СезонЗима                                                6.025e-01  5.811e-01
## СезонЛето                                               -8.082e-01  5.889e-01
## СезонОсень                                               1.388e+00  5.927e-01
##                                                         t value Pr(>|t|)    
## (Intercept)                                              27.378  < 2e-16 ***
## `Среднемесячная заработная плата`                         5.585 4.94e-07 ***
## `Число активных абонентов беспроводного доступа к сети`  -2.956  0.00434 ** 
## `Уровень безработицы населения`                           4.377 4.46e-05 ***
## СезонЗима                                                 1.037  0.30370    
## СезонЛето                                                -1.372  0.17467    
## СезонОсень                                                2.341  0.02231 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1.731 on 65 degrees of freedom
## Multiple R-squared:  0.6188, Adjusted R-squared:  0.5837 
## F-statistic: 17.59 on 6 and 65 DF,  p-value: 5.579e-12

В этом случае коэффициенты при количественных показателях, а также коэффициент при показателе Сезон-Осень, получились статистически значимыми на уровне значимости 5%. Проинтерпретируем влияние каждого фактора на зависимую переменную (при прочих равных условиях):

При увеличении среднемесячной заработной платы на 1 доля безналичных платежей также увеличивается на 0.00008194.
Число активных абонентов беспроводного интернета слабо отрицательно влияет на долю безналичных платежей: при увеличении на 1 зависимая переменная уменьшается на 0.0003714.
При увеличении уровня безработицы на 1 доля безналичных платежей увеличивается в среднем на 1.134.
Осенью доля безналичных платежей в среднем на 1.388 выше, чем весной. Для остальных сезонов статистически значимой разницы на выбранном нами уровне значимости = 0.05 не выявлено.

R-квадрат полученной модели выше, чем у предыдущей, и равен 0.6188. Это означает, что она может объяснять примерно 62% дисперсии зависимой переменной и обладает средней прогностической мощностью.

Логистическая регрессия:

# Добавим в датафрейм столбец, показыващие преобладание количества онлайн-заявок на кредит в офисе банка или онлайн, а также аналогичный столбец, где 1 соответсвует регионам с общим числом ипотечных сделок > 500:

top <- top %>% 
  mutate(is_online = ifelse(`Доля онлайн-заявок` > `Доля заявок в офисе банка`, 1, 0), .after = `Доля заявок в офисе банка`)

top$is_online <- as.factor(top$is_online)

top <- top %>% 
  mutate(hyp_love = ifelse(`Всего ипотечных сделок` %in% c("500 - 1 000", "> 1 000"), 1, 0), .after = `Всего ипотечных сделок`)

top$hyp_love <- as.factor(top$hyp_love)

Построим первую регрессионную модель, которая отображает зависимость склонности населения предпочитать онлайн-услуги банка от средней заработной платы в регионе, числа активных абонентов беспроводного наземного фиксированного доступа к сети Интернет и уровня безработицы населения.

# Так как зависимая переменная у нас может принимать только два значения - 0 и 1, с которыми мы для удобства связали склонность выбирать онлайн услуги, будем использовать логистическую регрессионную модель:
logmod <- glm(data = top,
              formula = top$is_online ~ `Среднемесячная заработная плата` + `Число активных абонентов беспроводного доступа к сети` + `Уровень безработицы населения` + `Сезон`, family = "binomial")

summary(logmod)

## 
## Call:
## glm(formula = top$is_online ~ `Среднемесячная заработная плата` + 
##     `Число активных абонентов беспроводного доступа к сети` + 
##     `Уровень безработицы населения` + 
##     Сезон, family = "binomial", data = top)
## 
## Deviance Residuals: 
##      Min        1Q    Median        3Q       Max  
## -1.83091  -0.55433  -0.24683  -0.00004   1.99553  
## 
## Coefficients:
##                                                           Estimate Std. Error
## (Intercept)                                             -7.228e+00  3.509e+00
## `Среднемесячная заработная плата`                        8.659e-05  3.353e-05
## `Число активных абонентов беспроводного доступа к сети`  5.246e-05  2.479e-04
## `Уровень безработицы населения`                          5.785e-01  4.626e-01
## СезонЗима                                               -1.600e+00  1.077e+00
## СезонЛето                                               -1.978e+01  2.263e+03
## СезонОсень                                              -7.198e-01  9.172e-01
##                                                         z value Pr(>|z|)   
## (Intercept)                                              -2.060  0.03943 * 
## `Среднемесячная заработная плата`                         2.583  0.00981 **
## `Число активных абонентов беспроводного доступа к сети`   0.212  0.83239   
## `Уровень безработицы населения`                           1.251  0.21108   
## СезонЗима                                                -1.486  0.13740   
## СезонЛето                                                -0.009  0.99302   
## СезонОсень                                               -0.785  0.43255   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 70.935  on 71  degrees of freedom
## Residual deviance: 44.304  on 65  degrees of freedom
## AIC: 58.304
## 
## Number of Fisher Scoring iterations: 18

exp(coef(logmod))

##                                             (Intercept) 
##                                            7.262714e-04 
##                       `Среднемесячная заработная плата` 
##                                            1.000087e+00 
## `Число активных абонентов беспроводного доступа к сети` 
##                                            1.000052e+00 
##                         `Уровень безработицы населения` 
##                                            1.783306e+00 
##                                               СезонЗима 
##                                            2.018236e-01 
##                                               СезонЛето 
##                                            2.561146e-09 
##                                              СезонОсень 
##                                            4.868326e-01

Выдача R показывает, что статистически значимым на уровне 5% оказался только коэффициент при одном факторе - среднемесячной заработной плате. Однако, его влияние незначительно (при увеличении средней зарплаты на 1 рубль вероятность склонности выбирать онлайн-услуги увеличивается в среднем в 1.000087 раз при прочих равных). Чтобы было более показательно, попробуем изменить столбец Средняя заработная плата таким образом, чтобы значения отображались в тысячах рублей.

top <- top %>% 
  mutate(`Средняя заработная плата тысячи` = `Среднемесячная заработная плата` / 1000, .after = `Среднемесячная заработная плата`)

# Продублируем модель с новой размерностью и посмотрим, как это повлияет на коэффициенты при показателях:
logmod1000 <- glm(data = top,
              formula = top$is_online ~ `Средняя заработная плата тысячи` + `Число активных абонентов беспроводного доступа к сети` + `Уровень безработицы населения` + `Сезон`, family = "binomial")

summary(logmod1000)

## 
## Call:
## glm(formula = top$is_online ~ `Средняя заработная плата тысячи` + 
##     `Число активных абонентов беспроводного доступа к сети` + 
##     `Уровень безработицы населения` + 
##     Сезон, family = "binomial", data = top)
## 
## Deviance Residuals: 
##      Min        1Q    Median        3Q       Max  
## -1.83091  -0.55433  -0.24683  -0.00004   1.99553  
## 
## Coefficients:
##                                                           Estimate Std. Error
## (Intercept)                                             -7.228e+00  3.509e+00
## `Средняя заработная плата тысячи`                        8.659e-02  3.353e-02
## `Число активных абонентов беспроводного доступа к сети`  5.246e-05  2.479e-04
## `Уровень безработицы населения`                          5.785e-01  4.626e-01
## СезонЗима                                               -1.600e+00  1.077e+00
## СезонЛето                                               -1.978e+01  2.263e+03
## СезонОсень                                              -7.198e-01  9.172e-01
##                                                         z value Pr(>|z|)   
## (Intercept)                                              -2.060  0.03943 * 
## `Средняя заработная плата тысячи`                         2.583  0.00981 **
## `Число активных абонентов беспроводного доступа к сети`   0.212  0.83239   
## `Уровень безработицы населения`                           1.251  0.21108   
## СезонЗима                                                -1.486  0.13740   
## СезонЛето                                                -0.009  0.99302   
## СезонОсень                                               -0.785  0.43255   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 70.935  on 71  degrees of freedom
## Residual deviance: 44.304  on 65  degrees of freedom
## AIC: 58.304
## 
## Number of Fisher Scoring iterations: 18

exp(coef(logmod1000))

##                                             (Intercept) 
##                                            7.262714e-04 
##                       `Средняя заработная плата тысячи` 
##                                            1.090449e+00 
## `Число активных абонентов беспроводного доступа к сети` 
##                                            1.000052e+00 
##                         `Уровень безработицы населения` 
##                                            1.783306e+00 
##                                               СезонЗима 
##                                            2.018236e-01 
##                                               СезонЛето 
##                                            2.561146e-09 
##                                              СезонОсень 
##                                            4.868326e-01

Мы получили вполне ожидаемый результат: статистически значимым остался только коэффициент при показателе среднемесячной заработной платы, при этом сам коэффициент немного поменялся, а его интерпретация стала более наглядной: при увеличении средней заработной платы на 1 тысячу рублей вероятность того, что потребитель предпочтет оформлять заявку на кредит онлайн возрастает в среднем в 1.09 раз при прочих равных условиях.

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

logmod2 <- glm(data = top, 
               formula = hyp_love ~ `Средняя заработная плата тысячи` +
                 `Число активных абонентов беспроводного доступа к сети` +
                 `Уровень безработицы населения` +
                 `Сезон`, family = "binomial")

summary(logmod2)

## 
## Call:
## glm(formula = hyp_love ~ `Средняя заработная плата тысячи` + 
##     `Число активных абонентов беспроводного доступа к сети` + 
##     `Уровень безработицы населения` + 
##     Сезон, family = "binomial", data = top)
## 
## Deviance Residuals: 
##     Min       1Q   Median       3Q      Max  
## -2.2149  -0.2453   0.2015   0.4761   1.6768  
## 
## Coefficients:
##                                                           Estimate Std. Error
## (Intercept)                                              2.6584777  2.8820906
## `Средняя заработная плата тысячи`                       -0.0748965  0.0337003
## `Число активных абонентов беспроводного доступа к сети`  0.0007899  0.0003221
## `Уровень безработицы населения`                         -0.1322649  0.4184972
## СезонЗима                                                0.1286347  0.9424217
## СезонЛето                                                1.0362722  1.0489430
## СезонОсень                                               2.4828064  1.3588832
##                                                         z value Pr(>|z|)  
## (Intercept)                                               0.922   0.3563  
## `Средняя заработная плата тысячи`                        -2.222   0.0263 *
## `Число активных абонентов беспроводного доступа к сети`   2.452   0.0142 *
## `Уровень безработицы населения`                          -0.316   0.7520  
## СезонЗима                                                 0.136   0.8914  
## СезонЛето                                                 0.988   0.3232  
## СезонОсень                                                1.827   0.0677 .
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 88.632  on 71  degrees of freedom
## Residual deviance: 47.110  on 65  degrees of freedom
## AIC: 61.11
## 
## Number of Fisher Scoring iterations: 6

exp(coef(logmod2))

##                                             (Intercept) 
##                                              14.2745423 
##                       `Средняя заработная плата тысячи` 
##                                               0.9278395 
## `Число активных абонентов беспроводного доступа к сети` 
##                                               1.0007902 
##                         `Уровень безработицы населения` 
##                                               0.8761089 
##                                               СезонЗима 
##                                               1.1372746 
##                                               СезонЛето 
##                                               2.8186899 
##                                              СезонОсень 
##                                              11.9748238

Выдача R показывает, что статистически значимыми на уровне 5% являются показатели Средняя заработная плата и Число активных абонентов беспроводного доступа к сети. При увеличении средней заработной платы на 1 тысячу рублей склонность потребителя оформлять ипотеку увеличивается в 0.93 раза (или уменьшается на 7%), при условии, что значения других показателей не изменяются. Число активных абонентов беспроводного интернета положительно связано с зависимой переменной: при увеличении первого на 1 склонность оформлять ипотеку увеличивается в 1.0008 раз (или ~ в 1.08 раз при увеличении на 100).

Оценка качества полученных моделей:

Оценка качества множественной линейной регрессии:

Сначала поработаем с моделями множественной линейной регрессии. В предыдущем пункте мы уже оценили R-квадрат полученных моделей, теперь протестируем их более детально: проверим, что выполняется следующий ряд условий: * В модели нет независимых переменных, которые сильно скоррелированы друг с другом (отсутствие мультиколлинеарности); * Остатки в модели распределены случайно, нет выраженной зависимости между остатками модели и независимыми переменными (проблема гетероскедастичности); * В данных нет влиятельных наблюдений, искажающих плоскость, которую описывает уравнение регрессии.

Выявление мультиколлинеарности:

Для этого построим корреляционную матрицу для независимых переменных:

# Вначале преобразуем факторный показатель Сезон:
top <- top %>% 
  mutate(Season_numeric = as.numeric(Сезон), .after = Сезон)
  
mltcor <- top %>% 
  select(`Среднемесячная заработная плата`, `Число активных абонентов беспроводного доступа к сети`, `Уровень безработицы населения`, Season_numeric)
cor(mltcor)

##                                                       Среднемесячная заработная плата
## Среднемесячная заработная плата                                            1.00000000
## Число активных абонентов беспроводного доступа к сети                     -0.50399291
## Уровень безработицы населения                                             -0.38481489
## Season_numeric                                                            -0.01735107
##                                                       Число активных абонентов беспроводного доступа к сети
## Среднемесячная заработная плата                                                                -0.503992914
## Число активных абонентов беспроводного доступа к сети                                           1.000000000
## Уровень безработицы населения                                                                  -0.008317665
## Season_numeric                                                                                  0.045851236
##                                                       Уровень безработицы населения
## Среднемесячная заработная плата                                        -0.384814885
## Число активных абонентов беспроводного доступа к сети                  -0.008317665
## Уровень безработицы населения                                           1.000000000
## Season_numeric                                                          0.255818292
##                                                       Season_numeric
## Среднемесячная заработная плата                          -0.01735107
## Число активных абонентов беспроводного доступа к сети     0.04585124
## Уровень безработицы населения                             0.25581829
## Season_numeric                                            1.00000000

Можно отметить корреляцию средней силы между числом активных абонентов беспроводного доступа к сети и среднемесячной заработной платой, а также среднемесячной заработной платой и уровнем безработицы населения. Максимальное значение по модулю коэффициента корреляции = 0.5. На мой взгляд, это допустимое значение, не требующее корректировки модели. Если бы мы встретили значения по модулю выше 0.7-0.8, то это уже могло бы свидетельствовать о возможных проблемах с качеством получаемых оценок. Вторая модель содержит такой же набор независимых переменных, поэтому повторять этот тест для нее мы не будем.

Выявление гетероскедастичности:

Для каждой модели построим диаграмму рассеивания вида предсказанные значения - остатки:

# Добавим в основной датафрейм предсказанные значения и остатки моделей lmod и lmod2:
top <- top %>% 
  mutate(lmod_predicted = lmod$fitted.values)
top <- top %>% 
  mutate(lmod_residuals = lmod$residuals)
top <-  top %>% 
  mutate(lmod2_predicted = lmod2$fitted.values)
top <- top %>% 
  mutate(lmod2_residuals = lmod2$residuals)

# Построим график для первой модели:
top %>% 
  ggplot(aes(x = lmod_predicted, y = lmod_residuals)) +
  geom_point(color = "steelblue", size = 2) +
  geom_hline(yintercept = 0, color = "tomato") +
  labs(title = "Диаграмма рассеивания для предсказанных значений и остатков lmod", 
       x = "Fitted values",
       y = "Residuals") +
  theme_minimal()

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

coeftest(lmod, vcov. = vcovHC(lmod, type = "HC0"))

## 
## t test of coefficients:
## 
##                                                            Estimate  Std. Error
## (Intercept)                                              5.6059e+01  7.2262e+00
## `Среднемесячная заработная плата`                       -1.8474e-05  4.4399e-05
## `Число активных абонентов беспроводного доступа к сети`  4.8906e-04  4.8575e-04
## `Уровень безработицы населения`                          8.5165e-01  1.0268e+00
## СезонЗима                                                9.7489e+00  2.5113e+00
## СезонЛето                                                1.4797e+01  2.5264e+00
## СезонОсень                                               1.4572e+01  2.7000e+00
##                                                         t value  Pr(>|t|)    
## (Intercept)                                              7.7577 7.864e-11 ***
## `Среднемесячная заработная плата`                       -0.4161 0.6787122    
## `Число активных абонентов беспроводного доступа к сети`  1.0068 0.3177582    
## `Уровень безработицы населения`                          0.8294 0.4099237    
## СезонЗима                                                3.8821 0.0002448 ***
## СезонЛето                                                5.8569 1.704e-07 ***
## СезонОсень                                               5.3971 1.020e-06 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Мы получили немного другие значения стандартных ошибок и p-value, при этом значимые факторы остались прежними, значения коэффициентов тоже (при округлении до сотых) не изменились. Таким образом, интерпретация модели остается прежней.

Для второй модели:

top %>% 
  ggplot(aes(x = lmod2_predicted, y = lmod2_residuals)) +
  geom_point(color = "steelblue", size = 2) +
  geom_hline(yintercept = 0, color = "tomato") +
  labs(title = "Диаграмма рассеивания для предсказанных значений и остатков lmod2", 
       x = "Fitted values",
       y = "Residuals") +
  theme_minimal()

Здесь также присутствует небольшая, менее выраженная, чем в первом случае, неоднородность (график имеет веретеновидную форму). Пересчитаем коэффциенты модели и сравним их с исходными:

coeftest(lmod2, vcov. = vcovHC(lmod2, type = "HC0"))

## 
## t test of coefficients:
## 
##                                                            Estimate  Std. Error
## (Intercept)                                              4.8561e+01  1.4638e+00
## `Среднемесячная заработная плата`                        8.1936e-05  1.3099e-05
## `Число активных абонентов беспроводного доступа к сети` -3.7144e-04  1.1263e-04
## `Уровень безработицы населения`                          1.1335e+00  2.0000e-01
## СезонЗима                                                6.0245e-01  5.5513e-01
## СезонЛето                                               -8.0817e-01  4.9594e-01
## СезонОсень                                               1.3876e+00  5.7528e-01
##                                                         t value  Pr(>|t|)    
## (Intercept)                                             33.1746 < 2.2e-16 ***
## `Среднемесячная заработная плата`                        6.2552 3.507e-08 ***
## `Число активных абонентов беспроводного доступа к сети` -3.2978  0.001583 ** 
## `Уровень безработицы населения`                          5.6676 3.578e-07 ***
## СезонЗима                                                1.0853  0.281818    
## СезонЛето                                               -1.6296  0.108027    
## СезонОсень                                               2.4120  0.018698 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Значимыми остались прежние коэффициенты, сами значения коэффициентов не изменились (при округлении до сотых). Соотвественно, интерпретация регрессионной модели остается прежней.

Выявление влиятельных наблюдений:

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

plot(lmod, 5)

Обычно явно влиятельные наблюдения располагаются на этом графике в правом верхнем углу и отделяются красной пунктирной линией. В нашем случае таких наблюдений не выявлено.

Повторим для второй модели:

plot(lmod2, 5)

Во втором случае также не выявлено влиятельных наблюдений.

Таким образом, мы провели ключевые тесты оценки качества и можем сделать вывод о том, что построенные нами линейные модели являются качественными, однако обладают не очень высокой прогностической силой.

Оценка качества логистических моделей:

Так как полученные нами модели являются логистическими (зависимая переменная может принимать значение 0 либо 1), оценим, насколько точно каждая модель предсказывает положительный и отрицательный исходы. Для этого сравним полученные на предыдущем этапе значения столбцов is_online и hyp_love с предсказанными значениями для каждой модели.

## Вычислим чувствительность и специфичность для первой модели:
top <- top %>% 
  mutate(is_online_predicted = ifelse(logmod$fitted.values > mean(logmod$fitted.values), 1, 0), .after = is_online)

# Верно предсказанные 0 (True Negative):
TN <- top %>% 
  filter(is_online == 0, is_online_predicted == 0) %>% 
  nrow()
# Ложно предсказанные 1 (False Positive):
FP <- top %>% 
  filter(is_online == 0, is_online_predicted == 1) %>% 
  nrow()
# Ложно предсказанные 0 (False Negative):
FN <- top %>% 
  filter(is_online == 1, is_online_predicted == 0) %>% 
  nrow()
# Верно предсказанные 1 (True Positive):
TP <- top %>% 
  filter(is_online == 1, is_online_predicted == 1) %>% 
  nrow()

# Специфичность модели:
specif <- TN / (TN + FP)

# Чувствительность:
sens <- TP / (TP + FN)

# Чувствительность и специфичность полученной модели = 0.785 и 0.793. Другими словами наша модель нашла 78.5% из всех истинно положительных исходов. Также в 79.3% случаев, модель правильно предсказала 0 среди общего числа истинно отрицательных исходов.

Построим ROC-кривую, еще один инструмент, который поможет наглядно отобразить качество исследуемой нами модели.

proc <- roc(top$is_online ~ logmod$fitted.values)

## Setting levels: control = 0, case = 1

## Setting direction: controls < cases

proc

## 
## Call:
## roc.formula(formula = top$is_online ~ logmod$fitted.values)
## 
## Data: logmod$fitted.values in 58 controls (top$is_online 0) < 14 cases (top$is_online 1).
## Area under the curve: 0.8793

plot(proc)
text(x = 1, y = 0.9, "AUC = 0.88")

Площадь под ROC-кривой (AUC) модет принимать значения от 0 до 1. Высокие значения AUC, близкие к 1, соотвествуют высокому качеству регрессионной модели. В нашем случае AUC = 0.88, в связи с чем можно говорить о достаточно высокой предсказательной силе модели.

Повторим то же самое для второй модели.

## Вычислим чувствительность и специфичность для второй модели:
top <- top %>% 
  mutate(hyp_love_predicted = ifelse(logmod2$fitted.values > mean(logmod2$fitted.values), 1, 0), .after = hyp_love)

# Верно предсказанные 0 (True Negative):
TN2 <- top %>% 
  filter(hyp_love == 0, hyp_love_predicted == 0) %>% 
  nrow()
# Ложно предсказанные 1 (False Positive):
FP2 <- top %>% 
  filter(hyp_love == 0, hyp_love_predicted == 1) %>% 
  nrow()
# Ложно предсказанные 0 (False Negative):
FN2 <- top %>% 
  filter(hyp_love == 1, hyp_love_predicted == 0) %>% 
  nrow()
# Верно предсказанные 1 (True Positive):
TP2 <- top %>% 
  filter(hyp_love == 1, hyp_love_predicted == 1) %>% 
  nrow()

# Специфичность модели:
specif2 <- TN2 / (TN2 + FP2)

# Чувствительность:
sens2 <- TP2 / (TP2 + FN2)

# Чувствительность и специфичность второй модели logmod2 оказались еще выше = 0.84 и 0.82 соответственно. Интерпретировать эти значения можно следующим образом: в 84% случаев модель правильно угадывает 1 и в 82% случаев правильно предсказывает 0.

ROC-кривая и AUC:

proc2 <- roc(top$hyp_love ~ logmod2$fitted.values)

## Setting levels: control = 0, case = 1

## Setting direction: controls < cases

proc2

## 
## Call:
## roc.formula(formula = top$hyp_love ~ logmod2$fitted.values)
## 
## Data: logmod2$fitted.values in 22 controls (top$hyp_love 0) < 50 cases (top$hyp_love 1).
## Area under the curve: 0.9118

plot(proc2)
text(x = 1.05, y = 0.9, "AUC = 0.91")

Значение AUC = 0.91, получилось еще выше, чем у предыдущей модели.

Таким образом, мы провели оценку качества для двух логистических регрессионных моделей двумя способами: рассчитав чувствительность и специфичность, а также значение AUC для ROC-кривой.

Для модели, предсказывающей склонность населения предпочитать онлайн-услуги банка:

Чувствительность = 0.785
Специфичность = 0.793
AUC = 0.88

Модель, предсказывающая склонность населения оформлять ипотеку:

Чувствительность = 0.84
Специфичность = 0.82
AUC = 0.91

Вторая модель характеризуется более высокими показателями.

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

Представление результатов:

Библиотека Stargazer позволяет в удобном виде отображать характеристики сразу нескольких регрессионных моделей.

Выведем описание полученных линейных моделей:

stargazer(lmod, lmod2, type = "html", title = "Зависимость потребительской активности от средней заработной платы в регионе, числа активных абонентов беспроводного интернета, уровня безработицы населения и времени года", style = )

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

	Dependent variable:

	`Индекс потребительской активности`	`Доля безналичных платежей`
	(1)	(2)

`Среднемесячная заработная плата`	-0.00002	0.0001^***
	(0.0001)	(0.00001)

`Число активных абонентов беспроводного доступа к сети`	0.0005	-0.0004^***
	(0.0005)	(0.0001)

`Уровень безработицы населения`	0.852	1.134^***
	(1.008)	(0.259)

СезонЗима	9.749^***	0.602
	(2.262)	(0.581)

СезонЛето	14.797^***	-0.808
	(2.293)	(0.589)

СезонОсень	14.572^***	1.388^**
	(2.307)	(0.593)

Constant	56.059^***	48.561^***
	(6.905)	(1.774)


Observations	72	72
R²	0.497	0.619
Adjusted R²	0.450	0.584
Residual Std. Error (df = 65)	6.740	1.731
F Statistic (df = 6; 65)	10.701^***	17.589^***

Note:	p<0.1; p<0.05; p<0.01

Аналогично для логистических моделей:

stargazer(logmod1000, logmod2, type = "html", title = "Зависимость склонности населения предпочитать онлайн-услуги банка и брать ипотеку от средней заработной платы, числа активных абонентов интернета и уровня безработицы населения")

Зависимость склонности населения предпочитать онлайн-услуги банка и брать ипотеку от средней заработной платы, числа активных абонентов интернета и уровня безработицы населения

	Dependent variable:

	is_online	hyp_love
	(1)	(2)

`Средняя заработная плата тысячи`	0.087^***	-0.075^**
	(0.034)	(0.034)

`Число активных абонентов беспроводного доступа к сети`	0.0001	0.001^**
	(0.0002)	(0.0003)

`Уровень безработицы населения`	0.578	-0.132
	(0.463)	(0.418)

СезонЗима	-1.600	0.129
	(1.077)	(0.942)

СезонЛето	-19.783	1.036
	(2,262.531)	(1.049)

СезонОсень	-0.720	2.483^*
	(0.917)	(1.359)

Constant	-7.228^**	2.658
	(3.509)	(2.882)


Observations	72	72
Log Likelihood	-22.152	-23.555
Akaike Inf. Crit.	58.304	61.110

Note:	p<0.1; p<0.05; p<0.01

Основные выводы:

Регионы России сильно различаются по показателям СберИндекс (доля безналичных платежей, индекс потребительской активности и количество внутренних туристов). Весной 2020 года в большинстве регионов имело место падение индекса потребительской активности и количества внутренних туристов, что, вероятно, было связано с эпидемиологической ситуацией в стране и противоэпидемическим режимом.
Такие регионы как Воронежская, Калининградская, Кировская, Нижегородская, Самарская области и Камчатский край характеризуются наиболее высокими значениями показателей СберИндекс и могут быть перспективными кандидатами для открытия филиалов TrustUs.
Максимальные значения показателей СберИндекс в целом по стране зафиксированы в течение двух периодов: февраль-март и июль-октябрь.
Не обнаружено статистически значимой связи (α = 5%) между долей безналичных платежей и количеством внутренних туристов.
Общее количество одобренных заявок на кредит сильно положительно коррелирует с общим числом ипотечных сделок.
Больше всего ипотечных сделок совершается в осенние месяцы, при этом на долю онлайн-заявок время года не влияет.
Ключевыми факторами, определяющими спрос населения на услуги банков, являются сезон, среднемесячная заработная плата и число активных абонентов беспроводного интернета.

Аналитический отчет

Плеханов Н.

2024-08-16