
简介AP聚类Affinity Propagation算法无需预先指定聚类数量通过消息传递自动确定簇结构在数据挖掘与模式识别等领域应用广泛。该Python实现完整覆盖算法核心流程适合学习聚类原理、开展无监督学习实验或在实际项目中快速接入AP聚类的开发者使用。压缩包内仅含1个文件apcluster.py体积约1KB代码简洁集中便于逐行阅读和二次修改。目前已有781人浏览学习。文件实现了从数据读取、距离矩阵计算到职责与可用性消息迭代更新、聚类中心判定及结果输出的完整环节并兼顾了输入检查与稀疏矩阵计算等细节可帮助读者直观理解AP算法从公式到工程落地的全过程。通过动手运行和调整参数还能进一步掌握责任矩阵R与可用性矩阵A的收敛机制为后续扩展更大规模数据集或对比其他聚类算法打下坚实基础。1. 一个不知道分几类的问题为什么先说AP聚类先看一个实际场景你拿到一批用户行为特征想分群做运营但没人告诉你该分成几类。跑k-means之前你得先定K手肘图在簇之间有重叠时给出的K值模棱两可即便定下K初始中心选偏一次收敛结果就翻一次脸。AP聚类Affinity Propagation亲和传播聚类在这个问题上更省心的原因在于它把所有样本点都当作候选的簇代表点exemplar通过样本两两之间传递“我推荐谁”“谁愿意接收我”两类消息自动收敛出聚类数量和每个样本的归属。标题里的apcluster.zip、apcluster这类命名落到Python环境里核心就是把这套消息传递迭代写成可运行代码。如果你手头数据量在几千条以内又想省去反复试探K的过程这篇文章会从算法判定逻辑讲到参数调整和结果验证全程给出可复现的命令。2. 先看懂AP聚类的核心机制两条消息的交替收敛2.1 exemplar是从样本里选出来的不是算出来的中心k-means的簇中心是均值向量这个词有两个隐含问题第一它是个“合成点”不一定落在真实样本上第二数据形状只要偏离团状高斯分布均值中心就开始失真。AP聚类换了一个思路——簇中心必须从原始样本里选选出来的那位就叫exemplar代表点。这个设计带来的直接收益是业务可解释性。用户分群完成后不需要再描述“第3类用户坐标均值是xxx”直接拿着代表点的画像说“第3类用户就是这个样本的样子”这对运营、推荐和风控场景都友好得多。代价是算法的比较逻辑变了不再计算点到假设中心的距离而是比较“候选代表k对样本i的吸引程度”和“候选代表k愿意接收多少成员”。为了做到这一点每个样本都同时扮演两种角色既是被分配的普通成员又是潜在的簇代表。最终结果是每个样本指向唯一一个exemplar并且这个指向关系不能成环只会收敛成一棵棵以exemplar为根的星形结构。2.2 responsibility样本对候选代表的“竞争后吸引力”第一类消息叫responsibility记作 r(i,k)含义是“样本i认为候选k适不适合当自己代表”。它的计算公式依赖相似度矩阵相似度 s(i,k) 在欧氏空间里通常取负的距离平方即 s(i,k) -||xi - xk||²值越大表示越相近。r(i,k) 的递推公式是r(i,k) s(i,k) - max_{k≠k} { a(i,k) s(i,k) }这个式子的意思是把k和所有其他候选k′放在一起比k有多少“超出竞争对手”的优势。注意右边的max要在去掉k之后取否则k自己的得分也被拿来和自己比结果会失真。如果r(i,k)大于0说明在排除k之后k依然是i最有吸引力的选择如果小于0说明存在至少一个候选比k更有资格。全部样本、全部候选计算一遍后每个样本i能看出“谁最值得跟随”。2.3 availability候选代表反过来决定“收不收你”光是样本觉得谁好还不够候选代表得有意愿接收。第二类消息availability记作 a(i,k)表示候选k对样本i的接收意愿它分两种写法。对 i≠k 的情况a(i,k) min( 0, r(k,k) sum_{i≠i,k} max(0, r(i,k)) )对 k自己a(k,k) sum_{i≠k} max(0, r(i,k))r(k,k)是k自荐当代表的决心后面那串求和是所有其他样本给k的正向responsibility之和——有多少人在推举k。这个总和越高说明k当代表的支持声越大。但a(i,k)被截断到0以下的原因也很直观一个候选代表的容量有限它如果同时答应所有样本簇就会膨胀成整块数据。当竞争激烈时它对某个特定样本的接收意愿就得降下来把这个名额让给更需要的样本。responsibility负责横向比“谁的吸引力最强”availability负责纵向比“谁推举我、我容不容得下”两条消息交错更新。更新一轮后要做阻尼混合避免相邻两轮结果跳变r_new (1 - λ) · r_calculated λ · r_old a_new (1 - λ) · a_calculated λ · a_oldλ就是sklearn里的damping参数取值必须在0.5到1之间。λ越接近1新旧更新之间的步伐越慢越不容易振荡越接近0.5更新越快但数据分布复杂时越容易陷入来回横跳。对比项responsibility r(i,k)availability a(i,k)传递方向样本i → 候选k候选k → 样本i回答的问题k对i来说够不够格k愿不愿意接收i主要依赖s(i,k)与其他候选的比较结果r(k,k)和其他样本对k的推荐正值含义k领先其他候选k有余力接纳i更新时需要排除排除k自身排除i和k2.4 收敛条件和“谁当代表”的最后判定迭代不会无限进行。每轮更新完算法会检查所有样本当前的exemplar指派是否连续convergence_iter轮都没有变化如果没到收敛条件但轮数先达到max_iter迭代也会停止。最终判定某个样本k是否真的成为exemplar看的是对角线上 r(k,k)a(k,k) 是否大于0也就是“自己推举自己”和“别人推举自己”的总分。所有exemplar加上它们吸引到的样本就组成最终的聚类结果。这里有一个实际运行中容易遇到的现象两个样本互相指定对方当exemplar形成互指环。遇到这种情况先怀疑是不是相似度矩阵选得不对或preference设得太极端导致算法陷在局部稳定状态。后续章节会在参数调节部分给出具体解法。3. Python里把AP聚类跑通先调sklearn再还原公式3.1 用scikit-learn跑通一次AP聚类的最小命令先确认Python环境能import numpy和scikit-learn各平台安装方式不同这里不展开。装依赖的命令是pip install scikit-learn numpy matplotlib下面是最小可运行代码用make_blobs造300个二维样本点里面真实簇数是4看AP聚类在不知道簇数的情况下能找回几个簇。from sklearn.cluster import AffinityPropagation from sklearn.datasets import make_blobs import numpy as np X, _ make_blobs(n_samples300, centers4, cluster_std0.8, random_state42) ap AffinityPropagation( damping0.9, preferenceNone, # None表示用相似度矩阵中位数 max_iter200, convergence_iter15, random_state42 ) y_pred ap.fit_predict(X) unique, counts np.unique(y_pred, return_countsTrue) print(聚类数:, len(unique)) print(每个簇的样本量:, counts) print(代表点索引:, ap.cluster_centers_indices_) print(实际迭代轮数:, ap.n_iter_)fit_predict在内部完成了相似度矩阵计算和消息迭代。preferenceNone表示使用相似度矩阵的中位数作为初始自荐值这个默认值在中型数据集上通常会给出偏多的簇比如通常会产生20个上下的小簇对300个点来说很容易超过真实簇数4。damping0.9是先把振荡风险压住避免第一轮结果不可复现。random_state42的作用是让结果稳定可复现不过在数据本身足以区分簇时AP聚类对随机种子并不敏感只有在preference恰好在多个候选点之间并列时才容易出现微小差异。把输出的聚类数和centers4对比会发现大概率落在4到10之间。这正是AP聚类的特点它不需要你给K但它给什么K取决于preference怎么设。如何把K压到目标区间看第4章。3.2 照公式手写核心循环确认自己真懂了再调参sklearn的实现为了效率和内存做了大量优化读源码容易被细节干扰。想确认自己对算法的理解照公式写一个无优化的参考实现就够了。下面这段代码严格按2.2和2.3的递推公式执行没有任何技巧性加速。import numpy as np def ap_reference(X, damping0.9, max_iter200, convergence_iter15, preferenceNone): n len(X) S -((X[:, None, :] - X[None, :, :]) ** 2).sum(axis2) if preference is not None: np.fill_diagonal(S, preference) R np.zeros((n, n)) A np.zeros((n, n)) center_history [] for it in range(max_iter): R_prev R.copy() A_prev A.copy() # 更新responsibility for i in range(n): A_plus_S A_prev[i] S[i] for k in range(n): others np.delete(A_plus_S, k) # 严格排除k自己 R[i, k] S[i, k] - others.max() # 更新availability for k in range(n): pos np.maximum(R[:, k], 0) sum_pos pos.sum() - pos[k] # 排除k本人 A[k, k] sum_pos for i in range(n): if i k: continue A[i, k] min(0.0, R[k, k] sum_pos - pos[i]) # 阻尼混合 R damping * R_prev (1 - damping) * R A damping * A_prev (1 - damping) * A # 连续convergence_iter轮exemplar集合不变则收敛 scores R.diagonal() A.diagonal() centers frozenset(np.where(scores 0)[0].tolist()) center_history.append(centers) if len(center_history) convergence_iter: recent center_history[-convergence_iter:] if len(set(recent)) 1: break scores R.diagonal() A.diagonal() labels (A S).argmax(axis1) return labels, np.where(scores 0)[0]这段代码有三个关键细节值得单独说明。第一更新R时用np.delete把k从竞争中剔除了这一步直接对应公式里的“k′≠k”。网上流传的若干简化版直接对整行取最大值在k本身就是最大竞争者时会把r(i,k)压低最终形成更多碎片簇。第二更新A时先算出所有样本对k的正向responsibility总和再依次减掉k自己和当前样本i的贡献和公式里的排除逻辑完全一致。第三收敛判断比较的是每一轮exemplar集合是否完全相同不是R或A矩阵数值是否接近因为聚类任务关心的是成员指派关系稳定。需要顺手提醒的是这个参考实现的复杂度大约O(n³)300个样本没问题超过1000个点就会慢到难以等待。生产环境请继续用sklearn参考实现只用于理解逻辑和验证自己的推导。3.3 相似度矩阵默认是负欧氏平方也能换成预计算矩阵AP聚类一切迭代都发生在相似度矩阵上。sklearn在fit_predict(X)时内部计算的是负欧氏距离平方多数数值型特征场景都适用。但如果特征本身是文本向量、点击序列或其他需要特殊度量的数据可以自己算相似度矩阵用affinityprecomputed传进去。from sklearn.metrics.pairwise import cosine_similarity S_cos cosine_similarity(X) # 余弦相似度值域[-1,1] p np.percentile(S_cos, 25) # 比中位数更严苛的自荐门槛 ap AffinityPropagation(affinityprecomputed, preferencep, damping0.9) y_pred ap.fit_predict(S_cos)提示affinityprecomputed时preference参数会直接写入相似度矩阵对角线。你也可以不传preference此时sklearn会用矩阵中位数填充对角线。很多人在这一步踩坑以为precomputed模式可以完全不管对角线结果收敛出的聚类数莫名其妙偏多。另外手写相似度矩阵时注意大n下的内存占用。上面参考实现用广播生成(n, n, dim)三阶数组样本过万时内存会直接爆掉这种场景应改用sklearn.metrics.pairwise.euclidean_distances或按块计算。4. AP聚类调参preference、damping、max_iter与一个自适应循环4.1 preference是决定聚类数的总开关preference表示每个样本“自荐当代表”的先验得分它直接写进相似度矩阵对角线。preference越大样本越容易自立为代表聚类数就越多preference越小代表名额越稀缺聚类数就越少。这是整个AP聚类里最值得花时间调的参数。sklearn默认取值是中位数这在大几百样本的数据上一般会给出几十个簇远多于实际业务期望。想少分几类就把preference调小常见做法是从相似度矩阵的10%到25%分位数开始尝试。比如想分成5到10簇可以先拿到相似度矩阵算一下最小值、10%分位数、中位数各是多少然后从10%分位数附近起步观察输出的聚类数。preference还可以做成非标量给每个样本不同的自荐分。比如在风控场景中已知某些样本明显是“典型代表”可以把这些样本的对角线分数提高引导算法优先选它们当代表。这是k-means实现不了的能力。4.2 damping控制振荡调参不能拉的旋钮damping取值在0.5到1之间默认0.5。它的本质是上一轮消息和本轮消息的混合比例damping0.9的含义是保留90%的旧值、加入10%的新值。数据本身分离度好时0.5也能稳定收敛一旦簇之间有重叠、或preference设得较极端消息就可能在几个候选点之间来回跳表现为max_iter耗尽且聚类数不稳定。遇到振荡先把damping提到0.9还不行就0.99。代价是收敛速度变慢原来50轮能收敛的问题可能要300轮。因此调大damping时一定要同步把max_iter放开否则新的报错就是“迭代轮数耗尽”。4.3 max_iter和convergence_iter判断收敛的标准别用默认值裸跑这两个参数决定算法什么时候停。convergence_iter15的含义是连续15轮exemplar指派完全不变才宣告收敛。这个阈值在数据干净时很宽松在噪声数据上15轮可能过于草率——有时中间出现一次偶然波动刚好打断了连续15轮的记录又要重新数。遇到这类情况把convergence_iter提到30结果会更稳。判断当前参数跑没跑好直接看ap.n_iter_。如果输出等于max_iter说明撞到了迭代上限优先检查damping是否过低而不是无脑加大max_iter。加大max_iter只是给振荡更长的时间去跳问题本身没有解决。下面是几个参数的速查表参数sklearn默认作用调参方向preference相似度矩阵中位数样本自荐当代表的门槛想少聚类就调小想多聚类就调大damping0.5R与A的阻尼混合比例振荡时提到0.9~0.99max_iter200最大迭代轮数高damping时扩到500~1000convergence_iter15连续多少轮exemplar不变视为收敛噪声数据可提到304.4 按目标聚类数自动搜索preference的参考脚本既然preference和聚类数呈单调负相关就可以用二分查找把聚类数压进目标区间。下面这个循环每次对半分preference区间跑完看聚类数是多了还是少了最多30轮能收敛到目标区间。def search_preference(X, target_k(5, 15), damping0.9): S -((X[:, None, :] - X[None, :, :]) ** 2).sum(axis2) lo, hi S.min(), np.median(S) # 下界几乎不分簇上界是默认行为 for _ in range(30): p (lo hi) / 2 ap AffinityPropagation(affinityprecomputed, preferencep, dampingdamping, max_iter500, random_state0) ap.fit(S) k len(ap.cluster_centers_indices_) if target_k[0] k target_k[1]: return p, ap if k target_k[1]: hi p # 簇太多说明preference偏大往下压 else: lo p # 簇太少说明preference偏小往上抬 return lo, None这里二分区间选在“相似度最小值”和“中位数”之间是因为最小值附近几乎不会有样本自荐成功聚类数接近1中位数是默认行为聚类数偏多。实际调用时如果返回的ap是None说明30轮内没有命中目标区间可以采用最后一次的lo作为preference继续跑或者放宽target_k。这段脚本在样本量几千、特征维度几十的场景下运行时间可接受样本上到几万单次拟合就会变慢需要在第5章给出的采样方案下使用。5. 先用轮廓系数验证AP聚类质量再用KNN把大样本扩出来5.1 聚类质量验证轮廓系数和ARIAP聚类调完参数不等于结果能直接用先用指标验证一次。轮廓系数不需要真实标签适合做无监督评估如果数据有真实类别用Adjusted Rand Index更直接。from sklearn.metrics import silhouette_score, adjusted_rand_score sil silhouette_score(X, y_pred) print(f轮廓系数: {sil:.3f}) if y_true is not None: ari adjusted_rand_score(y_true, y_pred) print(fARI: {ari:.3f})轮廓系数接近1说明簇内紧凑、簇间分离接近0或负值则说明簇边界模糊。AP聚类在这种评估下通常比k-means更能体现非凸形状的数据结构但这不代表它可以无视数据预处理特征量纲差异大的场景必须先标准化否则相似度矩阵会被量级大的特征主导。5.2 大样本的处理随机子集做AP再用KNN外推AP聚类每次迭代都在更新n×n的消息矩阵时间和空间复杂度都是O(n²)样本到两万以上会非常吃力。常见的落地方案是先随机抽样一至两千个样本做AP聚类得到代表点和簇结构再用KNN把剩余样本归类。from sklearn.neighbors import KNeighborsClassifier idx np.random.choice(len(X), size2000, replaceFalse) sample X[idx] ap_sample AffinityPropagation(damping0.9, max_iter300).fit(sample) sample_labels ap_sample.labels_ knn KNeighborsClassifier(n_neighbors5).fit(sample, sample_labels) full_labels knn.predict(X)抽样规模建议控制在2000以内AP部分能秒级收敛KNN只做最近邻查找样本量再大也能扛。这个方案牺牲了一小部分聚类精度换来了对全量数据的可伸缩性也是AP聚类在生产环境中最常见的用法。如果业务要求每个簇都必须有真实代表点抽样得到的exemplar仍然来自原始样本天然满足这个约束。把这段KNN外推逻辑包成一个函数配合5.1的轮廓系数脚本就是一个可以直接用于线下分析和线上打标的AP聚类流水线。本文还有配套的精品资源点击获取