AP聚类算法实现三维数据点分类
基于亲和传播(Affinity Propagation, AP)聚类算法的三维数据点分类MATLAB实现。AP聚类是一种无监督学习算法，无需预先指定聚类数量，通过数据点之间的"消息传递"自动确定聚类中心。
function ap_clustering_3d()
% 生成三维测试数据
[data, true_labels] = generate_3d_data();
% 计算相似度矩阵
S = compute_similarity(data);
% 执行AP聚类
[idx, centers, exemplars] = apcluster(S);
% 可视化聚类结果
visualize_clusters(data, idx, centers, true_labels);
% 评估聚类性能
evaluate_clustering(data, idx, true_labels, exemplars);
end
function [data, true_labels] = generate_3d_data()
% 生成三维测试数据
rng(42); % 设置随机种子保证可重复性
% 定义4个高斯分布的参数
num_points = 300; % 总点数
num_clusters = 4; % 真实聚类数量
% 聚类中心
centers = [1 1 1;
-1 -1 1;
1 -1 -1;
-1 1 -1];
% 生成数据点
data = [];
true_labels = [];
for i = 1:num_clusters
cluster_points = centers(i,:) + 0.3 * randn(num_points/num_clusters, 3);
data = [data; cluster_points];
true_labels = [true_labels; i * ones(num_points/num_clusters, 1)];
end
% 添加一些噪声点
noise_points = 4 * (rand(10, 3) - 2);
data = [data; noise_points];
true_labels = [true_labels; zeros(size(noise_points, 1), 1)];
% 绘制原始数据
figure('Name', '原始三维数据', 'Position', [100, 100, 800, 600]);
scatter3(data(:,1), data(:,2), data(:,3), 30, true_labels, 'filled');
title('原始三维数据（带真实标签）');
xlabel('X'); ylabel('Y'); zlabel('Z');
colormap(jet(num_clusters+1));
colorbar('Ticks', 0:num_clusters, 'TickLabels', ['噪声'; cellstr(num2str((1:num_clusters)'))]);
grid on; axis equal;
view(30, 30); % 设置视角
end
function S = compute_similarity(data)
% 计算相似度矩阵（负的欧氏距离平方）
n = size(data, 1);
S = zeros(n, n);
for i = 1:n
for j = 1:n
if i ~= j
% 计算欧氏距离平方
dist_sq = sum((data(i,:) - data(j,:)).^2);
% 相似度为负的欧氏距离平方
S(i,j) = -dist_sq;
else
% 对角线设置为偏好值（preference）
% 使用中位数作为偏好值
S(i,j) = median(S(:));
end
end
end
% 使用中位数作为偏好值
median_val = median(S(S < 0)); % 只考虑负值部分
S(1:n+1:end) = median_val;
end
function [idx, centers, exemplars] = apcluster(S)
% AP聚类算法实现
% 输入: S - 相似度矩阵
% 输出:
% idx - 每个点的聚类标签
% centers - 聚类中心索引
% exemplars - 每个点的代表点索引
% 算法参数
max_iter = 1000; % 最大迭代次数
conv_iter = 50; % 收敛检查迭代次数
lambda = 0.9; % 阻尼系数
n = size(S, 1);
% 初始化信息矩阵
R = zeros(n, n); % 责任矩阵
A = zeros(n, n); % 可用性矩阵
% 存储收敛状态
conv_history = zeros(1, max_iter);
% 迭代更新
for iter = 1:max_iter
% 保存旧的R和A用于收敛检查
R_old = R;
A_old = A;
% 更新责任矩阵 (Responsibility)
AS = A + S;
[Y, I] = max(AS, [], 2);
for i = 1:n
AS(i, I(i)) = -Inf;
end
Y2 = max(AS, [], 2);
R_new = S - repmat(Y, 1, n);
for i = 1:n
R_new(i, I(i)) = S(i, I(i)) - Y2(i);
end
% 应用阻尼
R = (1 - lambda) * R_new + lambda * R_old;
% 更新可用性矩阵 (Availability)
Rp = max(R, 0);
Rp(1:n+1:end) = R(1:n+1:end); % 对角线特殊处理
A_new = repmat(sum(Rp, 1), n, 1) - Rp;
dA = diag(A_new);
A_new = min(A_new, 0);
A_new(1:n+1:end) = dA;
% 应用阻尼
A = (1 - lambda) * A_new + lambda * A_old;
% 计算收敛指标
conv_history(iter) = mean(abs(R(:) - R_old(:))) + mean(abs(A(:) - A_old(:)));
% 检查收敛（连续conv_iter次变化小于阈值）
if iter > conv_iter
recent_changes = conv_history(iter-conv_iter+1:iter);
if all(recent_changes < 1e-6)
fprintf('AP聚类在 %d 次迭代后收敛\n', iter);
break;
end
end
end
if iter == max_iter
fprintf('达到最大迭代次数 %d\n', max_iter);
end
% 确定聚类中心和分配点
E = A + R; % 判断矩阵
exemplars = diag(E) > 0; % 代表点（聚类中心）
centers = find(exemplars); % 聚类中心索引
% 分配每个点到其聚类中心
[~, idx] = max(S(:, centers), [], 2);
% 处理未分配点（理论上AP会分配所有点）
unassigned = find(isnan(idx));
if ~isempty(unassigned)
fprintf('警告: %d 个点未被分配\n', length(unassigned));
idx(unassigned) = length(centers) + 1; % 分配到新聚类
end
fprintf('发现 %d 个聚类中心\n', length(centers));
end
function visualize_clusters(data, idx, centers, true_labels)
% 可视化聚类结果
figure('Name', 'AP聚类结果', 'Position', [100, 100, 1200, 500]);
% 子图1：聚类结果
subplot(1, 2, 1);
cluster_centers = data(centers, :);
unique_labels = unique(idx);
num_clusters = length(unique_labels);
% 为每个聚类分配颜色
colors = lines(num_clusters);
hold on;
for i = 1:num_clusters
cluster_points = data(idx == unique_labels(i), :);
scatter3(cluster_points(:,1), cluster_points(:,2), cluster_points(:,3), ...
30, colors(i,:), 'filled');
end
% 绘制聚类中心
scatter3(cluster_centers(:,1), cluster_centers(:,2), cluster_centers(:,3), ...
150, 'k', 'p', 'LineWidth', 2);
title(sprintf('AP聚类结果 (%d个聚类)', num_clusters));
xlabel('X'); ylabel('Y'); zlabel('Z');
grid on; axis equal;
view(30, 30);
% 添加图例
legend_str = arrayfun(@(x) sprintf('聚类 %d', x), 1:num_clusters, 'UniformOutput', false);
legend([legend_str, {'聚类中心'}], 'Location', 'best');
% 子图2：与真实标签对比
subplot(1, 2, 2);
scatter3(data(:,1), data(:,2), data(:,3), 30, idx, 'filled');
title('聚类分配');
xlabel('X'); ylabel('Y'); zlabel('Z');
colormap(lines(num_clusters));
colorbar('Ticks', 1:num_clusters, 'TickLabels', arrayfun(@num2str, 1:num_clusters, 'UniformOutput', false));
grid on; axis equal;
view(30, 30);
% 创建单独的聚类中心图
figure('Name', '聚类中心与数据分布', 'Position', [100, 100, 800, 600]);
hold on;
% 绘制所有数据点（灰色）
scatter3(data(:,1), data(:,2), data(:,3), 20, [0.7, 0.7, 0.7], 'filled');
% 高亮聚类中心
scatter3(cluster_centers(:,1), cluster_centers(:,2), cluster_centers(:,3), ...
200, 'r', 'p', 'LineWidth', 2);
% 连接聚类中心
for i = 1:size(cluster_centers, 1)
for j = i+1:size(cluster_centers, 1)
line([cluster_centers(i,1), cluster_centers(j,1)], ...
[cluster_centers(i,2), cluster_centers(j,2)], ...
[cluster_centers(i,3), cluster_centers(j,3)], ...
'Color', 'k', 'LineStyle', '--', 'LineWidth', 1);
end
end
title('聚类中心与数据分布');
xlabel('X'); ylabel('Y'); zlabel('Z');
legend('数据点', '聚类中心', 'Location', 'best');
grid on; axis equal;
view(30, 30);
end
function evaluate_clustering(data, idx, true_labels, exemplars)
% 评估聚类性能
% 计算轮廓系数
s = silhouette(data, idx);
avg_silhouette = mean(s);
% 计算Calinski-Harabasz指数
eval = evalclusters(data, idx, 'CalinskiHarabasz');
% 计算戴维斯-布尔丁指数
db = evalclusters(data, idx, 'DaviesBouldin');
% 显示评估指标
fprintf('\n聚类性能评估:\n');
fprintf('平均轮廓系数: %.4f\n', avg_silhouette);
fprintf('Calinski-Harabasz指数: %.4f\n', eval.CriterionValues);
fprintf('戴维斯-布尔丁指数: %.4f\n', db.CriterionValues);
% 显示聚类中心信息
centers = find(exemplars);
cluster_centers = data(centers, :);
fprintf('\n聚类中心坐标:\n');
disp(cluster_centers);
% 绘制轮廓系数图
figure('Name', '聚类质量评估', 'Position', [100, 100, 800, 400]);
subplot(1, 2, 1);
silhouette(data, idx);
title('轮廓系数');
% 绘制聚类大小分布
subplot(1, 2, 2);
cluster_counts = histcounts(idx, 'BinMethod', 'integers');
bar(cluster_counts);
title('聚类大小分布');
xlabel('聚类标签');
ylabel('点数');
grid on;
% 计算并显示每个聚类的统计信息
unique_labels = unique(idx);
num_clusters = length(unique_labels);
fprintf('\n聚类统计信息:\n');
fprintf('%-10s %-10s %-10s %-10s\n', '聚类ID', '点数', '平均轮廓', '中心坐标');
for i = 1:num_clusters
cluster_points = data(idx == unique_labels(i), :);
cluster_size = size(cluster_points, 1);
cluster_silhouette = mean(s(idx == unique_labels(i)));
if i <= size(cluster_centers, 1)
center_coord = sprintf('(%.2f, %.2f, %.2f)', ...
cluster_centers(i,1), ...
cluster_centers(i,2), ...
cluster_centers(i,3));
else
center_coord = 'N/A';
end
fprintf('%-10d %-10d %-10.4f %-10s\n', ...
unique_labels(i), cluster_size, cluster_silhouette, center_coord);
end
end
% 运行主函数
ap_clustering_3d();
AP聚类算法详解
1. 算法核心思想
亲和传播(AP)聚类是一种基于"消息传递"的聚类算法，其主要特点包括：
- 无需预先指定聚类数量：算法自动确定最佳聚类数
- 以实际数据点作为聚类中心：聚类中心称为"exemplars"
- 处理非凸形状聚类：适用于各种形状的数据分布
- 对初始点不敏感：结果稳定，不受随机初始化的影响
2. 核心数学原理
AP聚类通过两种消息传递进行迭代：
责任(Responsibility)消息 R(i,k)：
- 表示点i对点k作为其代表的累积证据
- 计算公式：$R(i,k) = S(i,k) - \max_{k' \neq k} {A(i,k') + S(i,k')}$
可用性(Availability)消息 A(i,k)：
- 表示点k适合作为点i代表的累积证据
- 计算公式：$A(i,k) = \min \left(0, R(k,k) + \sum_{i' \notin {i,k}} \max(0, R(i',k)) \right)$
其中S(i,k)是相似度矩阵，通常定义为负的欧氏距离平方：$S(i,k) = -|x_i - x_k|^2$
参考代码 ap聚类算法实现三维数据点的分类 youwenfan.com/contentalb/78315.html
3. 算法实现关键步骤
数据生成：
- 创建4个三维高斯分布作为基础聚类
- 添加随机噪声点模拟真实场景
相似度矩阵计算：
- 使用负的欧氏距离平方作为相似度度量
- 对角线设置为偏好值(preference)，控制聚类数量
消息传递迭代：
- 交替更新责任矩阵R和可用性矩阵A
- 应用阻尼因子(λ=0.9)确保收敛
- 连续50次迭代变化小于阈值则判定收敛
聚类确定：
- 结合A和R矩阵确定聚类中心(exemplars)
- 将每个点分配到相似度最高的聚类中心
结果可视化：
- 三维散点图展示聚类结果
- 突出显示聚类中心
- 对比真实标签与聚类结果
性能评估：
- 轮廓系数：衡量聚类内聚度和分离度
- Calinski-Harabasz指数：类间方差与类内方差之比
- Davies-Bouldin指数：类内距离与类间距离之比
4. 参数调整建议
偏好值(preference)：
- 控制聚类数量，值越大聚类越多
- 默认使用相似度矩阵的中位数
- 可尝试：
median(S)、min(S)或max(S)
阻尼因子(lambda)：
- 范围0.5-1.0，防止振荡
- 默认0.9，收敛慢时可增大到0.95
收敛条件：
conv_iter：连续迭代次数- 变化阈值：1e-6（可根据数据规模调整）
5. 实际应用扩展
大规模数据优化：
- 使用稀疏相似度矩阵
- 采用近似最近邻(ANN)算法
不同相似度度量：
- 余弦相似度：
S(i,j) = dot(x_i,x_j)/(norm(x_i)*norm(x_j)) - 曼哈顿距离：
S(i,j) = -sum(abs(x_i - x_j))
- 余弦相似度：
处理分类变量：
- 使用Jaccard相似度或汉明距离
动态聚类：
- 增量式AP聚类处理流数据
- 结合时间衰减因子处理时序数据
graph TD
A[生成三维测试数据] --> B[计算相似度矩阵]
B --> C[AP聚类迭代]
C --> D[确定聚类中心]
D --> E[分配点到聚类]
E --> F[可视化结果]
F --> G[评估聚类质量]
该实现完整展示了AP聚类在三维数据分类中的应用，适用于各种三维点云数据处理场景，如：
- 三维物体识别与分割
- 地理空间数据分析
- 生物医学图像处理
- 分子结构分析
- 运动轨迹聚类
