K均值 ClusterR 编程示例

什么是 Cluster 分析？

Cluster 分析属于 无监督学习。聚类是一组具有相似特征的数据。我们可以说，聚类分析更多的是发现而不是预测。机器在数据中寻找相似性。例如，您可以将聚类分析用于以下应用程序：

• 客户细分：寻找客户群体之间的相似性
• 股票市场聚类：根据表现对股票进行分组
• 通过分组降低数据集的维度ping 具有相似值的观测结果

Cluster分析实施起来并不太困难，并且对于业务来说很有意义且可操作。

监督学习和无监督学习之间最显著的区别在于结果。无监督学习会创建一个新变量，即标签，而监督学习则会预测结果。机器可以帮助从业者根据密切相关性标记数据。分析师需要利用这些组并为其命名。

让我们举一个例子来理解聚类的概念。为简单起见，我们从两个维度进行分析。您有客户总支出及其年龄的数据。为了改善广告效果，营销团队希望向客户发送更有针对性的电子邮件。

在下图中，您可以绘制总支出和顾客年龄。

library(ggplot2)
df <- data.frame(age = c(18, 21, 22, 24, 26, 26, 27, 30, 31, 35, 39, 40, 41, 42, 44, 46, 47, 48, 49, 54),
    spend = c(10, 11, 22, 15, 12, 13, 14, 33, 39, 37, 44, 27, 29, 20, 28, 21, 30, 31, 23, 24)
)
ggplot(df, aes(x = age, y = spend)) +
    geom_point()

此时可见一种模式

1. 左下角是购买力较低的年轻人
2. 中上阶层是指有工作且能够负担得起更多支出的人
3. 最后，预算较低的老年人。

在上图中，您可以手动对观测值进行聚类，并定义三个组中的每一个。这个例子相当简单，而且非常直观。如果将新的观测值附加到数据集中，您可以在圆圈内标记它们。您可以根据我们的判断来定义圆圈。相反，您可以使用 机器学习 客观地对数据进行分组。

在本教程中，您将学习如何使用 k均值 算法。

K均值算法

K均值聚类算法是应用最广泛的聚类方法。该算法起源于20世纪50年代，此后经过多次改进。

该算法试图通过最小化观测值之间的距离来寻找组，这称为 局部最优 解决方案。距离是根据观测点的坐标来测量的。例如，在二维空间中，坐标就是 x 和 y。

该算法的工作原理如下：

• 步骤 1：在特征空间中随机选择 k 个初始质心
• 步骤 2：将每个观测值分配到其最近的聚类中心， 重心这将产生 k 个组
• 第三步： Shift 初始质心到组内坐标的平均值。
• 步骤 4：根据新的质心最小化距离。创建新的边界。因此，观测值将从一个组移动到另一个组
• 重复，直到没有观察改变组

K均值算法通常衡量两个观测值x和y之间的欧氏距离：

还有其他度量方法可用，例如曼哈顿距离或闵可夫斯基距离。请注意，每次运行 K-means 算法时，返回的分组都会不同。回想一下，初始猜测是随机的，算法会不断计算距离，直到组内达到同质性。也就是说，K-means 对初始选择非常敏感，除非观测值和分组数量都很少，否则几乎不可能得到相同的聚类结果。

选择聚类数量

k均值算法的另一个难点在于聚类数量的选择。你可以设置较大的k值（即较大的聚类数），以提高组内同质性，但这会带来风险。 过度拟合 数据过拟合意味着模型在新数据上的性能急剧下降，因为模型记住的是特定样本中的噪声，而不是潜在的模式。

聚类的数量取决于数据集的性质、行业、业务等。不过，选择合适的聚类数量有一个经验法则：

其中 n 是数据集中的观测值数量。

在实践中，花时间寻找最符合业务需求的 k 值是值得的。

我们将使用个人电脑价格数据集进行聚类分析。该数据集包含 6259 个观测值和 10 个特征。该数据集观察了 1993 年至 1995 年美国 486 台个人电脑的价格。变量包括价格、速度、内存、屏幕、CD 等。

您将按如下方式进行：

• 导入日期
• 训练模型
• 评估模型

导入日期

K-means算法不适用于因子变量，因为它基于距离计算，而离散标签之间没有实际意义的距离。请删除三个分类列（cd、multi、premium）以及行索引X。此数据集不包含缺失值。

library(dplyr)
PATH <-"https://raw.githubusercontent.com/guru99-edu/R-Programming/master/computers.csv"
df <- read.csv(PATH) %>%
select(-c(X, cd, multi, premium))
glimpse(df)

输出

