k均值算法
算法介绍
给定样本集
D
=
{
x
1
,
x
2
,
.
.
.
,
x
m
}
D=\{x_1,x_2,...,x_m\}
D={x1,x2,...,xm},k均值算法针对聚类算法所得簇划分
C
=
{
C
1
,
C
2
,
.
.
.
,
C
k
}
\mathcal{C}=\{C_1,C_2,...,C_k\}
C={C1,C2,...,Ck}最小化平方误差,即:
E
=
∑
i
=
1
k
∑
x
∈
C
i
∣
∣
x
−
μ
i
∣
∣
2
2
(1)
E=\sum_{i=1}^k\sum_{x\in C_i}||x-\mu_i||_2^2\tag{1}
E=i=1∑kx∈Ci∑∣∣x−μi∣∣22(1)
其中 μ i = 1 ∣ C i ∣ ∑ x ∈ C i x \mu_i=\frac{1}{|C_i|}\sum_{x \in C_i}x μi=∣Ci∣1∑x∈Cix 是簇 C i C_i Ci的均值向量。由(1)式可知 E E E的值越小,其簇内样本的距离越小(相似度越高)。
它涉及到对所有可能的簇划分进行组合和比较。要找到全局最小值 E E E,需要尝试所有可能的组合,而随着数据量和簇数的增加,搜索空间呈指数增长。虽然可以使用启发式方法近似地求解这个问题,但找到确切的最优解对于大规模问题来说是计算上非常昂贵的,因此被认为是NP难问题。故k均值算法采用了贪心策略,通过迭代优化来近似的求解式(1)。k均值算法流程图如下图所示:
实验分析
数据集如下表所示:
读入数据集
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
data = pd.read_csv('data/4.0.csv')
k均值算法
# k值
k = 3
# 随机初始化簇中心
np.random.seed(0)
initial_centers_indices = np.random.choice(data.index, size=k, replace=False)
initial_centers = data.loc[initial_centers_indices]
# 迭代次数
max_iterations = 100
for iteration in range(max_iterations):
# 分配样本到簇
clusters = {i: [] for i in range(k)}
for index, row in data.iterrows():
# 计算该个样本点与其余所有均值向量的欧式距离
distances = [np.linalg.norm(row - center) for _, center in initial_centers.iterrows()]
assigned_cluster = np.argmin(distances) # 返回与某个均值向量最小的位置
clusters[assigned_cluster].append(index) # 将其划入到这个均值向量所对应的簇中
# 更新簇的均值向量
new_centers = pd.DataFrame([data.loc[cluster].mean() for cluster in clusters.values()])
# 检查簇中心变化是否小于某个阈值(收敛条件)
if np.allclose(initial_centers.values, new_centers.values, atol=1e-5):
break
initial_centers = new_centers
# 打印最终的簇划分
for i, (cluster_idx, cluster_data) in enumerate(clusters.items()):
print(f'Cluster {i + 1}:')
print(data.loc[cluster_data])
print('-' * 20)
输出结果:
Cluster 1:
Density Sugar inclusion rate
5 0.403 0.237
6 0.481 0.149
7 0.437 0.211
9 0.243 0.267
10 0.245 0.057
11 0.343 0.099
14 0.360 0.370
17 0.359 0.188
18 0.339 0.241
19 0.282 0.257
22 0.483 0.312
--------------------
Cluster 2:
Density Sugar inclusion rate
0 0.697 0.460
1 0.744 0.376
3 0.608 0.318
21 0.714 0.346
23 0.478 0.437
24 0.525 0.369
25 0.751 0.489
26 0.532 0.472
27 0.473 0.376
28 0.725 0.445
29 0.446 0.459
--------------------
Cluster 3:
Density Sugar inclusion rate
2 0.634 0.264
4 0.556 0.215
8 0.666 0.091
12 0.639 0.161
13 0.657 0.198
15 0.593 0.042
16 0.719 0.103
20 0.748 0.232
--------------------
绘制分类图像:
# 绘制原始数据点和簇中心,用不同颜色标记不同簇
for i, (cluster_idx, cluster_data) in enumerate(clusters.items()):
plt.scatter(data.loc[cluster_data]['Density'], data.loc[cluster_data]['Sugar inclusion rate'],
label=f'Cluster {i + 1}', alpha=0.7, s=50)
plt.scatter(initial_centers['Density'], initial_centers['Sugar inclusion rate'],
marker='X', s=200, label='Cluster Centers', c='black')
plt.title('K-Means Clustering with Decision Boundaries')
plt.xlabel('Density')
plt.ylabel('Sugar Inclusion Rate')
plt.legend()
plt.show()
学习向量量化(LVQ)
算法介绍
学习向量量化是一种用于模式分类和聚类的监督学习算法,试图通过找到一组原型向量来刻画聚类结构。
其算法流程图如下图所示:
若最近的原型向量
p
i
⋆
p_{i^\star}
pi⋆与
x
j
x_j
xj的类别标记相同 ,则令
p
i
⋆
p_{i^\star}
pi⋆向
x
j
x_j
xj的方向靠拢,如第7行所示:
p
′
=
p
i
⋆
+
η
⋅
(
x
j
−
p
i
⋆
)
(2)
p^\prime=p_{i^\star}+\eta\cdot (x_j-p_{i^\star}) \tag{2}
p′=pi⋆+η⋅(xj−pi⋆)(2)
故更新后的
p
′
p^\prime
p′和
x
j
x_j
xj之间的距离为:
∣
∣
p
′
−
x
j
∣
∣
2
=
∣
∣
p
i
⋆
+
η
⋅
(
x
j
−
p
i
⋆
)
−
x
j
∣
∣
2
=
(
1
−
η
)
⋅
∣
∣
p
i
⋆
−
x
j
∣
∣
2
(3)
\begin{aligned} ||p^\prime-x_j||_2&=||p_{i^\star}+\eta\cdot (x_j-p_{i^\star})-x_j||_2\\ &=(1-\eta)\cdot ||p_{i^\star}-x_j||_2 \end{aligned} \tag{3}
∣∣p′−xj∣∣2=∣∣pi⋆+η⋅(xj−pi⋆)−xj∣∣2=(1−η)⋅∣∣pi⋆−xj∣∣2(3)
令 η ∈ ( 0 , 1 ) \eta\in(0,1) η∈(0,1),则 1 − η < 1 1-\eta<1 1−η<1使得 ∣ ∣ p ′ − x j ∣ ∣ 2 < ∣ ∣ x j − p i ⋆ ∣ ∣ 2 ||p^\prime-x_j||_2<||x_j-p_{i^\star}||_2 ∣∣p′−xj∣∣2<∣∣xj−pi⋆∣∣2故原型向量 p i ⋆ p_{i^\star} pi⋆在更新为 p ′ p^\prime p′之后更接近 x j x_j xj。
最终学得一组原型向量
{
p
1
,
p
2
,
.
.
.
,
p
q
}
\{p_1,p_2,...,p_q\}
{p1,p2,...,pq}后即可实现对样本空间
χ
\chi
χ的簇划分,对于每个原型向量
p
i
p_i
pi定义了一个区域
R
i
R_i
Ri,在该区域中每个样本与
p
i
p_i
pi的距离不大于它与其他原型向量
p
i
′
p_{i^\prime}
pi′(
i
≠
i
′
i\ne i^\prime
i=i′)的距离,即:
R
i
=
{
x
∈
χ
∣
∣
∣
x
−
p
i
∣
∣
2
≤
∣
∣
x
−
p
i
′
∣
∣
2
,
i
≠
i
′
}
(4)
R_i=\{x\in \chi \mid ||x-p_i||_2\leq ||x-p_{i^\prime}||_2, i\ne i^\prime\} \tag{4}
Ri={x∈χ∣∣∣x−pi∣∣2≤∣∣x−pi′∣∣2,i=i′}(4)
由此形成了一个对样本空间 χ \chi χ的簇划分 { R 1 , R 2 , . . . , R q } \{R_1,R_2,...,R_q\} {R1,R2,...,Rq},该划分称为Voronoi剖分。
实验分析
数据集如下所示:
读入数据集:
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
df = pd.read_csv('data/4.0a.csv')
创建LQV模型:
# 定义学习率
learning_rate = 0.01
# 定义LVQ类
# 定义学习率
learning_rate = 0.01
# 定义LVQ类
class LVQ:
def __init__(self, input_size, num_prototypes):
self.prototypes = np.random.rand(num_prototypes, input_size)
self.prototype_labels = np.random.choice([0, 1], num_prototypes) # 随机设置预设的类别标记
def train(self, data, labels, epochs):
for epoch in range(epochs):
for i in range(len(data)):
current_data = data[i]
current_label = labels[i]
# 计算距离样本xj与原型向量pi的距离并返回距离最近的那个原型向量
winner_index = np.argmin(np.linalg.norm(self.prototypes - current_data, axis=1))
# 获取预设的类别标记
prototype_label = self.prototype_labels[winner_index]
if current_label == prototype_label: # 若标记相同,则使其更接近xj
self.prototypes[winner_index] += learning_rate * (current_data - self.prototypes[winner_index])
else: #若不同,则远离xj
self.prototypes[winner_index] -= learning_rate * (current_data - self.prototypes[winner_index])
提取数据集,并训练LVQ模型:
# 提取数据集
data = df[['Density', 'Sugar inclusion rate']].values
labels = df['label'].values
# 创建LVQ模型
lvq_model = LVQ(input_size=2, num_prototypes=5)
# 训练LVQ模型
lvq_model.train(data, labels, epochs=1000)
绘制原型向量坐标:
# 绘制数据点
plt.scatter(data[:, 0], data[:, 1], c=labels, cmap='viridis', label='Data Points')
# 绘制原型向量
plt.scatter(lvq_model.prototypes[:, 0], lvq_model.prototypes[:, 1], marker='X', s=100, c='red', label='Prototypes')
# 添加图例
plt.legend()
# 添加标签和标题
plt.xlabel('Density')
plt.ylabel('Sugar inclusion rate')
plt.title('LVQ Classification')
# 显示图像
plt.show()
高斯混合聚类
算法介绍
高斯混合聚类(Gaussian Mixture Model,简称GMM)是一种基于概率分布的聚类算法。它假设数据是由多个高斯分布组成的混合体,每个高斯分布对应一个聚类。GMM通过最大化似然函数来估计模型参数,包括每个高斯分布的均值、协方差矩阵和权重。
多元高斯分布概率密度函数为:
p
(
x
)
=
1
(
2
π
)
n
2
∣
Σ
∣
1
2
e
−
1
2
(
x
−
μ
)
T
Σ
−
1
(
x
−
μ
)
(5)
p(x)=\frac{1}{(2 \pi)^{\frac{n}{2}}|\Sigma|^{\frac{1}{2}}}e^{-\frac{1}{2}(x-\mu)^T\Sigma^{-1}(x-\mu)} \tag{5}
p(x)=(2π)2n∣Σ∣211e−21(x−μ)TΣ−1(x−μ)(5)
记为
x
∼
N
(
μ
,
Σ
)
x\sim\mathcal{N}(\mu,\Sigma)
x∼N(μ,Σ)。
其中
Σ
=
[
c
o
v
(
x
1
,
x
1
)
c
o
v
(
x
1
,
x
2
)
⋯
c
o
v
(
x
1
,
x
n
)
c
o
v
(
x
2
,
x
1
)
c
o
v
(
x
2
,
x
2
)
⋯
c
o
v
(
x
2
,
x
n
)
⋮
⋮
⋱
⋮
c
o
v
(
x
n
,
x
1
)
c
o
v
(
x
n
,
x
2
)
⋯
c
o
v
(
x
n
,
x
n
)
]
\Sigma=\begin{bmatrix} cov(x_1,x_1) & cov(x_1,x_2) & \cdots & cov(x_1,x_n) \\ cov(x_2,x_1) & cov(x_2,x_2) & \cdots & cov(x_2,x_n) \\ \vdots & \vdots & \ddots & \vdots \\ cov(x_n,x_1) & cov(x_n,x_2) & \cdots & cov(x_n,x_n) \\ \end{bmatrix}\
Σ=
cov(x1,x1)cov(x2,x1)⋮cov(xn,x1)cov(x1,x2)cov(x2,x2)⋮cov(xn,x2)⋯⋯⋱⋯cov(x1,xn)cov(x2,xn)⋮cov(xn,xn)
是
n
×
n
n\times n
n×n的协方差矩阵。
∣
Σ
∣
|\Sigma|
∣Σ∣是
Σ
\Sigma
Σ的行列式;
Σ
−
1
\Sigma^{-1}
Σ−1是
Σ
\Sigma
Σ的逆矩阵;
μ
\mu
μ是均值向量。
为了明确显示高斯分布的参数依赖关系,将概率密度函数记为 p ( x ∣ μ , Σ ) p(x\mid\mu,\Sigma) p(x∣μ,Σ)。
故我们可以定义高斯混合分布:
p
M
(
x
)
=
∑
i
=
1
k
α
i
⋅
(
x
∣
μ
i
,
Σ
i
)
p_{\mathcal{M}}(x)=\sum_{i=1}^k\alpha_i\cdotp(x\mid\mu_i,\Sigma_i)
pM(x)=i=1∑kαi⋅(x∣μi,Σi)
其中 α i > 0 \alpha_i>0 αi>0为权重(混合系数),使得 p M ( x ) p_{\mathcal{M}}(x) pM(x)符合概率密度函数的定义,即 ∫ − ∞ + ∞ p M ( x ) d x = 1 \int_{-\infty}^{+\infty}p_{\mathcal{M}}(x)dx=1 ∫−∞+∞pM(x)dx=1。该分布共由 k k k个混合成分组成,每个混合部分对应一个高斯分布。
若训练集
D
=
{
x
1
,
x
2
,
.
.
.
,
x
m
}
D=\{x_1,x_2,...,x_m\}
D={x1,x2,...,xm},令随机变量
z
j
∈
{
1
,
2
,
.
.
.
,
k
}
z_j\in\{1,2,...,k\}
zj∈{1,2,...,k}表示生成样本
x
j
x_j
xj的高斯混合成分(
z
j
=
i
z_j=i
zj=i表示生成样本
x
j
x_j
xj的高斯混合成分属于第
j
j
j个簇),故样本
x
j
x_j
xj的高斯混合成分属于第
j
j
j个簇的概率为:
p
M
(
z
j
=
i
∣
x
j
)
=
p
M
(
z
j
=
i
,
x
j
)
p
M
(
x
j
)
=
P
(
z
j
=
i
)
⋅
p
M
(
x
j
∣
z
j
=
i
)
p
M
(
x
j
)
=
α
i
⋅
p
(
x
j
∣
μ
i
,
Σ
i
)
∑
l
=
1
k
α
l
⋅
p
(
x
j
∣
μ
l
,
Σ
l
)
(6)
\begin{aligned} p_{\mathcal{M}}(z_j=i\mid x_j)&=\frac{p_{\mathcal{M}}(z_j=i, x_j)}{p_{\mathcal{M}}(x_j)}\\ &=\frac{P(z_j=i)\cdot p_{\mathcal{M}}(x_j\mid z_j=i)}{p_{\mathcal{M}}(x_j)}\\ &=\frac{\alpha_i\cdot p(x_j\mid\mu_i,\Sigma_i)}{\sum_{l=1}^k\alpha_l\cdot p(x_j\mid\mu_l,\Sigma_l)} \end{aligned} \tag{6}
pM(zj=i∣xj)=pM(xj)pM(zj=i,xj)=pM(xj)P(zj=i)⋅pM(xj∣zj=i)=∑l=1kαl⋅p(xj∣μl,Σl)αi⋅p(xj∣μi,Σi)(6)
我们将
p
M
(
z
j
=
i
∣
x
j
)
p_{\mathcal{M}}(z_j=i\mid x_j)
pM(zj=i∣xj)记为
γ
j
i
(
i
=
1
,
2
,
.
.
.
,
k
)
\gamma_{ji}(i=1,2,...,k)
γji(i=1,2,...,k)。
当确定了式(5)时,高斯混合聚类将样本集
D
D
D划分为
k
k
k个簇
C
=
{
C
1
,
C
2
,
.
.
.
,
C
k
}
\mathcal{C}=\{C_1,C_2,...,C_k\}
C={C1,C2,...,Ck},每个样本
x
j
x_j
xj的簇标记
λ
j
\lambda_j
λj如下确定:
λ
j
=
arg max
i
∈
{
1
,
2
,
.
.
.
,
k
}
γ
j
i
(7)
\lambda_j=\underset{i \in\{1,2,...,k\}}{\text{arg max}} \ \gamma_{ji} \tag{7}
λj=i∈{1,2,...,k}arg max γji(7)
即找到样本 x j x_j xj概率最大的那个标签。
那么对于式(5),模型参数
{
(
α
i
,
μ
i
,
Σ
i
)
∣
,
1
≤
i
≤
k
}
\{(\alpha_i,\mu_i,\Sigma_i)\mid ,1\leq i\leq k\}
{(αi,μi,Σi)∣,1≤i≤k}采用极大似然估计法,最大化对数似然函数:
L
L
(
D
)
=
ln
(
∑
j
=
1
m
p
M
(
x
j
)
)
=
∑
j
=
1
m
ln
(
∑
i
=
1
k
α
i
⋅
(
x
∣
μ
i
,
Σ
i
)
)
(8)
\begin{aligned} LL(D)&=\ln(\sum_{j=1}^mp_{\mathcal{M}}(x_j))\\ &=\sum_{j=1}^m\ln(\sum_{i=1}^k\alpha_i\cdotp(x\mid\mu_i,\Sigma_i)) \end{aligned} \tag{8}
LL(D)=ln(j=1∑mpM(xj))=j=1∑mln(i=1∑kαi⋅(x∣μi,Σi))(8)
令
∂
L
L
(
D
)
∂
μ
i
=
0
\frac{\partial LL(D)}{\partial \mu_i}=0
∂μi∂LL(D)=0有:
μ
i
=
∑
j
=
1
m
γ
j
i
x
j
∑
j
=
1
m
γ
j
i
(9)
\mu_i=\frac{\sum_{j=1}^m\gamma_{ji}x_j}{\sum_{j=1}^m \gamma_{ji}} \tag{9}
μi=∑j=1mγji∑j=1mγjixj(9)
令
∂
L
L
(
D
)
∂
Σ
i
=
0
\frac{\partial LL(D)}{\partial \Sigma_i}=0
∂Σi∂LL(D)=0有:
Σ
i
=
∑
j
=
1
m
γ
j
i
(
x
j
−
μ
i
)
(
x
j
−
μ
i
)
T
∑
j
=
1
m
γ
j
i
(10)
\Sigma_i=\frac{\sum_{j=1}^m\gamma_{ji}(x_j-\mu_i)(x_j-\mu_i)^T}{\sum_{j=1}^m\gamma_{ji}} \tag{10}
Σi=∑j=1mγji∑j=1mγji(xj−μi)(xj−μi)T(10)
对于混合系数
α
i
\alpha_i
αi,不仅要最大化
L
L
(
D
)
LL(D)
LL(D),还要满足
∑
i
=
1
k
α
i
=
1
\sum_{i=1}^k\alpha_i=1
∑i=1kαi=1,考虑拉格朗日乘子法有:
L
L
(
D
)
+
λ
(
∑
i
=
1
k
α
i
−
1
)
(11)
LL(D)+\lambda(\sum_{i=1}^k\alpha_i-1)\tag{11}
LL(D)+λ(i=1∑kαi−1)(11)
令 ∂ L L ( D ) ∂ α i = 0 \frac{\partial LL(D)}{\partial \alpha_i}=0 ∂αi∂LL(D)=0和 ∂ L L ( D ) ∂ λ = 0 \frac{\partial LL(D)}{\partial \lambda}=0 ∂λ∂LL(D)=0有:
{ ∑ j = 1 m p ( x j ∣ μ i , Σ i ) ∑ l = 1 k α l ⋅ p ( x j ∣ μ l , Σ l ) + λ = 0 ∑ i = 1 k α i = 1 (12) \begin{cases} \sum_{j=1}^m \frac{p\left(\boldsymbol{x}_j \mid \boldsymbol{\mu}_i, \boldsymbol{\Sigma}_i\right)}{\sum_{l=1}^k \alpha_l \cdot p\left(\boldsymbol{x}_j \mid \boldsymbol{\mu}_l, \boldsymbol{\Sigma}_l\right)}+\lambda=0\\ \sum_{i=1}^k\alpha_i=1 \end{cases} \tag{12} {∑j=1m∑l=1kαl⋅p(xj∣μl,Σl)p(xj∣μi,Σi)+λ=0∑i=1kαi=1(12)
解出:
{
λ
=
−
m
α
i
=
1
m
∑
j
=
1
m
γ
j
i
(13)
\begin{cases} \lambda=-m\\ \alpha_i=\frac{1}{m}\sum_{j=1}^m\gamma_{ji} \end{cases}\tag{13}
{λ=−mαi=m1∑j=1mγji(13)
高斯混合聚类算法流程如下图所示:
实验分析
数据集如下表所示:
读入数据集:
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
data = pd.read_csv('data/4.0.csv')
定义多元高斯分布函数:
# 选择需要聚类的特征
X = data[['Density', 'Sugar inclusion rate']].values
# 定义高斯分布函数
def gaussian(x, mean, cov):
exponent = -0.5 * np.dot(np.dot((x - mean).T, np.linalg.inv(cov)), (x - mean))
return np.exp(exponent) / (2 * np.pi * np.sqrt(np.linalg.det(cov)))
初始化参数:
# 初始化参数
num_clusters = 3
max_iterations = 500
tolerance = 1e-4
# 初始化均值、协方差矩阵和权重
means = np.random.rand(num_clusters, X.shape[1])
covariances = [np.identity(X.shape[1]) for _ in range(num_clusters)]
weights = np.ones(num_clusters) / num_clusters
EM算法
# EM算法
for iteration in range(max_iterations):
# E 步骤
responsibilities = np.zeros((X.shape[0], num_clusters))
for i in range(X.shape[0]):
for j in range(num_clusters):
responsibilities[i, j] = weights[j] * gaussian(X[i], means[j], covariances[j])
responsibilities[i, :] /= np.sum(responsibilities[i, :])
# M 步骤
Nk = np.sum(responsibilities, axis=0)
weights = Nk / X.shape[0]
means = np.dot(responsibilities.T, X) / Nk[:, None]
for j in range(num_clusters):
# 计算新的协方差矩阵
covariances[j] = np.dot((responsibilities[:, j] * (X - means[j]).T), (X - means[j])) / Nk[j]
画出可视化结果:
# 根据聚类结果可视化数据
labels = np.argmax(responsibilities, axis=1)
data['Cluster'] = labels
plt.scatter(data['Density'], data['Sugar inclusion rate'], c=data['Cluster'], cmap='viridis')
plt.xlabel('Density')
plt.ylabel('Sugar inclusion rate')
plt.title('Gaussian Mixture Clustering (Manual Implementation)')
plt.show()
总结
特性 | k-means | LVQ | 高斯混合聚类 |
---|---|---|---|
优点 | - 简单、易于理解和实现 - 计算效率较高 | - 具有在线学习能力 - 可以保留原始数据的拓扑结构 - 可以处理动态数据 | - 对于复杂的数据分布有较好的拟合能力 - 软聚类:对每个数据点都给出其属于每个簇的概率 - 每个簇的形状由协方差矩阵决定,因此能够适应各种形状的簇 |
缺点 | - 对初始聚类中心敏感 - 对噪声和异常值敏感 - 不能处理非球形簇 | - 对初始权重和原型向量的选择敏感 - 需要调整学习率和邻域大小 - 需要事先确定原型向量的数量和类别 | - 计算复杂度较高 - 需要选择合适的簇数目 - 对于高维数据,样本数量较少时可能过拟合 - 初始参数的选择对结果影响较大 |