如何按site ggplot2上色(轮廓图)

Question

我有一个轮廓图，我想要根据站点而不是默认的集群编号进行着色。之前有人写了一些代码，我又开始尝试做plot How can I change the color to a variable other than cluster number in fviz_silhouette，不幸的是，他们的代码只会给我产生一个错误。我相信这是因为当把pamspec放到一个列表中时，所有的东西都变成了pamspec$data，而不是像site这样的单个组件。有人能帮上忙吗？

轮廓图

regforamcountsall<-read.csv("regionalforamscountsall.csv",header=TRUE)

library(ggplot2)
library(plyr)
library(cluster)
library(factoextra)
library(ggfortify)
library(dplyr)

spec<-regforamcountsall
str(spec)

#load in data and scale
spec
head(spec)

#visualise data and use clustering number from Av silhouette
pamspec<-pam(spec,4,keep.diss = TRUE)
plot(pamspec)
spec<-cbind(pamspec$clustering)

#create plot using cluster number

fviz_silhouette(pamspec, label=TRUE)+theme_classic()

#attempt to create new plot

myfviz_silhouette <- function (sil.obj, var.col, label = FALSE, print.summary = TRUE, ...) {
  if (inherits(sil.obj, c("eclust", "hcut", "pam", "clara", 
                          "fanny"))) {
    df <- as.data.frame(sil.obj$silinfo$widths, stringsAsFactors = TRUE)
  }
  else if (inherits(sil.obj, "silhouette")) 
    df <- as.data.frame(sil.obj[, 1:3], stringsAsFactors = TRUE)
  else stop("Don't support an oject of class ", class(sil.obj))
  df <- df[order(df$cluster, -df$sil_width), ]
  if (!is.null(rownames(df))) 
    df$name <- factor(rownames(df), levels = rownames(df))
  else df$name <- as.factor(1:nrow(df))
  df$cluster <- as.factor(df$cluster)
  df$var_col <- var.col
  mapping <- aes_string(x = "name", y = "sil_width", color = "var_col", 
                        fill = "var_col")
  p <- ggplot(df, mapping) + geom_bar(stat = "identity") + 
    labs(y = "Silhouette width Si", x = "", title = paste0("Clusters silhouette plot ", 
                                                           "\n Average silhouette width: ", round(mean(df$sil_width), 
                                                                                                  2))) + ggplot2::ylim(c(NA, 1)) + geom_hline(yintercept = mean(df$sil_width), 
                                                                                                                                              linetype = "dashed", color = "red")
  p <- ggpubr::ggpar(p, ...)
  if (!label) 
    p <- p + theme(axis.text.x = element_blank(), axis.ticks.x = element_blank())
  else if (label) 
    p <- p + theme(axis.text.x = element_text(angle = 45))
  ave <- tapply(df$sil_width, df$cluster, mean)
  n <- tapply(df$cluster, df$cluster, length)
  sil.sum <- data.frame(cluster = names(ave), size = n, ave.sil.width = round(ave, 
                                                                              2), stringsAsFactors = TRUE)
  if (print.summary) 
    print(sil.sum)
  p
}

library(factoextra)
library(cluster)
pamspec <- pam(pamspec, 4, keep.diss = TRUE) 

color_var <- pamspec$Site
myfviz_silhouette(pamspec, color_var, label=TRUE) +
  theme_classic()

#Error in myfviz_silhouette(pamspec, Site, label = TRUE) : 
#  object 'Site' not found

数据

structure(list(Sample = structure(c(2L, 23L, 30L, 24L, 25L, 28L
), .Label = c("LG1", "LG120", "LG130", "LG135", "LG160", "LG170", 
"LG185", "LG2", "LG225", "LG230", "LG240", "LG245", "LG255", 
"LG260", "LG275", "LG280", "LG285", "LG290", "LG295", "LG3", 
"LG305", "LG315", "LG32", "LG36", "LG38", "LG4", "LG48", "LG5", 
"LG60", "LG7", "LSP010", "LSP020", "LSP030", "LSP040", "LSP050", 
"LSP060", "LSP070", "LSP080", "LSP089", "LSP100", "LSP110", "LSP120", 
"LSP130", "LSP140", "LSP150", "LSP160", "LSP165", "ST-2LG0", 
"ST-2LG100", "ST-2LG120", "ST-2LG140", "ST-2LG160", "ST-2LG190", 
"ST-2LG40", "ST-2LG60", "ST-2LG80", "T3LB11.301", "T3LB12.05", 
"T3LB12.844", "T3LB13.87", "T3LB14.51", "T3LB14.63", "T3LB15.321", 
"T3LB15.59", "T3LB15.95", "T3LB16.69", "T3LB18.226", "T3LB19.762", 
"T3LB21.078", "T3LB26.256", "T3LB28.57", "T3LB28.84", "T3LB29.03", 
"T3LB31.056", "T3LB31.365", "T3LB7.008", "T3LB7.18", "T3LB7.303", 
"T3LB7.5", "T3LB7.9", "T3LB8.73", "T3LB9.45", "WAP 0 ST-2", "WAP 10 ST-2", 
"WAP 110 ST1", "WAP 120 ST-1", "WAP 122 ST-1", "WAP 125 ST1", 
"WAP 130 ST1", "WAP 135 ST-1", "WAP 140 ST-1", "WAP 144 ST-1", 
"WAP 150 ST-1 ", "WAP 155 ST-1", "WAP 159 ST1", "WAP 160 ST-1", 
"WAP 170 ST-1", "WAP 175 ST 1", "WAP 180 ST-1", "WAP 190 ST-1", 
"WAP 200 ST-1", "WAP 210 ST-1", "WAP 230 ST-1", "WAP 240 ST-1", 
"WAP 25 ST-2", "WAP 40 ST-2", "WAP 45 ST-2", "WAP 5  ST-2", "WAP 50 ST-2", 
"WAP 55 ST-2", "WAP 60 ST-1", "WAP 60 ST-2"), class = "factor"), 
    T.salsa = c(44.4, 75.4, 68, 51.9, 56, 0), H.wilberti = c(9.4, 
    1.8, 8.2, 10.3, 28.8, 0), Textularia = c(0.4, 0, 0, 0, 0, 
    0), T.irregularis = c(0, 3.6, 2.5, 0, 0.4, 0), P.ipohalina = c(0, 
    0, 0, 0, 0, 0), J.macrescens = c(5.4, 17.4, 7.4, 13.4, 6.6, 
    45.5), T.inflata = c(40.4, 1.8, 7.4, 24.4, 8.2, 53), S.lobata = c(0, 
    0, 0, 0, 0, 0), M.fusca = c(0, 0, 6.6, 0, 0, 1.5), A.agglutinans = c(0, 
    0, 0, 0, 0, 0), A.exiguus = c(0, 0, 0, 0, 0, 0), P.hyperhalina = c(0, 
    0, 0, 0, 0, 0), H.germanica = c(0, 0, 0, 0, 0, 0), Quinqueloculina = c(0, 
    0, 0, 0, 0, 0), Elphidium = c(0, 0, 0, 0, 0, 0), SWLI = c(190.4972376, 
    190.3867403, 190.2762431, 188.3977901, 186.8508287, 185.4143646
    ), Site = c(1, 1, 1, 1, 1, 1)), row.names = c(NA, 6L), class = "data.frame")

Answer 1

您遇到的问题在pamspec$Site中，因为它不存在

pamspec <- pam(regforamcountsall, 4, keep.diss = TRUE) 

color_var <- as.factor(pamspec$clustering)
myfviz_silhouette(pamspec, var.col = color_var, label=TRUE) +
  theme_classic()

​

​