企业官网建设流程全解析

做聚类最让人纠结的从来不是算法本身，而是动手之前那个灵魂拷问：“你这堆数据，到底应该分几类？”K-means要先拍一个K，K值拍不好后面全是白搭；层次聚类虽然不用K，但树状图怎么看怎么像玄学。近邻传播聚类算法（AP算法）就聪明在这：它不需要你事先给定聚类数目，也不需要初始化聚类中心，只要把样本两两的相似度算好，剩下的事交给消息传递机制自己长出来。这篇东西我会从原理、Matlab手写实现到实测算例完整过一遍，适合刚接触AP算法、想赶紧跑通代码并且理解背后逻辑的人。

我最早接触AP算法是做人脸聚类那会儿，样本量几千个，类别数想提前估都估不准，K-means换了好几个K都出不来稳定的结果。后来换成AP算法，跑出来的exemplar（聚类中心）是真实样本，不是均值合成出来的“平均脸”，这对后续分析特别友好。这几年用它处理过不少数据，也踩过不少坑，下面把能直接复用的经验写出来。

1. 聚类不用猜K：AP算法到底解决了什么痛点

1.1 传统聚类的三个“心累”时刻

先说K-means。你拿一份数据过来，第一件事就是定K。看轮廓系数？画肘部法则？跑一遍要时间，画完图发现肘部拐得不明显，好几个K看起来都像手肘。就算你咬牙定了K，初始聚类中心是随机给的，同样的K跑三遍，可能出来三套不同的结果，初值敏感这个问题到哪都躲不掉。

再说层次聚类。它不用预先指定K，但你需要从树状图里找一个高度去切，这本质上还是在变相选类别数。而且层次聚类算完一张距离矩阵后不停合并，复杂度O(N^2·logN)，数据量一上来消耗很大。

还有一类场景更麻烦：聚类中心本身需要是“有意义的真实样本”。比如图像检索系统里要把相似图片归堆，然后从每堆里选一张作为代表图缓存起来，你肯定希望这张代表图是真实存在的图片，而不是K-means算出来的一个“平均像素模板”。这些需求加起来，就逼着我去找一条更省心的路，AP算法刚好接得住。

1.2 AP算法的核心思想：让样本自己去“竞选代表”

AP算法的全称是Affinity Propagation，中文叫近邻传播聚类。它最核心的概念叫exemplar，也就是聚类中心。注意，这个中心不是靠计算均值得到的“虚拟点”，而是从样本里选出来的一个真实数据点。

算法把所有样本点同时当成两种角色：一种是“选民”，每个点都要找一个自己认为最合适的代表；另一种是“候选人”，每个点都可以站出来说“我可以当大家的代表”。关键的是，算法不提前指定谁当候选人，也不限定候选人数，大家通过一轮一轮的消息传递，自发地形成共识：一部分点出来当聚类中心，剩下的点各自归附到某个中心。所以聚类数目是自动涌现的，不是人为拍脑袋定的。

打个比方，班里要选几个课代表，老师不说选几个，只让大家互相评价“你觉得谁合适、谁对你有吸引力”。每个人把自己的真实想法传出去，经过几轮信息交换，最后自然会有几个威信高的同学成了课代表，其他同学围拢过去形成几个小组。AP算法干的差不多就是这件事，只是它把“看谁顺眼”量化成了相似度。

这套机制的好处很明显：

• 不用预设聚类数目，类别数从数据里长出来；
• 不用初始化聚类中心，不存在随机初值导致结果漂移的问题；
• 聚类中心是真实样本，可解释性强，可以直接拿去做后续分析。

2. 吸引度与归属度：消息传递机制完全拆解

2.1 两个矩阵、两种消息，各管什么事

AP算法的消息传递围着两个矩阵转：吸引度（Responsibility）矩阵R和归属度（Availability）矩阵A，维度都是N×N，N是样本数。

R(i, k)表示的是样本i认为“样本k适合当我的聚类中心”的证据强度，可以理解为i对k的“粉丝值”。这个值越大，说明i越认可k。注意它是单方向的，i对k的认可不等于k对i的认可。

A(i, k)表示的是样本k根据所有其他样本对它的认可情况，给出“我愿意作为样本i的聚类中心”的证据强度。可以反向理解为“k愿意收下i当自己人”的程度。这个值越大，说明k作为i的中心的“底气”越足。

迭代过程中，R和A交替更新、互相促进：R告诉A“我多喜欢k”，A再告诉R“k有多愿意接纳大家”。两个矩阵收敛稳定后，每个样本i取A(i,k)+R(i,k)最大的那个k作为自己的聚类中心，如果恰好有样本满足“自己选自己”，那它就是最终的exemplar，也就是聚类中心。

