Скрипт на R для работы

 # ввожу данные из файла
table_1 <- read.csv(file = (description = "i:/10pl.csv"),
                    sep = ";", h=F, dec = ",")
table_1[is.na(table_1 == "NA")] <- 0

# создаю заготовки матриц
v_I_J <- NULL # для индекса Жаккара
v_I_J <- as.data.frame(v_I_J[length(table_1), length(table_1)])
v_K_R <- NULL # для меры Ружички
v_K_R <- as.data.frame(v_K_R[length(table_1), length(table_1)])

# начинаю цикл
for(i in 1:length(table_1)) {
  wp_1 <- table_1[, i] # считываю первую выборку
  for(j in 1:length(table_1)) {
    wp_2 <- table_1[, j] # считываю вторую выборку
    
    # задаю счётные переменные
    a <- 0 # число видов, общих для двух списков
    b <- 0 # число видов, присутствующих только во втором списке
    c <- 0 # число видов, присутствующих только в первом списке
    d <- 0 # число видов, отсутствующих в обоих списках, но присутствующих в системе
    
    v_min <- 0 # сумма минимальных обилий из каждой пары
    v_max <- 0 # сумма максимальных обилий из каждой пары
    
    mn <- 0 # сумма произведений обилий
    
    sq_1 <- 0 # сумма квадратов обилий на первой площадке
    sq_2 <- 0 # сумма квадратов обилий на второй площадке
    
    # выполняю подсчёт основных параметров
    for(k in 1:length(wp_1)) {
      
      # считаю числа видов
      ifelse(wp_1[k] != 0 & wp_2[k] != 0, a <- a + 1, a <- a)
      ifelse(wp_1[k] != 0 & wp_2[k] == 0, b <- b + 1, b <- b)
      ifelse(wp_1[k] == 0 & wp_2[k] != 0, c <- c + 1, c <- c)
      ifelse(wp_1[k] == 0 & wp_2[k] == 0, d <- d + 1, d <- d)
      
      # считаю суммарные обилия
      v_min <- v_min + min(wp_1[k], wp_2[k])
      v_max <- v_max + max(wp_1[k], wp_2[k])
      
      # считаю сумму произведений обилий
      mn <- mn + wp_1[k] * wp_2[k]
      
      # считаю сумму квадратов обилий
      sq_1 <- sq_1 + wp_1[k]^2
      sq_2 <- sq_2 + wp_2[k]^2
    }
    
    # подсчитываю число видов в каждой выборке
    S_1 <- length(wp_1) - length(which(wp_1 == 0))
    S_2 <- length(wp_2) - length(which(wp_2 == 0))
    
    # подсчитываю суммарное обилие по каждой выборке
    v_S_1 <- sum(wp_1)
    v_S_2 <- sum(wp_2)
    
    # рассчитываю индексы и заполняю матрицы
    I_J <- a / (a + b + c) # считаю индекс Жаккара
    v_I_J[i, j] <- I_J # заношу его в соответствующую ячейку
    K_R <- v_min / v_max # считаю меру Ружички
    v_K_R[i, j] <- K_R # заношу её в соответствующую ячейку
  }
}

# рассчитываю матрицу расстояний через евклидово расстояние
table_J <- dist(scale(v_I_J), method = "euclidean") # для индекса Жаккара
table_R <- dist(scale(v_K_R), method = "euclidean") # для меры Ружички

# провожу кластеризацию для обоих индексов
hc_J1 <- hclust(table_J) # провожу кластеризацию методом полной связи
hc_R1 <- hclust(table_R) # провожу кластеризацию методом полной связи
par(mfrow=c(1, 2)) # делю лист на две части по горизонтали
plot(hc_J1, cex = 0.9, hang = -1) # строю дендрограмму для Жаккара
plot(hc_R1, cex = 0.9, hang = -1) # строю дендрограмму для Ружички

# проверяю качество кластеризации
library(pvclust) # запускаю нужную библиотеку
par(mfrow=c(1, 2)) # делю лист на две колонки

# для индекса Жаккара
# проверяю способы нахождения расстояния в матрице
table_pv_1J <- pvclust(t(v_I_J), method.dist = "euclidean",
                      method.hclust = "complete", nboot = 1000)
