# подключение библиотек
library (factoextra)
library (cluster)
library(ggplot2)
library(dplyr)
library(reshape2)
library(RColorBrewer)

data <- read.csv("Kvartiry_Ufa.csv", sep = ",")

# вычисляем корреляцию признаков
cormat <- round(cor(data),2)
# изменяем форму и удлиняем корреляционную матрицу
melted_cormat <- melt(cormat)
# оставим в матрице только ее верхний треугольник
get_upper_tri <- function(cormat){
  cormat[lower.tri(cormat)]<- NA
  return(cormat)
}
upper_tri <- get_upper_tri(cormat)
melted_cormat <- melt(upper_tri, na.rm = TRUE)
# строим хитмэп
ggheatmap <- ggplot(data = melted_cormat, aes(Var2, Var1, fill = value))+
  geom_tile(color = "white")+
  scale_fill_gradient2(low = "blue", high = "red", mid = "white",
                       midpoint = 0, limit = c(-1,1), space = "Lab",
                       name="Pearson\nCorrelation") +
  theme_minimal()+
  theme(axis.text.x = element_text(angle = 45, vjust = 1,
                                   size = 12, hjust = 1))+
  coord_fixed()

# расположим значения на хитмэпе
ggheatmap +
  geom_text(aes(Var2, Var1, label = value), color = "black", size = 4) +
  theme(
    axis.title.x = element_blank(),
    axis.title.y = element_blank(),
    panel.grid.major = element_blank(),
    panel.border = element_blank(),
    panel.background = element_blank(),
    axis.ticks = element_blank(),
    legend.justification = c(1, 0),
    legend.position = c(0.6, 0.7),
    legend.direction = "horizontal")+
  guides(fill = guide_colorbar(barwidth = 7, barheight = 1,
                               title.position = "top", title.hjust = 0.5))

# исходя из матрицы корреляции необходимо избавиться от living area
drop <- c("living_area")
data = data[,!(names(data) %in% drop)]

# стандартизируем значения для корректной кластеризации
data_scaled <- scale(data)

# определим оптимальное количество кластеров методом локтя
fviz_nbclust(data_scaled, kmeans, method = "wss")
# определим оптимальное количество кластеров методом силуэта
# устанавливаем максимальное количество кластеров для тестирования
k.max <- 10

# используем функцию fviz_nbclust с методом силуэта
silhouette_nb <- fviz_nbclust(iris[, -5], kmeans, method = "silhouette", k.max = k.max)
# достаем оптимальное значение
n_clust <- silhouette_nb$optimal

# определяем диапазон значений nstart для проверки
nstart_range <- seq(1, 20, by = 1)

# создаем пустой вектор для хранения средней ширины силуэта
avg_silhouette <- vector(length = length(nstart_range))

# перебираем каждое значение nstart и вычисляем среднюю ширину силуэта
for (i in seq_along(nstart_range)) {
  set.seed(123) # установили начальное значение для воспроизводимости
  km <- kmeans(iris[, -5], centers = n_clust, nstart = nstart_range[i])
  avg_silhouette[i] <- mean(silhouette(km$cluster, dist(iris[, -5])))
}

# отображаем среднюю ширину силуэта для каждого значения nstart
plot(nstart_range, avg_silhouette, type = "b", xlab = "nstart", ylab = "Average silhouette width")

# находим оптимальное значение nstart методом силуэта
optimal_nstart <- nstart_range[which.max(avg_silhouette)]

# выведем оптимальное значение nstart
cat("Optimal nstart value:", optimal_nstart)

# проведем кластеризацию методом kmeans
km <- kmeans(data_scaled, n_clust, nstart = optimal_nstart)
# сохраним вектор с кластером в отдельной переменной
grp <- km$cluster

# отобразим кластеризацию с использованием метода главных компонентов
fviz_cluster(km, data_scaled,
             palette = "Set2", ggtheme = theme_minimal())

# отобразим диаграмму рассеяния признаков
pairs(data, col = brewer.pal(n = n_clust, name = "Set2")[grp], upper.panel=NULL)

# необходимо избавиться от незначительных признаков
drop <- c("floors_total", "ceiling_height", "rooms", "studio")
data_vs = data[,!(names(data) %in% drop)]

# отобразим диаграмму рассеяния признаков с учетом кластеризации
pairs(data_vs, col = brewer.pal(n = n_clust, name = "Set2")[grp], upper.panel=NULL)

# посмотрим на расположение центроидов кластеров
km$centers
km_centroids <- aggregate(data, by=list(cluster=km$cluster), mean)
km_centroids
# загрузим датасет с описанием центров кластеров
write.csv(km_centroids, "pr_km_centers.csv", row.names=FALSE)