漫谈 Clustering (2): k-Mysql Tutorial-php.cn

原文：http://blog.pluskid.org/?p=40 上一次我们了解了一个最基本的 clustering 办法 k-means ，这次要说的 k-medoids 算法，其实从名字上就可以看出来，和 k-means 肯定是非常相的。事实也确实如此，k-medoids 可以算是 k-means 的一个变种。 k-medoids 和

原文：http://blog.pluskid.org/?p=40

上一次我们了解了一个最基本的 clustering 办法 k-means ，这次要说的 k-medoids 算法，其实从名字上就可以看出来，和 k-means 肯定是非常相似的。事实也确实如此，k-medoids 可以算是 k-means 的一个变种。

k-medoids 和 k-means 不一样的地方在于中心点的选取，在 k-means 中，我们将中心点取为当前 cluster 中所有数据点的平均值：

漫谈 Clustering (2): k

Rough Collie

并且我们已经证明在固定了各个数据点的 assignment 的情况下，这样选取的中心点能够把目标函数漫谈 Clustering (2): k 最小化。然而在 k-medoids 中，我们将中心点的选取限制在当前 cluster 所包含的数据点的集合中。换句话说，在 k-medoids 算法中，我们将从当前 cluster 中选取这样一个点——它到其他所有（当前 cluster 中的）点的距离之和最小——作为中心点。k-means 和 k-medoids 之间的差异就类似于一个数据样本的均值 (mean) 和中位数 (median) 之间的差异：前者的取值范围可以是连续空间中的任意值，而后者只能在给样本给定的那些点里面选。那么，这样做的好处是什么呢？
一个最直接的理由就是 k-means 对数据的要求太高了，它使用欧氏距离描述数据点之间的差异 (dissimilarity) ，从而可以直接通过求均值来计算中心点。这要求数据点处在一个欧氏空间之中。

然而并不是所有的数据都能满足这样的要求，对于数值类型的特征，比如身高，可以很自然地用这样的方式来处理，但是类别 (categorical) 类型的特征就不行了。举一个简单的例子，如果我现在要对犬进行聚类，并且希望直接在所有犬组成的空间中进行，k-means 就无能为力了，因为欧氏距离漫谈 Clustering (2): k 在这里不能用了：一只Samoyed 减去一只 Rough Collie 然后在平方一下？天知道那是什么！再加上一只 German Shepherd Dog 然后求一下平均值？根本没法算，k-means 在这里寸步难行！

在 k-medoids 中，我们把原来的目标函数漫谈 Clustering (2): k 中的欧氏距离改为一个任意的 dissimilarity measure 函数：

<img  src="/inc/test.jsp?url=http%3A%2F%2Fblog.pluskid.org%2Flatexrender%2Fpictures%2Fc2f42fa0d2b5b49f31e8a7459af89a4e.png&refer=http%3A%2F%2Fblog.csdn.net%2Fzhazhiqiang%2Farticle%2Fdetails%2F19554235" alt="漫谈 Clustering (2): k" >

Copy after login

最常见的方式是构造一个 dissimilarity matrix 漫谈 Clustering (2): k 来代表，其中的元素表示第只狗和第只狗之间的差异程度，例如，两只 Samoyed 之间的差异可以设为 0 ，一只 German Shepherd Dog 和一只 Rough Collie 之间的差异是 0.7，和一只 Miniature Schnauzer 之间的差异是 1 ，等等。

除此之外，由于中心点是在已有的数据点里面选取的，因此相对于 k-means 来说，不容易受到那些由于误差之类的原因产生的 Outlier 的影响，更加 robust 一些。

扯了这么多，还是直接来看看 k-medoids 的效果好了，由于 k-medoids 对数据的要求比 k-means 要低，所以 k-means 能处理的情况自然 k-medoids 也能处理，为了能先睹为快，我们偷一下懒，直接在中的 k-means 代码的基础上稍作一点修改，还用同样的例子。将代码的 45 到 47 行改成下面这样：