2.2 更新公式的直觉推导：为什么r是“相似度减最强劲敌”

R的更新公式长这样：

r(i,k) = s(i,k) − max_{j≠k} { a(i,j) + s(i,j) }

这里的s(i,k)是样本i和样本k的相似度，一般在算法运行前由你算好，默认越大代表越相似。整个公式能翻译成人话：i对k的认可度，等于“我和k本身有多像”，再减去“我在其他候选中心那里能拿到的最大诱惑”。如果k跟i很亲近，但还有另一个候选点比k更讨i欢心，那么这个“最强劲敌”的A+S值就会很高，把r(i,k)压下去。

这里有三个细节值得说：

• A(i,j)+S(i,j)可以理解为“i对j的综合好感度”，前半段是j对其他人的欢迎程度，后半段是i和j的直接相似度；
• max的下标排除j=k，也就是不跟k自己比，而要去比较其他选项；
• 当k=i时，这个公式自动变成r(i,i)=s(i,i)−max_{j≠i}{a(i,j)+s(i,j)}，s(i,i)其实就是偏好值p，它的大小直接决定了i“自我推荐当中心”的底气。

A的更新公式对应两条：

a(i,k) = min{0, r(k,k) + Σ_{j≠i,k} max(0, r(j,k))}（i≠k）

a(k,k) = Σ_{j≠k} max(0, r(j,k))

第一条的意思是：k想收i当成员，前提是k自己愿意当中心（r(k,k)），而且有不少其他点也认可k（那些正的r(j,k)），但最终要“克制”一点，用min(0, ...)限制住，避免出现一个样本疯狂簇拥很多中心的混乱局面。第二条则完全放开，把所有正向的认可加总，作为k自我定位为中心的依据。A更新完成后，k要真能成为一个中心，必须是“自己信得过自己，并且大家也愿意跟它”。

2.3 阻尼系数λ：防止震荡的“刹车片”

两个矩阵这么来回传消息，容易刹不住车，表现就是R和A来回震荡、迟迟不收敛。AP算法专门加了一个阻尼系数λ（典型值0.9）来做平滑：

R_new = (1−λ)·R_calc + λ·R_old

A_new = (1−λ)·A_calc + λ·A_old

公式摆出来其实就是在旧值和新计算值之间做加权平均。λ越接近1，更新越平稳，收敛越慢；λ越接近0，更新冲劲越大，但越容易震荡。我实际跑数据的经验是：当迭代到最大次数还不收敛时，最优先的操作就是把λ从0.9慢慢提到0.95、0.98，而不是去改数据或者改偏好值，大部分震荡问题都能解决。

3. Matlab手写AP聚类：不用工具箱也能跑

网上能搜到Frey教授公布的AP工具箱，功能很全，但我一直建议想弄清楚AP算法的人至少手写一遍核心循环。原因很简单：工具箱是个黑盒，你调参时不知道里面发生了什么；自己写一遍，才知道哪些参数带来了哪些影响，遇到问题也更好定位。好在AP的Matlab实现并不复杂，核心迭代也就二三十行。

3.1 输入准备：相似度矩阵和偏好参数

算法的第一份输入是相似度矩阵S，N×N。S(i, j)越大表示样本i和j越相似。最常用的构造方式是负平方欧氏距离：

S(i,j) = −‖x_i − x_j‖²

负号保证“距离越小，相似度越大”。也有人用余弦相似度、高斯核相似度，关键是所有相似度要同向比较，越大越相似，这个必须全局一致。

第二份输入是偏好值p，也就是矩阵对角线上的值S(k,k)。p表示每个样本对自己成为聚类中心的“先天意愿”。p给得越大，越多的点想自立门户，聚类数目就越多；p给得越小，大家越谦让，聚类数目就越少。默认推荐是取相似度矩阵所有非对角线元素的中位数，这时产生的聚类数目通常比较合理。我见过很多人把这个参数忽略掉，实际上它才是AP算法的“旋钮”，后面第四章会专门讲怎么用它控制聚类数。

代码里我把数据规范化放在主脚本中做，算法本体不掺和，保证函数只负责迭代：

% 数据标准化：消除量纲影响 X = zscore(X);

注意，如果用负欧氏距离构造相似度矩阵，数据标准化这一步很关键。否则一个量纲很大的特征会主导整个距离计算，聚类结果基本就跟着那个特征走了。

