k-medoids聚类算法,即k-中心聚类算法,它是基于k-means聚类算法的改进。我们知道,k-means算法执行过程,首先需要随机选择初始质心,只有第一次随机选择的初始质心才是实际待聚类点集中的点,而后续将非质心点指派到对应的质心点后,重新计算得到的质心并非是待聚类点集中的点,而且如果某些非质心点是离群点的话,导致重新计算得到的质心可能偏离整个簇,为了解决这个问题,提出了改进的k-medoids聚类算法。
k-medoids聚类算法也是通过划分的方式来计算得到聚类结果,它使用绝对差值和(Sum of Absolute Differences,SAD)来衡量聚类结果的优劣,在n维欧几里德空间中,计算SAD的公式如下所示:
围绕中心点划分(Partitioning Around Medoids,PAM)的方法是比较常用的,使用PAM方法进行处理,可以指定一个最大迭代次数的参数,在迭代过程中基于贪心策略来选择使得聚类的质量最高的划分。PAM的方法,每次交换一个中心点和非质心点,进行指派以后,计算得到的SAD值越小,则聚类质量越好,如此不断地迭代,直到找到一个最好的划分。
在维基百科上,它给出的基于PAM方法计算聚类的过程,描述如下:
上面算法描述,应该是按照某种顺序的取中心点集合中的点,也从非中心点集合中取点,分别计算生成新的划分的代价,我们可以考虑,每次取点的时候,采用随机取点的策略,随机性越强越好,只要满足最终迭代终止的条件即可。通常,如果能够迭代所有情况,那么最终得到的划分一定是最优的划分,即聚类结果最好,这通常适用于比较小的点的集合,但是如果待聚类的点比较大,则需要通过限制迭代次数来终止迭代计算,从而得到一个能够满足实际需要的聚类结果。我们在下面的实现中,选择随机取中心点和非中心点进行交换,通过设置允许最大迭代次数这个参数值,来使聚类迭代计算停止。
首先,便于理解后面的实现过程,我们描述一下代码实现聚类过程的基本逻辑,如下所示:
我们实现的k-medoids聚类算法,需要指定2个聚类相关参数,另外一个参数是程序计算并行度,可以通过构造方法看到,代码如下所示:
public KMedoidsClustering(int k, int maxIterations, int parallism) { super(k, maxIterations, parallism); distanceCache = new DistanceCache(Integer.MAX_VALUE); executorService = Executors.newCachedThreadPool(new NamedThreadFactory("SEEKER")); latch = new CountDownLatch(parallism); }
上面代码中的参数含义如下:
聚类实现的核心代码如下所示:
@Override public void clustering() { // parse sample files FileUtils.read2DPointsFromFiles(allPoints, "[/t,;//s]+", inputFiles); // 从文件读取点数据,加入到集合allPoints中 LOG.info("Total points: count=" + allPoints.size()); ClusterHolder currentHolder = new ClusterHolder(); // 每一次迭代过程中的需要的数据结构都封装到ClusterHolde对象中 ClusterHolder previousHolder = null; currentHolder.medoids = initialCentroidsSelectionPolicy.select(k, allPoints); // 随机选择初始中心 LOG.info("Initial selected medoids: " + currentHolder.medoids); // start seeker threads for (int i = 0; i < parallism; i++) { // 启动parallism个线程,执行非中心点到中心点的指派 final NearestMedoidSeeker seeker = new NearestMedoidSeeker(seekerQueueSize); executorService.execute(seeker); seekers.add(seeker); } // ///////////////// // make iterations // ///////////////// boolean firstTimeToAssign = true; int numIterations = 0; double previousSAD = 0.0; double currentSAD = 0.0; try { while(!finallyCompleted) { try { LOG.debug("Current medoid set: " + currentHolder.medoids); if(firstTimeToAssign) { // 第一次处理时,只是根据随机选择的初始中心集合,和全部点的集合,指派给多个线程处理 assignNearestMedoids(currentHolder, true); firstTimeToAssign = false; } else { // 非第一次处理时,每次迭代得到的聚类结果,都是基于中心点进行分组的,处理逻辑稍微不同 assignNearestMedoids(currentHolder, false); } // merge result mergeMedoidAssignedResult(currentHolder); // 每个线程处理一部分,最后要合并多个线程分别处理的结果 LOG.debug("Merged result: " + currentHolder.medoidWithNearestPointSet); // compare cost for 2 iterations, we use SAD (sum of absolute differences) if(previousSAD == 0.0) { // first time compute SAD previousSAD = currentSAD; currentSAD = computeSAD(currentHolder); // 第一次计算SAD } else { RandomPoint randomPoint = selectNonCenterPointRandomly(currentHolder); // 随机选择一个非中心点 LOG.debug("Randomly selected: " + randomPoint); // compute current cost when using random point to substitute for the medoid currentSAD = computeSAD(currentHolder); // // 计算用随机选择非中心点替换一个中心点得到的SAD值 // compare SADs if(currentSAD - previousSAD < 0.0) { // 如果此次迭代得到的SAD值,比上次迭代计算得到SAD小,替换previousHolder和previousSAD,以保证最终算法终止后,该最小值对应的划分能够保留下来 previousHolder = currentHolder; previousSAD = currentSAD; } // construct new cluster holder currentHolder = constructNewHolder(currentHolder, randomPoint); // 根据随机选择的中心点,创建一个新的 ClusterHolde对象,用于下次迭代 } LOG.info("Iteration #" + (++numIterations) + ": previousSAD=" + previousSAD + ", currentSAD=" + currentSAD); if(numIterations > maxIterations) { // 如果达到指定的最大迭代次数,则终止 finallyCompleted = true; } } catch(Exception e) { Throwables.propagate(e); } finally { try { if(!finallyCompleted) { latch = new CountDownLatch(parallism); completeToAssignTask = false; } Thread.sleep(10); synchronized(signalLock) { signalLock.notifyAll(); } } catch (InterruptedException e) {} } } } finally { LOG.info("Shutdown executor service: " + executorService); executorService.shutdown(); } // finally result centerPointSet.addAll(previousHolder.medoids); // 处理最终的聚类结果 Iterator<Entry<CenterPoint, List<Point2D>>> iter = previousHolder.medoidWithNearestPointSet.entrySet().iterator(); while(iter.hasNext()) { Entry<CenterPoint, List<Point2D>> entry = iter.next(); int clusterId = entry.getKey().getId(); Set<ClusterPoint<Point2D>> set = Sets.newHashSet(); for(Point2D p : entry.getValue()) { set.add(new ClusterPoint2D(p, clusterId)); } clusteredPoints.put(clusterId, set); } }
通过上面代码及其注释,我们可以了解到基本的处理流程,首先看一下几个工具类ClusterHolder和RandomPoint:
private class ClusterHolder { /** snapshot of clustering result: medoids of clustering result, as well as non-medoid points */ private TreeMap<CenterPoint, List<Point2D>> medoidWithNearestPointSet; /** center point set represented by Point2D */ private Set<Point2D> centerPoints; /** center point set represented by CenterPoint */ private TreeSet<CenterPoint> medoids; public ClusterHolder() { super(); } } private class RandomPoint { /** medoid which the random point belongs to */ private final CenterPoint medoid; // 随机选择的中心点 /** a non-medoid point selected randomly */ private final Point2D point; // 随机选择的非中心点,该点被指派给上面的中心点medoid public RandomPoint(CenterPoint medoid, Point2D point) { super(); this.medoid = medoid; this.point = point; } @Override public String toString() { return "RandomPoint[medoid=" + medoid + ", point=" + point + "]"; } }
下面我们看一下上面代码调用的比较重要的方法的实现逻辑。
将非中心点指派到最近的中心点的计算,是调用assignNearestMedoids方法,该方法的代码实现,如下所示:
private void assignNearestMedoids(final ClusterHolder holder, boolean firstTimeToAssign) { LOG.debug("firstTimeToAssign=" + firstTimeToAssign); try { // assign tasks to seeker threads if(firstTimeToAssign) { // 第一次进行指派,因为还没有进行指派过,所以只有随机选择的一组中心点,和全部待聚类的点的集合 holder.centerPoints = Sets.newHashSet(); for(CenterPoint medoid : holder.medoids) { holder.centerPoints.add(medoid.toPoint()); // 构造ClusterHolder对象,将中心点加入到集合中 } LOG.debug("holder.centerPoints: " + holder.centerPoints); for(Point2D p : allPoints) { // 对全部待聚类的点作为任务,加入到每个线程的队列中,但是要排除已经被选择为中心点的点 LOG.debug("Assign point: " + p); if(!holder.centerPoints.contains(p)) { selectSeeker().q.put(new Task(holder.medoids, p)); } } } else { for(List<Point2D> points : holder.medoidWithNearestPointSet.values()) { // 如果笔试第一次进行指派,已经在构造ClusterHolder对象的时候,将随机选择的中心点和非中心点进行了交换,这里直接进行指派即可 for(Point2D p : points) { selectSeeker().q.put(new Task(holder.medoids, p)); } } } } catch(Exception e) { Throwables.propagate(e); } finally { try { completeToAssignTask = true; latch.await(); } catch (InterruptedException e) { } } }
上面代码调用selectSeeker()方法,获取到一个NearestMedoidSeeker线程,将待指派的点加入到其队列中,selectSeeker()方法实现代码如下所示:
private NearestMedoidSeeker selectSeeker() { int index = taskIndex++ % parallism; return seekers.get(index); }
下面,我们看一下NearestMedoidSeeker线程的实现,它也比较简单,实现了将队列q中的点取出,计算到该点最近的中心,然后指派给该中心,线程实现代码如下所示:
private class NearestMedoidSeeker implements Runnable { private final Log LOG = LogFactory.getLog(NearestMedoidSeeker.class); private final BlockingQueue<Task> q; private Map<CenterPoint, List<Point2D>> clusteringNearestPoints = Maps.newHashMap(); private int processedTasks = 0; public NearestMedoidSeeker(int qsize) { q = new LinkedBlockingQueue<Task>(qsize); } @Override public void run() { while(!finallyCompleted) { // 每一轮迭代,调用一次assign方法 try { assign(); Thread.sleep(200); } catch (Exception e) { e.printStackTrace(); } } } private void assign() throws InterruptedException { try { LOG.debug("Q size: " + q.size()); while(!(q.isEmpty() && completeToAssignTask)) { processedTasks++; final Task task = q.poll(); if(task != null) { final Point2D p1 = task.point; double minDistance = Double.MAX_VALUE; CenterPoint nearestMedoid = null; for(CenterPoint medoid : task.medoids) { final Point2D p2 = medoid.toPoint(); Double distance = distanceCache.computeDistance(p1, p2); // 计算非中心点p1到中心点p2的欧几里德距离 if(distance < minDistance) { minDistance = distance; nearestMedoid = medoid; } } LOG.debug("Nearest medoid seeked: point=" + p1 + ", medoid=" + nearestMedoid); List<Point2D> points = clusteringNearestPoints.get(nearestMedoid); if(points == null) { points = Lists.newArrayList(); clusteringNearestPoints.put(nearestMedoid, points); } points.add(p1); // 将非中心点p1,指派给到中心点的欧几里德距离最近的点 } else { Thread.sleep(150); } } } catch (Exception e) { e.printStackTrace(); } finally { latch.countDown(); LOG.debug("Point processed: processedTasks=" + processedTasks); synchronized(signalLock) { signalLock.wait(); } clusteringNearestPoints = Maps.newHashMap(); processedTasks = 0; } } }
每一轮指派,多个线程都计算得到一个非中心点指派的子集,最后还要将这些子集合并为一个全局的指派结果,合并的实现在mergeMedoidAssignedResult()方法中,代码如下所示:
private void mergeMedoidAssignedResult(ClusterHolder currentHolder) { currentHolder.medoidWithNearestPointSet = Maps.newTreeMap(); for(NearestMedoidSeeker seeker : seekers) { LOG.debug("seeker.clusteringNearestPoints: " + seeker.clusteringNearestPoints); Iterator<Entry<CenterPoint, List<Point2D>>> iter = seeker.clusteringNearestPoints.entrySet().iterator(); while(iter.hasNext()) { Entry<CenterPoint, List<Point2D>> entry = iter.next(); List<Point2D> set = currentHolder.medoidWithNearestPointSet.get(entry.getKey()); if(set == null) { set = Lists.newArrayList(); currentHolder.medoidWithNearestPointSet.put(entry.getKey(), set); } set.addAll(entry.getValue()); } } }
合并后的指派结果,都存放在ClusterHolder对象中,为下一轮迭代准备了数据。
随机选择一个中心点和非中心点,并交换,实现代码如下所示:
private RandomPoint selectNonCenterPointRandomly(ClusterHolder holder) { List<CenterPoint> medoids = new ArrayList<CenterPoint>(holder.medoidWithNearestPointSet.keySet()); CenterPoint selectedMedoid = medoids.get(random.nextInt(medoids.size())); // 随机选择一个中心点 List<Point2D> belongingPoints = holder.medoidWithNearestPointSet.get(selectedMedoid); Point2D point = belongingPoints.get(random.nextInt(belongingPoints.size())); // 随机选择一个非中心点 return new RandomPoint(selectedMedoid, point); // 返回这2个点 }
因为每一次迭代,我们都得到一个非中心点指派到最近的中心点的聚类结果集合,所以在设计随机选择中心点和非中心点进行交换时,我们首先从中心点集合中选择一个中心点,然后再从该中心点对应的非中心点的簇的集合中选择一个非中心点,当然也可以考虑用其他的方法,比如,待交换的中心点和非中心点不在同一个簇中。
我们处理的策略是,事后处理,也就是每次先实现非中心点和中心点的交换,再进行指派,计算SAD值,如果此轮迭代得到的SAD值比上一轮的大,则直接丢弃结果,将上一轮的指派结果作为最终候选结果,直到最后,保留着具有最小SAD值的指派结果。创建ClusterHolder对象,交换非中心点和中心点的实现逻辑,如下所示:
private ClusterHolder constructNewHolder(final ClusterHolder holder, RandomPoint randomPoint) { ClusterHolder newHolder = new ClusterHolder(); // collect center points with type Point2D for a holder object // from previous result of clustering procedure newHolder.centerPoints = Sets.newHashSet(); for(CenterPoint c : holder.medoidWithNearestPointSet.keySet()) { newHolder.centerPoints.add(c.toPoint()); } Point2D newPoint = randomPoint.point; CenterPoint oldMedoid = randomPoint.medoid; // create a new center point with type CenterPoint based on the randomly selected non-medoid point // and it's id is equal to the old medoid's CenterPoint newMedoid = new CenterPoint(oldMedoid.getId(), newPoint); // use new medoid above to substitute the old medoid newHolder.centerPoints.remove(oldMedoid.toPoint()); newHolder.centerPoints.add(newPoint); newHolder.medoids = Sets.newTreeSet(); newHolder.medoids.addAll(holder.medoidWithNearestPointSet.keySet()); newHolder.medoids.remove(oldMedoid); // remove old medoid from center point set of new holder object newHolder.medoids.add(newMedoid); // copy the holder's medoidWithNearestPointSet, and modify it newHolder.medoidWithNearestPointSet = Maps.newTreeMap(); newHolder.medoidWithNearestPointSet.putAll(holder.medoidWithNearestPointSet); List<Point2D> oldPoints = newHolder.medoidWithNearestPointSet.get(oldMedoid); oldPoints.remove(newPoint); // remove new randomly selected non-medoid point from previous result set of clustering oldPoints.add(oldMedoid.toPoint()); // add old medoid point to the non-medoid set newHolder.medoidWithNearestPointSet.put(newMedoid, oldPoints); return newHolder; }
为了保留上一次迭代指派的结果,这里不要修改holder对应的结果的集合,而是拷贝出一份,在拷贝的结果上交换中心点和非中心点。
我们将使用k-medoids算法与k-means算法聚类,分别对结果进行比较。其中k-means算法由于随机选择初始质心,每次执行聚类结果不同,而k-medoids算法聚类依赖于与迭代次数,所以我们选择maxIterations分别为300、1000、3000时,对比效果,如下图所示:
通过上图可以看出,k-medoids聚类算法,当maxIterations越大的时候,可能更加靠近最优解,聚类结果的质量越高,也就是SAD值越小。
本文基于 署名-非商业性使用-相同方式共享 4.0 许可协议发布,欢迎转载、使用、重新发布,但务必保留文章署名时延军(包含链接),不得用于商业目的,基于本文修改后的作品务必以相同的许可发布。如有任何疑问,请与我联系。