Copy after login

        <span><strong>for</strong></span> j <span><strong>in</strong></span> <span>range</span>(k):
            idx_j = (labels == j).nonzero()
            distj = distmat(X[idx_j], X[idx_j])
            distsum = ml.<span>sum</span>(distj, axis=<span>1</span>)
            icenter = distsum.argmin()
            centers[j] = X[idx_j[<span>0</span>][icenter]]

Copy after login

可以看到 k-medoids 在这个例子上也能得到很好的结果：

漫谈 Clustering (2): k

而且，同 k-means 一样，运气不好的时候也会陷入局部最优解中：

漫谈 Clustering (2): k

如果仔细看上面那段代码的话，就会发现，从 k-means 变到 k-medoids ，时间复杂度陡然增加了许多：在 k-means 中只要求一个平均值漫谈 Clustering (2): k 即可，而在 k-medoids 中则需要枚举每个点，并求出它到所有其他点的距离之和，复杂度为。

看完了直观的例子，让我们再来看一个稍微实际一点的例子好了：Document Clustering ——那个永恒不变的主题，不过我们这里要做的聚类并不是针对文档的主题，而是针对文档的语言。实验数据是从 Europarl 下载的包含 Danish、German、Greek、English、Spanish、Finnish、French、Italian、Dutch、Portuguese 和 Swedish 这些语言的文本数据集。

在 N-gram-based text categorization 这篇 paper 中描述了一种计算由不同语言写成的文档的相似度的方法。一个（以字符为单位的） N-gram 就相当于长度为 N 的一系列连续子串。例如，由 hello 产生的 3-gram 为：hel、ell 和 llo ，有时候还会在划分 N-gram 之前在开头和末尾加上空格（这里用下划线表示）：_he、hel、ell、llo、lo_ 和 o__ 。按照 Zipf’s law ：

The nth most common word in a human language text occurs with a frequency inversely proportional to n.

这里我们用 N-gram 来代替 word 。这样，我们从一个文档中可以得到一个 N-gram 的频率分布，按照频率排序一下，只保留频率最高的前 k 个（比如，300）N-gram，我们把叫做一个“Profile”。正常情况下，某一种语言（至少是西方国家的那些类英语的语言）写成的文档，不论主题或长短，通常得出来的 Profile 都差不多，亦即按照出现的频率排序所得到的各个 N-gram 的序号不会变化太大。这是非常好的一个性质：通常我们只要各个语言选取一篇（比较正常的，也不需要很长）文档构建出一个 Profile ，在拿到一篇未知文档的时候，只要和各个 Profile 比较一下，差异最小的那个 Profile 所对应的语言就可以认定是这篇未知文档的语言了——准确率很高，更可贵的是，所需要的训练数据非常少而且容易获得，训练出来的模型也是非常小的。

不过，我们这里且撇开分类（Classification）的问题，回到聚类（Clustering）上，按照前面的说法，在 k-medoids 聚类中，只需要定义好两个东西之间的距离（或者 dissimilarity ）就可以了，对于两个 Profile ，它们之间的 dissimilarity 可以很自然地定义为对应的 N-gram 的序号之差的绝对值，在 Python 中用下面这样一个类来表示：

Copy after login

