df <- read_csv("insurance_cost.csv")
skimr::skim(df)
| Name | df |
| Number of rows | 1338 |
| Number of columns | 7 |
| _______________________ | |
| Column type frequency: | |
| character | 3 |
| numeric | 4 |
| ________________________ | |
| Group variables | None |
Variable type: character
| skim_variable | n_missing | complete_rate | min | max | empty | n_unique | whitespace |
|---|---|---|---|---|---|---|---|
| sex | 0 | 1 | 4 | 6 | 0 | 2 | 0 |
| smoker | 0 | 1 | 2 | 3 | 0 | 2 | 0 |
| region | 0 | 1 | 9 | 9 | 0 | 4 | 0 |
Variable type: numeric
| skim_variable | n_missing | complete_rate | mean | sd | p0 | p25 | p50 | p75 | p100 | hist |
|---|---|---|---|---|---|---|---|---|---|---|
| age | 0 | 1 | 39.21 | 14.05 | 18.00 | 27.00 | 39.00 | 51.00 | 64.00 | ▇▅▅▆▆ |
| bmi | 0 | 1 | 30.66 | 6.10 | 15.96 | 26.30 | 30.40 | 34.69 | 53.13 | ▂▇▇▂▁ |
| children | 0 | 1 | 1.09 | 1.21 | 0.00 | 0.00 | 1.00 | 2.00 | 5.00 | ▇▂▂▁▁ |
| charges | 0 | 1 | 13270.42 | 12110.01 | 1121.87 | 4740.29 | 9382.03 | 16639.91 | 63770.43 | ▇▂▁▁▁ |
plot_ly(data = df,
x = ~ bmi,
y = ~ charges,
color = ~ smoker) %>%
layout(title = 'Отношение страховых застрат и индекса массы тела',
yaxis = list(title = 'Страховые затраты'),
xaxis = list(title = 'Индекс массы тела'))
plot <-
ggplot(df, aes(x = bmi, y = charges, color = smoker)) +
geom_point(alpha = 0.5) +
ggtitle('Отношение страховых застрат и индекса массы тела') +
ylab('Страховые затраты') +
xlab('Индекс массы тела')
ggplotly(plot)
library(corrplot)
df.corr <-
df %>%
select(where(is.numeric)) %>%
cor(method = "spearman")
corrplot(df.corr)
corrplot(df.corr, method = 'color')
library(caret)
dummy <- dummyVars(" ~ .", data = df, sep = "_")
df.dummy <- as_tibble(predict(dummy, newdata = df))
skimr::skim(df.dummy)
| Name | df.dummy |
| Number of rows | 1338 |
| Number of columns | 12 |
| _______________________ | |
| Column type frequency: | |
| numeric | 12 |
| ________________________ | |
| Group variables | None |
Variable type: numeric
| skim_variable | n_missing | complete_rate | mean | sd | p0 | p25 | p50 | p75 | p100 | hist |
|---|---|---|---|---|---|---|---|---|---|---|
| age | 0 | 1 | 39.21 | 14.05 | 18.00 | 27.00 | 39.00 | 51.00 | 64.00 | ▇▅▅▆▆ |
| sexfemale | 0 | 1 | 0.49 | 0.50 | 0.00 | 0.00 | 0.00 | 1.00 | 1.00 | ▇▁▁▁▇ |
| sexmale | 0 | 1 | 0.51 | 0.50 | 0.00 | 0.00 | 1.00 | 1.00 | 1.00 | ▇▁▁▁▇ |
| bmi | 0 | 1 | 30.66 | 6.10 | 15.96 | 26.30 | 30.40 | 34.69 | 53.13 | ▂▇▇▂▁ |
| children | 0 | 1 | 1.09 | 1.21 | 0.00 | 0.00 | 1.00 | 2.00 | 5.00 | ▇▂▂▁▁ |
| smokerno | 0 | 1 | 0.80 | 0.40 | 0.00 | 1.00 | 1.00 | 1.00 | 1.00 | ▂▁▁▁▇ |
| smokeryes | 0 | 1 | 0.20 | 0.40 | 0.00 | 0.00 | 0.00 | 0.00 | 1.00 | ▇▁▁▁▂ |
| regionnortheast | 0 | 1 | 0.24 | 0.43 | 0.00 | 0.00 | 0.00 | 0.00 | 1.00 | ▇▁▁▁▂ |
| regionnorthwest | 0 | 1 | 0.24 | 0.43 | 0.00 | 0.00 | 0.00 | 0.00 | 1.00 | ▇▁▁▁▂ |
| regionsoutheast | 0 | 1 | 0.27 | 0.45 | 0.00 | 0.00 | 0.00 | 1.00 | 1.00 | ▇▁▁▁▃ |
| regionsouthwest | 0 | 1 | 0.24 | 0.43 | 0.00 | 0.00 | 0.00 | 0.00 | 1.00 | ▇▁▁▁▂ |
| charges | 0 | 1 | 13270.42 | 12110.01 | 1121.87 | 4740.29 | 9382.03 | 16639.91 | 63770.43 | ▇▂▁▁▁ |
library(factoextra)
df.hc <-
df.dummy %>%
scale() %>%
dist(method = "euclidean") %>%
hclust(method = "ward.D2")
fviz_dend(df.hc, show_labels = F)
Раскрасим кластеры, волюнтаристски выбрав отсечение на 5 кластерах.
fviz_dend(df.hc, k = 5,
show_labels = F,
k_colors = "Set1",
color_labels_by_k = TRUE,
rect = TRUE)
Визуализируем PCA с идентифицированными выше кластерами.
grp <- cutree(df.hc, k = 5)
fviz_cluster(list(data = df.dummy, cluster = grp),
palette = "Set1",
ellipse.type = "convex",
show.clust.cent = FALSE,
labelsize = 0,
ggtheme = theme_minimal())
Попробуем перебор методов рассчета дистанции (все варианты параметра
method функции dist()) и методов кластеризации (рекомендованные
методы complete или ward.D2 функции hclust()). В целом
кластеризация очень неустойчива в зависимости от метода, как и было
упомянуто в лекции. Вот самое симпатичное, что получилось:
grp <-
df.dummy %>%
scale() %>%
dist(method = "canberra") %>%
hclust(method = "complete") %>%
cutree(k = 4)
fviz_cluster(list(data = df.dummy, cluster = grp),
palette = "Set1",
ellipse.type = "convex",
show.clust.cent = FALSE,
labelsize = 0,
ggtheme = theme_minimal())
Построим heat map с кластеризацией по первому использованому нами методу.
library(pheatmap)
pheatmap(scale(df.dummy))
library(FactoMineR)
library(ggbiplot)
df.pca <- prcomp(df.dummy, scale = T)
fviz_eig(df.pca)
fviz_pca_var(df.pca, col.var = "contrib")
fviz_pca_var(df.pca,
select.var = list(contrib = 4),
col.var = "contrib")
fviz_contrib(df.pca, choice = "var", axes = 1)
fviz_contrib(df.pca, choice = "var", axes = 2)
fviz_contrib(df.pca, choice = "var", axes = 3)
ggbiplot(df.pca, scale = 0, alpha = 0.1)
ggbiplot(df.pca,
scale = 0,
groups = as.factor(df$smoker),
ellipse = T,
alpha = 0.1)
ggbiplot(df.pca,
scale = 0,
groups = as.factor(df$sex),
ellipse = T,
alpha = 0.1)
Выводы:
- PC1 объясняет около 23% дисперсии. Наибольший вклад в PC1 вносит статус по курению и страховым затратам.
- PC2 объясняет около 16% дисперсии и почти полностью состоит из вариации по полу.
- PC3 объясняет около 13% дисперсии и наибольший вклад в нее вносят регион клиентов и BMI.
Поиграемся с новыми категориальными переменными.
df.2 <-
df %>%
mutate(age_group = case_when(age < 35 ~ "age: 18-34",
age >= 50 ~ "age: 50+",
age >= 35 & age < 50 ~ "age: 35-49")) %>%
mutate(obesity = case_when(bmi < 30 ~ "BMI < 30",
bmi >= 30 ~ "BMI >= 30")) %>%
mutate(smoker_sex = paste0(sex, ", ", smoker, " smoker"))
ggbiplot(df.pca,
scale=0,
groups = df.2$age_group,
ellipse = F,
alpha = 0.1)
ggbiplot(df.pca,
scale=0,
groups = df.2$obesity,
ellipse = F,
alpha = 0.1)
ggbiplot(df.pca,
scale=0,
groups = df.2$smoker_sex,
ellipse = T,
alpha = 0.1)
Попробуем сначала исключить некоторые из dummy-переменных. "Качество" PCA будем оценивать по проценту объясненной дисперсии (scree plot).
df.pca <-
df.dummy %>%
prcomp(scale = T)
p1 <-
fviz_eig(df.pca)+
ggtitle("Nothing excluded") +
ylim(0, 35)
df.pca <-
df.dummy %>%
select(-contains("smoker")) %>%
prcomp(scale = T)
p2 <-
fviz_eig(df.pca)+
ggtitle("Smoker dummies excluded") +
ylim(0, 35)
df.pca <-
df.dummy %>%
select(-contains("region")) %>%
prcomp(scale = T)
p3 <-
fviz_eig(df.pca)+
ggtitle("Region dummies excluded") +
ylim(0, 35)
df.pca <-
df.dummy %>%
select(-contains("sex")) %>%
prcomp(scale = T)
p4 <-
fviz_eig(df.pca)+
ggtitle("Sex dummies excluded") +
ylim(0, 35)
ggarrange(p1, p2, p3, p4, ncol = 2, nrow = 2)
Итого, больше всего повышает качество PCA удаление информации о регионах.
Теперь создадим новые dummy-переменные для новых факторных переменных, созданных в задании 10-11.
dummy <- dummyVars(" ~ .", data = df.2, sep = "_")
df.2.dummy <- as_tibble(predict(dummy, newdata = df.2))
df.pca <-
df.dummy %>%
prcomp(scale = T)
p1 <-
fviz_eig(df.pca)+
ggtitle("Nothing added") +
ylim(0, 35)
df.pca <-
df.2.dummy %>%
select(-contains("age_group"), -contains("obesity")) %>%
prcomp(scale = T)
p2 <-
fviz_eig(df.pca)+
ggtitle("Smoker_sex dummies added") +
ylim(0, 35)
df.pca <-
df.2.dummy %>%
select(-contains("smoker_sex"), -contains("obesity")) %>%
prcomp(scale = T)
p3 <-
fviz_eig(df.pca)+
ggtitle("Age_group dummies added") +
ylim(0, 35)
df.pca <-
df.2.dummy %>%
select(-contains("age_group"), -contains("smoker_sex")) %>%
prcomp(scale = T)
p4 <-
fviz_eig(df.pca)+
ggtitle("Obesity dummies added") +
ylim(0, 35)
ggarrange(p1, p2, p3, p4, ncol = 2, nrow = 2)
Итого, больше всего из рассмотренных вариантов повышает качество PCA добавление переменной smoker_sex, которая по сути есть комбинация переменных smoker и sex.
Теперь совместим оба подхода (добавление и исключение dummy-переменных).
dummy <- dummyVars(" ~ .", data = df.2, sep = "_")
df.2.dummy <- as_tibble(predict(dummy, newdata = df.2))
df.pca <-
df.dummy %>%
prcomp(scale = T)
p1 <-
fviz_eig(df.pca)+
ggtitle("Nothing added or excluded") +
ylim(0, 35)
df.pca <-
df.2.dummy %>%
select(-contains("age_group"), -contains("obesity"), -contains("region")) %>%
prcomp(scale = T)
p2 <-
fviz_eig(df.pca)+
ggtitle("Smoker_sex added, region excluded") +
ylim(0, 35)
ggarrange(p1, p2, ncol = 1, nrow = 2)
Вуаля! Вы великолепны. Суммарно PC1 и PC2 объясняли порядка 40% дисперсии, а после наших манипуляций это число повысилось до 60%.
Можно сделать такие выводы:
- Комбинация переменных smoker и sex лучше объясняет вариабельность данных, чем эти переменные по отдельности. Можно подозревать, что именно их комбинация - залог разброса в данных (что мы, в принципе, видели и по кластеризации по фактору smoker_sex).
- Переменная region вносит дополнительную вариабельность в данные. Поэтому если ее исключить, объяснять придется меньше, а потому качество PCA растет.





