3.2 核心迭代循环的实现要点

下面是我整理的一份脱离工具箱的手写版AP函数，可以直接复制保存成apcluster_manual.m使用。为了照顾可读性，尽量保留了计算过程的原有逻辑：

function [idx, centers] = apcluster_manual(S, p, maxits, convits, lam) % AP聚类算法手工实现（仿照Frey & Dueck 2007） % 输入： % S : N×N 相似度矩阵，S(i,j)越大表示i和j越相似 % p : 偏好值，标量或N×1向量，p越大聚类数越多 % maxits : 最大迭代次数，默认300 % convits : 收敛判定窗口，连续convits次聚类中心不变即停止，默认15 % lam : 阻尼系数，默认0.9 % 输出： % idx : N×1列向量，idx(i)是样本i所属聚类中心的样本编号 % centers : 聚类中心（exemplar）编号列表 N = size(S, 1); if isscalar(p) p = p * ones(N, 1); end S(1:N+1:end) = p; % 将偏好值写入对角线 R = zeros(N, N); A = zeros(N, N); stable_count = 0; stable_centers = []; centers = []; for iter = 1:maxits Rold = R; Aold = A; %% ---------- 更新吸引度 R ---------- % 技巧：每行先求最大值，再求第二大值 % 对大部分列直接用"行最大值"即可，只有最大值所在列要换成第二大值 B = A + S; [m1, k1] = max(B, [], 2); % 每行最大值及列索引 B2 = B; for i = 1:N B2(i, k1(i)) = -inf; % 遮住最大值位置 end [m2, ~] = max(B2, [], 2); % 每行第二大值 Rnew = S - repmat(m1, 1, N); for i = 1:N Rnew(i, k1(i)) = S(i, k1(i)) - m2(i); end R = (1 - lam) * Rnew + lam * Rold; %% ---------- 更新归属度 A ---------- Rpos = max(R, 0); % 只保留正向吸引度 sumR = sum(Rpos, 1); % 每列求和 Rdiag = diag(R); % 取对角线 Anew = min(0, repmat(Rdiag + sumR(:), 1, N) - Rpos); Anew(1:N+1:end) = sumR(:); % 对角线单独赋值 A = (1 - lam) * Anew + lam * Aold; %% ---------- 收敛判定 ---------- [~, idx] = max(A + R, [], 2); centers_now = find(idx == (1:N)'); % 选出自我认可的样本 if numel(centers_now) == numel(stable_centers) && ... isequal(sort(centers_now), sort(stable_centers)) stable_count = stable_count + 1; if stable_count >= convits centers = centers_now; break; end else stable_count = 1; stable_centers = centers_now; end end end

其中吸引度更新的“最大/第二大值”技巧值得单独说一下：标准公式里，每个r(i,k)都要排除j=k，再扫一遍最大值，复杂度是O(N³)，直接写循环会慢得你想哭。上面的写法是先把每行最大值m1算出来，除最大值所在的列外，其余列的r(i,k)直接用S(i,k)−m1；只有那一列需要换成第二大值m2。这样整段更新变成O(N²)，代码也简短得多。

3.3 收敛判断：看“中心集合”而不是“矩阵数值”

很多初学AP的人会把收敛条件写成“R矩阵和A矩阵变化很小”，这是一个性价比很低的做法，因为矩阵数值层面哪怕有一点点波动，最终聚类中心集合也可能早就稳定了。反过来也一样，矩阵“看起来”还在变，但中心点集合可能已经连续好几轮没变过。我用的判据是：每次迭代后把每个样本的归属中心idx算出来，再筛选出exemplar集合centers_now，这个集合连续convits轮完全不变，就判定收敛。这也和Frey原版工具箱的收敛标准保持一致。

需要注意，convits不要设成1或2，否则很容易在震荡初期误判为收敛。我自己一般固定15，最大迭代次数300。

4. 跑通第一个例子：三簇高斯数据聚类

4.1 完整可运行的主脚本

纸上谈兵没意思，直接构造一个三簇高斯混合数据来验证。每个簇50个样本，总共150个点，维度2，方便可视化：

%% 1. 生成三簇高斯数据 rng(42); X = [randn(50,2)*0.6 + repmat([ 1, 1], 50, 1); randn(50,2)*0.6 + repmat([-1, 1], 50, 1); randn(50,2)*0.6 + repmat([ 0, -1.2], 50, 1)]; %% 2. 标准化 X = zscore(X); %% 3. 相似度矩阵：负平方欧氏距离 S = -pdist2(X, X).^2; %% 4. 默认偏好：中位数 p = median(S(:)); %% 5. 运行AP聚类 [idx, centers] = apcluster_manual(S, p, 300, 15, 0.9); %% 6. 可视化 figure; gscatter(X(:,1), X(:,2), idx); hold on; plot(X(centers,1), X(centers,2), 'ks', ... 'MarkerSize', 14, 'LineWidth', 2, 'MarkerFaceColor', 'k'); title('AP聚类结果（黑色方块为自动选出的聚类中心）'); xlabel('x1'); ylabel('x2'); grid on;

我用pdist2是因为它能一次性算出两两距离矩阵，比嵌套循环快而且写法干净。跑完后你会看到三类簇基本被分开，黑色方块落在每个簇的“众望所归”位置。这里有一个值得注意的现象：三个簇给的都是高斯分布，理论上K=3是正确答案，但AP算法并没有刻意去找“三个类”，它只是通过消息传递自然地收敛出类似结构，聚类数目是结果，不是预设。如果你把簇之间的间隔调得再近一些，AP可能自动合并成两类，这是很正常的。

4.2 参数p与簇数量的动态关系

默认的中位数偏好只是一个通用起点，并不保证一定符合你期望的类别数。真正的调参手段是扫描p。p越大，每个点越愿意“自己当中心”，聚类数越多；p越小，样本越互相谦让，聚类数越少。这个趋势非常稳定，我每次拿到新数据都会先画一张“p-聚类数曲线”看看整体走势：

%% 扫描p，观察聚类数变化 pvals = linspace(min(S(:)), max(S(:)), 25); ncvals = zeros(size(pvals)); for i = 1:length(pvals) [idx, ~] = apcluster_manual(S, pvals(i), 300, 15, 0.9); ncvals(i) = numel(unique(idx)); end figure; plot(pvals, ncvals, 'o-', 'LineWidth', 1.5); xlabel('偏好值 p'); ylabel('聚类数目 K'); title('p-聚类数关系曲线'); grid on;

这条曲线几乎总是单调递减或单调递增（具体方向取决于p的定义），非常直观。跑一遍只需要几秒钟，却能让你对本份数据的“聚类敏感性”一目了然。我在实际项目中，会先跑这个扫描图，而不是直接拿中位数跑一次就交差。看完曲线，你脑海里对“p该往哪边调”就有了确切的判断。

4.3 用二分法自动找到你想要的簇数

如果你所在业务里聚类数目确实有明确要求，比如“我要5类”，不用手搓p试来试去，直接用二分法在扫描区间里找p。当K与p满足单调关系时，二分法几十次迭代内就能拿到一个满足要求或非常接近的p值：

%% 二分法：寻找聚类数接近targetK的偏好值 targetK = 5; plo = min(S(:)); phi = max(S(:)); for t = 1:50 pmid = (plo + phi) / 2; [idx, ~] = apcluster_manual(S, pmid, 300, 15, 0.9); K = numel(unique(idx)); if K > targetK phi = pmid; % p偏大，聚类数过多，下调上界 elseif K < targetK plo = pmid; % p偏小，聚类数过少，上调下界 else break; end end fprintf('找到的p = %.4f，对应聚类数 = %d\n', pmid, K);

这段代码需要注意，聚类结果本身存在一定随机性？其实AP是确定性算法，只要相似度矩阵和p相同，结果就是可复现的，所以二分是稳定的。唯一的麻烦是单调性可能存在局部平台，二分过程可能落在一个邻近值上，但一般差别不大。我个人经验是，如果你需要精确的K，二分法跑完后，再在得到的p附近左右微调几步，通常就能稳定命中。

5. 踩坑经验：震荡、内存爆炸和结果不稳

5.1 迭代不收敛：先调λ，再调convits

AP算法在数据类别比较少或者相似度矩阵区分度不强的时候，特别容易进入震荡。最明显的症状是跑到maxits=300次还没停止，警告级别的事件。这时候第一反应不是去乱改p，而是把λ从0.9往上提。我做过对比测试：一组数据在λ=0.9时跑了300次不收敛，改成0.95后大约80次收敛；改成0.98需要150次左右，但全程非常平稳。代价就是收敛变慢，所以λ不是越大越好，震荡不严重的话保持0.9就够了。

另一个相关经验是，检查一下是不是convits设太宽导致迟迟不触发停止。如果你设置convits=50，那即使已经很稳定了，也要连续50轮完全一致才停，这会让运行时间变得非常长。一般情况下convits=10~20即可。

5.2 相似度矩阵的尺度影响

AP算法对相似度的数值尺度比K-means敏感得多。原因是偏好值p是对角线上的一个“标尺”，而相似度矩阵的整体分布直接决定这个标尺的位置。举例来说，如果你把余弦相似度（取值范围−1到1）拿来直接用，然后又用中位数当p，出来的聚类数通常会很少；而换用负平方欧氏距离后，矩阵的中位数可能是−几或−几十，两者对应的类别数可能差很多。

所以不要换了一种相似度定义，还沿用之前的p值经验。每次换相似度度量，都要重新跑一遍扫描曲线，重新选定适合自己的偏好。对于距离类相似度，我习惯先把原始特征标准化到零均值单位方差，再算负平方欧氏距离，这样尺度和方向都比较稳。

5.3 数据量变大：内存、稀疏化和替代思路

AP算法的标准实现里需要存一个N×N的相似度矩阵和两个N×N的消息矩阵，内存随样本量平方增长。N=2000时，三个双精度矩阵就是96MB左右，还勉强能跑；N=10000时就是2.4GB，彻底吃不消。

我遇到几千样本量时一般不会硬上全矩阵，而是用稀疏策略：

• 只保留每个样本与最近t个邻居的相似度，其余置为−inf或0；
• 相似度矩阵用Matlab的sparse格式存储；
• 在R和A的更新过程中同样限制在邻居范围内，本质上是“局部消息传播”。

如果你不想自己改稀疏版，一个更省事的替代是先用KNN图做一次粗聚类，或者用子采样跑一版AP，把大概的聚类结构摸清楚，再对剩余样本做最近邻指派。这个思路在处理上百万样本的项目里特别实用，毕竟AP算法本身的定位本来就是中等规模数据——几千个样本跑起来很舒服，几万个样本就要认真规划资源了。我自己之前做微动面波数据的事件聚类时，样本量也在千级别，直接全矩阵AP完全够用，再大的数据集建议走稀疏化。

5.4 聚类结果“太离谱”时怎么诊断

有时跑完发现所有样本并成一类，或者每个样本都自称中心，这时候不要急着骂算法，一般先看两个点。

第一，p的取值范围是否偏离了相似度矩阵的分布。p如果远小于矩阵所有非对角元素，相当于大家都在说“我不想当中心”，最后自然合成一个巨大聚类；p如果接近甚至超过矩阵最大值，相当于人人都想当中心，结果就是几十上百个碎片簇。

第二，相似度是否构造错误。最常见的问题是把“小=好”的距离直接当相似度传入S，导致算法逻辑完全反着来。AP要求S越大代表越相似，如果你手里只有距离矩阵D，务必先做S=−D（或−D.^2）。

6. AP与K-means、DBSCAN的横向对比：什么场景真正适合AP

6.1 一张表看懂三种算法

这张表基本能说明问题。K-means胜在快和简单，但输入K和初值两大主观因素绕不开；DBSCAN好处是不需要K且支持任意形状，但密度参数eps对数据分布很敏感，密度不均时非常难调；AP则牺牲了一些规模和速度，换来了“不用预先指定聚类数目”和“中心是真实样本”两个硬优势。

6.2 适合AP算法的四个信号

结合这些年的使用体验，我一般这么判断要不要上AP：

• 样本量适中，千级到万级出头，太大就去想稀疏化方案；
• 类别数确实难以预估，行业内也没有标准答案；
• 你需要聚类中心作为后续决策依据，比如选代表图片、挑选典型用户画像；
• 两两之间天然可以定义“相似度”指标，而不是必须靠距离公式强行转换。

如果四个信号全中，AP基本就是最优解。如果只是给上百万条无标注日志快速分个桶，K-means明显更合适。

6.3 个人使用心得与建议

每次跑AP，我最不省心的其实不是算法，而是相似度矩阵的设计。AP本身对相似度质量非常依赖，相似度定义得糙，后面再怎么调p也救不回来。所以我建议拿到数据后别急着写代码，先花时间思考“什么样算相似”这件事，把这个问题想透了，AP的聚类质量通常不会差。

刚上手AP算法时，可以从一个小技巧切入：先跑一遍默认中位数偏好的结果，画出散点图，再看一眼当前的聚类数目。然后扫描一遍p-聚类数关系曲线，结合业务背景选一个合理的聚类粒度。整个过程不需要反复试K，也不用担心初始化漂移，确实比K-means省心得多。这份Matlab代码我用了很久，改造成Python或者R版本也不难——核心就是那两个消息矩阵的交替更新，理解了它，你也就在聚类这件事上省掉了一大半烦恼。