## Observations: 6,259
## Variables: 7
## $ price <int> 1499, 1795, 1595, 1849, 3295, 3695, 1720, 1995, 2225, 2...
## $ speed <int> 25, 33, 25, 25, 33, 66, 25, 50, 50, 50, 33, 66, 50, 25, ...
## $ hd <int> 80, 85, 170, 170, 340, 340, 170, 85, 210, 210, 170, 210...
## $ ram <int> 4, 2, 4, 8, 16, 16, 4, 2, 8, 4, 8, 8, 4, 8, 8, 4, 2, 4, ...
## $ screen <int> 14, 14, 15, 14, 14, 14, 14, 14, 14, 15, 15, 14, 14, 14, ...
## $ ads <int> 94, 94, 94, 94, 94, 94, 94, 94, 94, 94, 94, 94, 94, 94, ...
## $ trend  <int> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1...

从汇总统计信息可以看出，数据值较大。在使用任何基于距离的方法之前，一个好的做法是对数据进行标准化，使每个变量都具有一个均值。 均值为零，标准差为一如果没有它，价格（其数值高达数千）将主导距离计算，而屏幕尺寸将几乎无关紧要。

summary(df)

输出：

##      price          speed              hd              ram        
##  Min.   : 949   Min.   : 25.00   Min.   :  80.0   Min.   : 2.000  
##  1st Qu.:1794   1st Qu.: 33.00   1st Qu.: 214.0   1st Qu.: 4.000
##  Median :2144   Median : 50.00   Median : 340.0   Median : 8.000  
##  Mean   :2220   Mean   : 52.01   Mean   : 416.6   Mean   : 8.287  
##  3rd Qu.:2595   3rd Qu.: 66.00   3rd Qu.: 528.0   3rd Qu.: 8.000  
##  Max.   :5399   Max.   :100.00   Max.   :2100.0   Max.   :32.000  
##      screen           ads            trend      
##  Min.   :14.00   Min.   : 39.0   Min.   : 1.00  
##  1st Qu.:14.00   1st Qu.:162.5   1st Qu.:10.00  
##  Median :14.00   Median :246.0   Median :16.00  
##  Mean   :14.61   Mean   :221.3   Mean   :15.93  
##  3rd Qu.:15.00   3rd Qu.:275.0   3rd Qu.:21.50  
##  Max.   :17.00   Max.   :339.0   Max.   :35.00

您可以使用 dplyr 库的 scale() 函数重新调整变量。转换可减少异常值的影响，并允许将单个观察值与平均值进行比较。如果标准化值（或 分数) 很高，那么您可以确信该观察结果确实高于平均值（较大的 z 分数意味着该点在标准差方面远离平均值。z 分数为 2 表示该值与平均值相差 XNUMX 个标准差。请注意，z 分数遵循高斯分布，并围绕平均值对称。

# Note: speed is deliberately left out of the scaled data frame
rescale_df <- df %>%
    mutate(price_scal = scale(price),
    hd_scal = scale(hd),
    ram_scal = scale(ram),
    screen_scal = scale(screen),
    ads_scal = scale(ads),
    trend_scal = scale(trend)) %>%
select(-c(price, speed, hd, ram, screen, ads, trend))

R base 有一个函数可以运行 k 均值算法，k 均值的基本函数是：

kmeans(df, k)
arguments:
-df: dataset used to run the algorithm
-k: Number of clusters

训练模型

在图三中，您详细展示了算法的工作原理。您可以使用谢义辉（Yihui Xie）编写的动画包以图形方式查看每个步骤，谢义辉也是 R Markdown 的 knitr 的开发者。该包不在 conda 频道中，因此请从 CRAN 安装：

install.packages("animation")

加载库后，在 kmeans 后添加 .ani，然后 R 将绘制所有步骤。为了便于说明，您仅使用具有三个集群的重新缩放变量 hd 和 ram 运行该算法。

set.seed(2345)
library(animation)
kmeans.ani(rescale_df[2:3], 3)

Code 说明

• kmeans.ani(rescale_df[2:3], 3)：选择 rescale_df 数据集的第 2 列和第 3 列，并运行 k 设置为 3 的算法。绘制动画。

您可以按如下方式解读该动画：

• 步骤 1：R 随机选择三个点
• 步骤 2：计算欧几里德距离并绘制聚类。左下方有一个绿色聚类，右侧有一个黑色大聚类，中间有一个红色聚类。
• 步骤 3：计算质心，即聚类的平均值
• 重复直到集群中没有数据改变

该算法在七次迭代后收敛。您可以在包含五个集群的数据集中运行 k-mean 算法，并将其称为 pc_cluster。

pc_cluster <-kmeans(rescale_df, 5)

列表 pc_cluster 包含七个有用元素：