<span><strong>class</strong></span> Profile(<span>object</span>):
    <span><strong>def</strong></span> <span>__init__</span>(<span>self</span>, path, N=<span>3</span>, psize=<span>400</span>):
        <span>self</span>.N = N
        <span>self</span>.psize = psize
        <span>self</span>.build_profile(path)
 
    sep = <span>re</span>.<span>compile</span>(r<span>'<span><strong>\W</strong></span>+'</span>)
    <span><strong>def</strong></span> build_profile(<span>self</span>, path):
        grams = {}
        <span><strong>with</strong></span> <span>open</span>(path) <span><strong>as</strong></span> inf:
            <span><strong>for</strong></span> line <span><strong>in</strong></span> inf:
                <span><strong>for</strong></span> tok <span><strong>in</strong></span> <span>self</span>.sep.split(line):
                    <span><strong>for</strong></span> n <span><strong>in</strong></span> <span>range</span>(<span>self</span>.N):
                        <span>self</span>.feed_ngram(grams, tok, n+<span>1</span>)
        <span>self</span>.create_profile(grams.items())
 
    <span><strong>def</strong></span> create_profile(<span>self</span>, grams):
        <span><em># keep only the top most psize items</em></span>
        grams.sort(key=itemgetter(<span>1</span>), reverse=<span>True</span>)
        grams = grams[:<span>self</span>.psize]
 
        <span>self</span>.<span>profile</span> = <span>dict</span>()
        <span><strong>for</strong></span> i <span><strong>in</strong></span> <span>range</span>(<span>len</span>(grams)):
            <span>self</span>.<span>profile</span>[grams[i][<span>0</span>]] = i
 
    <span><strong>def</strong></span> <span>__getitem__</span>(<span>self</span>, key):
        idx = <span>self</span>.<span>profile</span>.get(key)
        <span><strong>if</strong></span> idx <span><strong>is</strong></span> <span>None</span>:
            <span><strong>return</strong></span> <span>len</span>(<span>self</span>.<span>profile</span>)
        <span><strong>return</strong></span> idx
 
    <span><strong>def</strong></span> dissimilarity(<span>self</span>, o):
        <span>dis</span> = <span>0</span>
        <span><strong>for</strong></span> tok <span><strong>in</strong></span> <span>self</span>.<span>profile</span>.keys():
            <span>dis</span> += <span>abs</span>(<span>self</span>[tok]-o[tok])
        <span><strong>for</strong></span> tok <span><strong>in</strong></span> o.<span>profile</span>.keys():
            <span>dis</span> += <span>abs</span>(<span>self</span>[tok]-o[tok])
        <span><strong>return</strong></span> <span>dis</span>
 
    <span><strong>def</strong></span> feed_ngram(<span>self</span>, grams, tok, n):
        <span><strong>if</strong></span> n <span>!</span>= <span>0</span>:
            tok = <span>'_'</span> + tok
        tok = tok + <span>'_'</span> <span>*</span> (n-<span>1</span>)
        <span><strong>for</strong></span> i <span><strong>in</strong></span> <span>range</span>(<span>len</span>(tok)-n+<span>1</span>):
            gram = tok[i:i+n]
            grams.setdefault(gram, <span>0</span>)
            grams[gram] += <span>1</span>

Copy after login

europarl 数据集共有 11 种语言的文档，每种语言包括大约 600 多个文档。我为这七千多个文档建立了 Profile 并构造出一个 7038×7038 的 dissimilarity matrix ，最后在这上面用 k-medoids 进行聚类。构造 dissimilarity matrix 的过程很慢，在我这里花了将近 10 个小时。相比之下，k-medoids 的过程在内存允许的情况下，采用向量化的方法来做实际上还是很快的，并且通常只要数次迭代就能收敛了。实际的 k-medoids 实现可以在 mltk 中找到，今后如果有时间的话，我会陆续地把一些相关的比较通用的代码放到那里面。

Hungarian algorithm 来求解。

我们这里有 11 种语言，全排列有 11! = 39916800 种情况，对于每一种排列，我们需要遍历一次 label list ，并数出真正的 label （语言）与聚类得出的结果相同的文档的个数，再除以总的文档个数，得到 accuracy 。假设每次遍历并求出 accuracy 只需要 1 毫秒的时间的话，总共也需要 11 个小时才能得到结果。看上去好像也不是特别恐怖，不过相比起来，用 Hungarian algorithm 的话，我们可以几乎瞬间得到结果。由于文章的篇幅已经很长了，就不在这里介绍具体的算法了，感兴趣的同学可以参考 Wikipedia ，这里我直接使用了一个现有的 Python 实现。

虽然这个实验非常折腾，不过最后的结果其实就是一个数字：accuracy ——在我这里达到了 88.97% ，证明 k-medoids 聚类和 N-gram Profile 识别语言这两种方法都是挺不错的。最后，如果有感兴趣的同学，代码可以从这里下载。需要最新版的 scipy， munkres.py 和 mltk 以及 Python 2.6 。