plot(table_pv_1J, print.num = FALSE)
table_pv_2J <- pvclust(t(v_I_J), method.dist = "manhattan",
                      method.hclust = "complete", nboot = 1000)
plot(table_pv_2J, print.num = FALSE)

# проверяю способы нахождения расстояния между кластерами
table_pv_3J <- pvclust(t(v_I_J), method.dist = "euclidean",
                      method.hclust = "complete", nboot = 1000)
plot(table_pv_3J, print.num = FALSE)
table_pv_4J <- pvclust(t(v_I_J), method.dist = "euclidean",
                      method.hclust = "average", nboot = 1000)
plot(table_pv_4J, print.num = FALSE)

# для меры Ружечки
# проверяю способы нахождения расстояния в матрице
table_pv_1R <- pvclust(t(v_K_R), method.dist = "euclidean",
                       method.hclust = "complete", nboot = 1000)
plot(table_pv_1R, print.num = FALSE)
table_pv_2R <- pvclust(t(v_K_R), method.dist = "manhattan",
                       method.hclust = "complete", nboot = 1000)
plot(table_pv_2R, print.num = FALSE)

# проверяю способы нахождения расстояния между кластерами
table_pv_3R <- pvclust(t(v_K_R), method.dist = "manhattan",
                       method.hclust = "complete", nboot = 1000)
plot(table_pv_3R, print.num = FALSE)
table_pv_4R <- pvclust(t(v_K_R), method.dist = "manhattan",
                       method.hclust = "average", nboot = 1000)
plot(table_pv_4R, print.num = FALSE)

# другой способ построения дендрограммы
par(mfrow=c(1, 1)) # возвращаю лист к целому

# для индекса Жаккара
hc_J2 <- as.dendrogram(hclust(table_J))
nodePar <- list(lab.cex = 0.7, pch = c(NA, 19), cex = 0.9, col = "blue")
plot(hc_J2, xlab = "Height", nodePar = nodePar, horiz = TRUE,
     edgePar = list(col = 2:3, lwd = 2 : 1))

# для меры Ружечки
hc_R2 <- as.dendrogram(hclust(table_R))
nodePar <- list(lab.cex = 0.7, pch = c(NA, 19), cex = 0.9, col = "blue")
plot(hc_R2, xlab = "Height", nodePar = nodePar, horiz = TRUE,
     edgePar = list(col = 2:3, lwd = 2 : 1))

# построение клубочковой дендрограммы
library(dendextend) # подключаю нужную библиотеку

m_J <- dist(scale(v_I_J), method = "euclidean") # считаю евклидово расстояние для индекса Жаккара
m_R <- dist(scale(v_K_R), method = "manhattan") # считаю манхеттенское расстояние для меры Ружички

hc_J3 <- hclust(m_J, method = "average") # провожу кластеризацию методом средней связи для индекса Жаккара
hc_R3 <- hclust(m_R, method = "complete") # провожу кластеризацию методом полной связи для меры Ружички

par(mfrow=c(1, 2)) # делю лист на две колонки
plot(hc_J3, cex = 0.9, hang = -1) # строю дендрограмму
plot(hc_R3, cex = 0.9, hang = -1) # строю дендрограмму

# на диаграмме типа «клубок» показываю цветом одинаковые ветви
tanglegram(hc_J3, hc_R3, common_subtrees_color_branches = TRUE)

# провожу неиерархическую кластеризацию
library(cluster) # подключаю нужную библиотеку

# для индекса Жаккара
non_ie_J <- fanny(v_I_J, 3) # провожу кластеризацию на 3 кластера
plot(non_ie_J, which=1) # строю диаграмму
head(data.frame(cl=c(1:10), non_ie_J$membership)) # вывожу данные

# для меры Ружички
non_ie_R <- fanny(v_K_R, 2) # провожу кластеризацию на 2 кластера
plot(non_ie_R, which=1) # строю диаграмму
head(data.frame(cl=c(1:10), non_ie_R$membership)) # вывожу данные

Последнее изменение: Понедельник, 6 декабря 2021, 08:57