• pc_cluster$cluster：分配给每个观测值的聚类
• pc_cluster$centers：集群中心
• pc_cluster$totss：总平方和
• pc_cluster$withinss：组内平方和，每个簇一个值
• pc_cluster$tot.withinss：内部总和
• pc_cluster$betweenss：总平方和减去组内平方和
• pc_cluster$size：每个集群内的观察数量

您将使用平方和（即 tot.withinss）的总和来计算最佳聚类数 k。寻找 k 确实是一项艰巨的任务。

如何使用肘部法则找到最优 k 值

选择最佳 k 的一种方法称为 肘法。该方法使用组内同质性或组内异质性来评估变异性。换句话说，您感兴趣的是每个聚类解释的方差百分比。您可以预期变异性会随着聚类数量的增加而增加，或者异质性会降低。我们的挑战是找到超出收益递减的 k。添加新聚类不会改善数据的变异性，因为剩下的信息很少。

在本教程中，我们使用异质性度量来找到这一点。簇内平方和总计是 kmean() 返回的列表中的 tot.withinss。

您可以构建肘部图并找到最佳 k，如下所示：

• 步骤 1：构建一个函数来计算聚类内平方和的总和
• 步骤 2：对一系列 k 值运行算法
• 步骤 3：使用算法结果创建数据框
• 步骤 4：绘制结果

步骤1） 构建一个函数来计算聚类内平方和的总和

您创建运行 k 均值算法的函数，并将总平方和存储在聚类内

kmean_withinss <- function(k) {
    cluster <- kmeans(rescale_df, k)
    return (cluster$tot.withinss)
}

Code 说明

• function(k)：设置函数中的参数数量
• kmeans(rescale_df, k): 针对此 k 值运行该算法
• return(clu​​ster$tot.withinss)：存储簇内平方和的总和

当 k 等于 2 时，测试该函数。

输出：

## Try with 2 cluster

kmean_withinss(2)

输出：

## [1] 27087.07

步骤2） 运行算法 n 次

您将使用 sapply() 函数在 k 范围内运行算法。此技术比创建循环并存储值更快。

# Set maximum cluster 
max_k <-20 
# Run algorithm over a range of k 
wss <- sapply(2:max_k, kmean_withinss)

Code 说明

• max_k <- 20：将 k 的最大值设置为 20。
• sapply(2：max_k，kmean_withinss)：在范围 2：max_k（即 2 到 20）内运行函数 kmean_withinss()。

步骤3） 使用算法的结果创建数据框

函数编写并测试完毕后，在 2 到 20 的范围内运行该函数，并存储每个 tot.withinss 值。

# Create a data frame to plot the graph
elbow <-data.frame(2:max_k, wss)

Code 说明

• data.frame(2:max_k, wss): 使用算法的输出创建一个数据框，存储在 wss 中

步骤4） 绘制结果

绘制图表来直观地了解肘点的位置

# Plot the graph with gglop
ggplot(elbow, aes(x = X2.max_k, y = wss)) +
    geom_point() +
    geom_line() +
    scale_x_continuous(breaks = seq(1, 20, by = 1))

从图中可以看出，最佳 k 为 7，其中曲线开始呈现收益递减的趋势。

一旦获得了最佳 k，即可重新运行 k 等于 7 的算法并评估聚类。

检查聚类

pc_cluster_2 <-kmeans(rescale_df, 7)

如前所述，您可以在 kmean() 返回的列表中访问其余有趣的信息。

pc_cluster_2$cluster
pc_cluster_2$centers
pc_cluster_2$size

评估是主观的，并且取决于聚类的用途。这里的目标是将具有相似规格的计算机分组。领域专家可以手动完成这项工作，但这会很慢且容易出错。K-means 算法可以实现分组。ping 客观地呈现结果，并由专家进行解释和命名。

作为事先评估，您可以检查集群的大小。

pc_cluster_2$size

输出：

## [1] 608 1596 1231  580 1003  699  542

最大的集群（编号 2）包含 1,596 个观测值，而最小的集群（编号 7）仅包含 542 台计算机。集群之间的同质性可能比较好，否则可能需要更精简的数据预处理。

通过中心组件，您可以更深入地了解数据。行代表聚类编号，列代表算法使用的变量。这些值是每个聚类在目标列上的平均得分。标准化使结果更易于解释。正值表示给定聚类的 z 分数高于总体均值。例如，聚类 4 的平均价格最高 (price_scal = 1.09)，而聚类 5 的平均价格最低 (-0.82)。

center <-pc_cluster_2$centers
center

输出：

##   price_scal    hd_scal     ram_scal screen_scal   ads_scal trend_scal
## 1 -0.6372457 -0.7097995 -0.691520682  -0.4401632  0.6780366 -0.3379751
## 2 -0.1323863  0.6299541  0.004786730   2.6419582 -0.8894946  1.2673184
## 3  0.8745816  0.2574164  0.513105797  -0.2003237  0.6734261 -0.3300536
## 4  1.0912296 -0.2401936  0.006526723   2.6419582  0.4704301 -0.4132057
## 5 -0.8155183  0.2814882 -0.307621003  -0.3205176 -0.9052979  1.2177279
## 6  0.8830191  2.1019454  2.168706085   0.4492922 -0.9035248  1.2069855
## 7  0.2215678 -0.7132577 -0.318050275  -0.3878782 -1.3206229 -1.5490909

您可以使用 ggplot 创建热图来帮助我们突出显示类别之间的差异。

ggplot 的默认颜色需要使用 RColorBrewer 库进行更改。您可以使用 conda 图书馆 以及在终端启动的代码：

conda install -c r r-rcolorbrewer

要创建热图，请按以下三个步骤进行：

• 使用中心值构建数据框，并使用聚类编号创建一个变量
• 使用 tidyr 库的 gather() 函数重塑数据。您需要将数据从宽转换为长。
• 用颜色创造调色板RampPalette() 函数

步骤1） 构建数据框

让我们创建重塑数据集

library(tidyr)

# create dataset with the cluster number

cluster <- c(1: 7)
center_df <- data.frame(cluster, center)

# Reshape the data

center_reshape <- gather(center_df, features, values, price_scal: trend_scal)
head(center_reshape)

输出：

##   cluster   features     values
## 1       1 price_scal -0.6372457
## 2       2 price_scal -0.1323863
## 3       3 price_scal  0.8745816
## 4       4 price_scal  1.0912296
## 5       5 price_scal -0.8155183
## 6       6 price_scal  0.8830191		

步骤2） 创建调色板

以下代码构建了热图使用的颜色调色板。

library(RColorBrewer)
# Create the palette
hm.palette <-colorRampPalette(rev(brewer.pal(10, 'RdYlGn')),space='Lab')

步骤 3）可视化

您可以绘制图表并查看聚类的样子。

# Plot the heat map
ggplot(data = center_reshape, aes(x = features, y = cluster, fill = values)) +
    scale_y_continuous(breaks = seq(1, 7, by = 1)) +
    geom_tile() +
    coord_equal() +
    scale_fill_gradientn(colours = hm.palette(90)) +
    theme_classic()

如何使用 set.seed() 和 nstart 获得可复现的 K-Means 结果

教程中提到，k-means算法每次运行都会返回不同的分组。这并非可以忽略的怪异现象，而是一个存在问题且有两种标准解决方法的缺陷，但上述示例中并未采用这两种方法。

1. 使用 set.seed() 确定起始点。 初始质心是随机抽取的，因此每次相同的调用都会产生不同的聚类结果。设置种子可以使整个分析过程可重复：

set.seed(123)
pc_cluster_2 <- kmeans(rescale_df, 7)

2. 运行多次启动，并使用 nstart 保留最佳启动结果。 单次随机起始点可能会收敛到较差的局部最优解。nstart 参数会从不同的随机起始点运行算法 n 次，并返回平方和最小的那个结果：

set.seed(123)
pc_cluster_2 <- kmeans(rescale_df, centers = 7, nstart = 25)

通常建议使用 25 个起始点，对于这种规模的数据集来说，这样做成本很低。如果没有 25 个起始点，肘部曲线本身就会变得不稳定，因为图上的每个点都来自一次运气好或运气不好的起始点。因此，前面定义的 kmean_withinss() 函数应该写成：

kmean_withinss <- function(k) {
    cluster <- kmeans(rescale_df, centers = k, nstart = 25)
    return (cluster$tot.withinss)
}

3. 如果算法发出警告，则增加 iter.max。 如果 R 报告“10 次迭代后仍未收敛”，请传递 iter.max = 50 以给予其更多空间。

K均值聚类与层次聚类 Cluster在 R 中

K-means 并不是唯一可用的聚类方法，选择哪种方法通常取决于数据的大小以及你是否已经知道你想要多少个组。

# Hierarchical alternative on the same scaled data
hc <- hclust(dist(rescale_df), method = "ward.D2")
plot(hc)
groups <- cutree(hc, k = 7)

由于该数据集包含 6,259 个观测值，层次聚类需要构建一个约 19.6 万个点对的距离矩阵，因此 k-means 算法是该数据集的实用选择。需要注意的是，k-means 算法假设聚类大小相近且形状近似球形；如果该假设不成立，DBSCAN 算法可以处理不规则形状并识别异常值，而不是强制将每个点都归入一个组。

K均值 ClusterR 语言中的函数参考

本教程中使用的所有函数如下所示：

常见问题