第六章——图像聚类

介绍聚类方法,展示如何利用它们对图像进行聚类,从而寻找相似的图像组。聚类可以用于识别、划分图像数据集,组织与导航。第三节会对聚类后的图像进行相似性可视化。

先来大致了解一下本章的聚类方法:

聚类方法定义思想优点缺点
K-means聚类将输入数据划分成K个簇反复提炼初始评估的类中心适用情形广泛不能保证得到最优结果;
需预先设定聚类数k
层次聚类(或凝聚式聚类)基于样本间成对距离建立一个简相似树利用树结构可以可视化数据间的关系,并显示簇是如何关联的;
一个好的特征向量可给出一个好的分离结果
对于给定的不同的阈值,可直接利用原来的树,不需要重新计算
对于实际需要的聚类簇,需选择一个合适的阈值
谱聚类对该谱矩阵进行特征分解得到的特征向量可用于降维,然后聚类仅需输入相似性矩阵,并用任意度量方式构建相似性矩阵

补充说明:K-means聚类和层次聚类特征向量需求平均,而谱聚类的特征向量无类别限制。

K-means聚类

  1. K-means反复提炼初始评估的类中心的步骤为:

    (1) 以随机或猜测的方式初始化类中心 u i u_i ui i = 1 … k i=1…k i=1k

    (2) 将每个数据点归并到离它距离最近的类中心所属的类 c i c_i ci

    (3) 对所有属于该类的数据点求平均,将平均值作为新的类中心;

    (4) 重复步骤(2)和步骤(3)直到收敛。

    目的是为了使得类内总方差最小:
    V = ∑ i = 1 k ∑ x j ∈ c i ( x i − μ i ) 2 V=\sum_{i=1}^{k} \sum_{x_{j} \in c_{i}}\left(x_{i}-\mu_{i}\right)^{2} V=i=1kxjci(xiμi)2
    该算法通常会初始化不同的类中心进行多次运算,然后选择方差V最小的结果。

    例如,先生成简单的二维正态分布数据,用k=2对数据进行聚类,直到选出方差最小的结果,将其可视化,画出这些数据点及最终的聚类中心。

    # coding=utf-8
    """
    Function:  figure 6.1
        An example of k-means clustering of 2D points
    """
    from pylab import *
    from scipy.cluster.vq import *
    
    # 添加中文字体支持
    from matplotlib.font_manager import FontProperties
    font = FontProperties(fname='SimSun.ttc', size=14)
    
    class1 = 1.5 * randn(100, 2)
    class2 = randn(100, 2) + array([5, 5])
    features = vstack((class1, class2))
    centroids, variance = kmeans(features, 2)
    code, distance = vq(features, centroids)
    figure()
    ndx = where(code == 0)[0]
    plot(features[ndx, 0], features[ndx, 1], 'g*')
    ndx = where(code == 1)[0]
    plot(features[ndx, 0], features[ndx, 1], 'y.')
    plot(centroids[:, 0], centroids[:, 1], 'ro')
    
    title(u'2维数据点聚类', fontproperties=font)
    axis('off')
    show()
    

在这里插入图片描述

  1. 图像聚类

    用K-means对字体图像进行聚类,文件selectedfontimages.zip 包含 66 幅来自该字体数据集 fontinages 的图像,用投影系数作为每幅图像的向量描述符,在主成分上对图像进行投影,然后用下面的方法聚类:

    import imtools
    import pickle
    from pylab import *
    from scipy.cluster.vq import *
    from PIL import Image
    # 获取 selectedfontimages 文件下图像文件名,并保存在列表中
    imlist = imtools.get_imlist('selectedfontimages/a_selected_thumbs')
    imnbr = len(imlist)
    # 载入模型文件
    with open('fontimages/font_pca_modes.pkl','rb') as f:
        immean = pickle.load(f)
        V = pickle.load(f)
    # 创建矩阵,存储所有拉成一组形式后的图像
    immatrix = array([array(Image.open(im)).flatten()
        for im in imlist],'f')
    # 投影到前 40 个主成分上
    immean = immean.flatten()
    projected = array([dot(V[:40],immatrix[i]-immean) for i in range(imnbr)])
    # 进行 k-means 聚类
    projected = whiten(projected)
    centroids,distortion = kmeans(projected,4)
    code,distance = vq(projected,centroids)
    # 绘制聚类簇
    for k in range(4):
        ind = where(code==k)[0]
        figure()
        gray()
        for i in range(minimum(len(ind),40)):
            subplot(4,10,i+1)
            imshow(immatrix[ind[i]].reshape((25,25)))
            axis('off')
    show()
    

    运行结果如下:

在这里插入图片描述

如果设定聚类数k=4,并进行归一化,可视化聚类后的结果为:

在这里插入图片描述

代码为:

 # -*- coding: utf-8 -*-
import imtools, pca
from PIL import Image, ImageDraw
from pylab import *

imlist = imtools.get_imlist('selectedfontimages/a_selected_thumbs')
imnbr = len(imlist)

# Load images, run PCA.
immatrix = array([array(Image.open(im)).flatten() for im in imlist], 'f')
V, S, immean = pca.pca(immatrix)

# Project on 2 PCs.
projected = array([dot(V[[0, 1]], immatrix[i] - immean) for i in range(imnbr)])  # P131 Fig6-3左图
#projected = array([dot(V[[1, 2]], immatrix[i] - immean) for i in range(imnbr)])  # P131 Fig6-3右图

# height and width
h, w = 1200, 1200

# create a new image with a white background
img = Image.new('RGB', (w, h), (255, 255, 255))
draw = ImageDraw.Draw(img)

# draw axis
draw.line((0, h/2, w, h/2), fill=(255, 0, 0))
draw.line((w/2, 0, w/2, h), fill=(255, 0, 0))

# scale coordinates to fit
scale = abs(projected).max(0)
scaled = floor(array([(p/scale) * (w/2 - 20, h/2 - 20) + (w/2, h/2)
                      for p in projected])).astype(int)

# paste thumbnail of each image
for i in range(imnbr):
  nodeim = Image.open(imlist[i])
  nodeim.thumbnail((25, 25))
  ns = nodeim.size
  box = (scaled[i][0] - ns[0] // 2, scaled[i][1] - ns[1] // 2,
         scaled[i][0] + ns[0] // 2 + 1, scaled[i][1] + ns[1] // 2 + 1)
  img.paste(nodeim, box)

img.show()
img.save('pca_font.png')
  1. 像素聚类

    下面是对单幅图像中的像素而非全部图像进行聚类的例子,仅会在RGB三通道的像素值上运用K-means进行聚类。

    在这里采用最近邻插值法,以便在类间进行变换时不需要引入新的像素值。
    在这里插入图片描述

    但如果图像某些区域没有空间一致性,则很难将它们分开,如下方图中小男孩和草坪的图。

    在这里插入图片描述

    代码为:

    # coding=utf-8
    """
    Function: figure 6.4
        Clustering of pixels based on their color value using k-means.
    """
    from scipy.cluster.vq import *
    from scipy.misc import imresize
    from pylab import *
    from PIL import Image
    
    # 添加中文字体支持
    from matplotlib.font_manager import FontProperties
    font = FontProperties(fname='SimSun.ttc', size=14)
    
    steps = 50  # image is divided in steps*steps region
    infile = 'boy_on_hill.jpg'
    im = array(Image.open(infile))
    dx = im.shape[0] / steps
    dy = im.shape[1] / steps
    # compute color features for each region
    features = []
    for x in range(steps):
        for y in range(steps):
            R = mean(im[int(x * dx):int((x + 1) * dx), int(y * dy):int((y + 1) * dy), 0])
            G = mean(im[int(x * dx):int((x + 1) * dx), int(y * dy):int((y + 1) * dy), 1])
            B = mean(im[int(x * dx):int((x + 1) * dx), int(y * dy):int((y + 1) * dy), 2])
            features.append([R, G, B])
    features = array(features, 'f')     # make into array
    # cluster
    centroids, variance = kmeans(features, 3)
    code, distance = vq(features, centroids)
    # create image with cluster labels
    codeim = code.reshape(steps, steps)
    codeim = imresize(codeim, im.shape[:2], 'nearest')
    
    figure()
    ax1 = subplot(121)
    title(u'原图', fontproperties=font)
    #ax1.set_title('Image')
    axis('off')
    imshow(im)
    
    ax2 = subplot(122)
    title(u'聚类后的图像', fontproperties=font)
    #ax2.set_title('Image after clustering')
    axis('off')
    gray()
    imshow(codeim)
    
    show()
    

层次聚类

  1. 建立简相似性树步骤为:

    首先将特征向量距离最近的两个样本归并为一组,并在树中创建一个“平均”节点,将这两个距离最近的样本作为该“平均”节点下的子节点;然后在剩下的包含任意平均节点的样本中寻找下一个最近的对,重复进行前面的操作。

    以下面的例子观察该聚类过程:

    先创建一些二维数据点,在对这些数据点进行聚类,设阈值为5,从列表中提取这些聚类簇,并在控制台打印出来。

    运行结果为:

    在这里插入图片描述

    代码为:

    import hcluster
    
    class1 = 1.5 * randn(100, 2)
    class2 = randn(100, 2) + array([5, 5])
    features = vstack((class1, class2))
    tree = hcluster.hcluster(features)
    clusters = tree.extract_clusters(5)
    print('number of clusters', len(clusters))
    for c in clusters:
        print(c.get_cluster_elements())
    
  2. 图像聚类

    看一个基于图像颜色信息对图像进行聚类的例子。

    下载sunsets.zip,用颜色直方图作为每幅图像的特征向量,在包含日落的文件夹中运行下面代码:

    import os
    from PIL import Image
    import hcluster
    # create a list of images
    from matplotlib.pyplot import *
    from numpy import *
    
    path = 'sunsets\\flickr-sunsets-small\\'
    imlist = [os.path.join(path, f) for f in os.listdir(path) if f.endswith('.jpg')]
    # extract feature vector (8 bins per color channel)
    features = zeros([len(imlist), 512])
    for i, f in enumerate(imlist):
        im = array(Image.open(f))
        # multi-dimensional histogram
        h, edges = histogramdd(im.reshape(-1, 3), 8, normed=True, range=[(0, 255), (0, 255), (0, 255)])
        features[i] = h.flatten()
    tree = hcluster.hcluster(features)
    
    # visualize clusters with some (arbitrary) threshold
    clusters = tree.extract_clusters(0.23 * tree.distance)
    # plot images for clusters with more than 3 elements
    for c in clusters:
        elements = c.get_cluster_elements()
        nbr_elements = len(elements)
        if nbr_elements > 3:
            figure()
            for p in range(minimum(nbr_elements,20)):
                subplot(4, 5, p + 1)
                im = array(Image.open(imlist[elements[p]]))
                imshow(im)
                axis('off')
    show()
    hcluster.draw_dendrogram(tree,imlist,filename='sunset.pdf')
    

    运行结果为:

    在这里插入图片描述

    在这里插入图片描述

    它能够将日落划分为2个簇。

    为了可视化聚类树,同时画出了树状图,如下图所示:

    在这里插入图片描述

    也可对之前的手写字母a进行聚类,运行结果为:

    在这里插入图片描述

谱聚类

  1. 相似矩阵定义:相似矩阵又称为亲和矩阵或距离矩阵,是一个n×n的矩阵,矩阵每个元素表示两两之间的相似性分数。它的特征向量没有类别限制,只要有一个“距离”或“相似性”的概念即可。

    谱聚类的过程:

    创建拉普拉斯矩阵 L = I − D − 1 / 2 S D − 1 / 2 L=I-D^{-1/2}SD^{-1/2} L=ID1/2SD1/2,对于相似性矩阵中的元素 s i j s_{ij} sij要求 s i j ≥ 0 s_{ij}≥0 sij0计算 L 的特征向量,并使用 k 个最大特征值对应的 k 个特征向量,构建出一个特征向量集,从而找到聚类簇。其中,
    D − 1 / 2 = [ 1 d 1 1 d 2 ⋱ 1 d n ] D^{-1 / 2}=\left[\begin{array}{llll} \frac{1}{\sqrt{d_{1}}} & & & \\ & \frac{1}{\sqrt{d_{2}}} & & \\ & & \ddots & \\ & & & \frac{1}{\sqrt{d_{n}}} \end{array}\right] D1/2=d1 1d2 1dn 1
    创建一个矩阵,该矩阵的各列是由之前求出的 k 个特征向量构成,每一行可以看做一个新的特征向量,长度为 k。这些新的特征向量可以用诸如 K-means 方法进行聚类,生成最终的聚类簇。

  2. 举例,我们使用两两间的欧氏距离创建矩阵S,并对k个特征向量用常规的K-means进行聚类,最后绘制出运行后的聚类簇,如下图所示:

    在这里插入图片描述

    代码如下:

    # -*- coding: utf-8 -*-
    import imtools, pca
    from PIL import Image, ImageDraw
    from pylab import *
    from scipy.cluster.vq import *
    
    imlist = imtools.get_imlist('selectedfontimages/a_selected_thumbs')
    imnbr = len(imlist)
    
    # Load images, run PCA.
    immatrix = array([array(Image.open(im)).flatten() for im in imlist], 'f')
    V, S, immean = pca.pca(immatrix)
    
    # Project on 2 PCs.
    projected = array([dot(V[[0, 1]], immatrix[i] - immean) for i in range(imnbr)])  # P131 Fig6-3左图
    # projected = array([dot(V[[1, 2]], immatrix[i] - immean) for i in range(imnbr)])  # P131 Fig6-3右图
    
    n = len(projected)
    # 计算距离矩阵
    S = array([[sqrt(sum((projected[i] - projected[j]) ** 2))
                for i in range(n)] for j in range(n)], 'f')
    # 创建拉普拉斯矩阵
    rowsum = sum(S, axis=0)
    D = diag(1 / sqrt(rowsum))
    I = identity(n)
    L = I - dot(D, dot(S, D))
    # 计算矩阵 L 的特征向量
    U, sigma, V = linalg.svd(L)
    k = 5
    # 从矩阵 L 的前k个特征向量(eigenvector)中创建特征向量(feature vector) # 叠加特征向量作为数组的列
    features = array(V[:k]).T
    # k-means 聚类
    features = whiten(features)
    centroids, distortion = kmeans(features, k)
    code, distance = vq(features, centroids)
    # 绘制聚类簇
    for c in range(k):
        ind = where(code == c)[0]
        figure()
        gray()
        for i in range(minimum(len(ind), 39)):
            im = Image.open(imlist[ind[i]])
            subplot(4, 10, i + 1)
            imshow(array(im))
            axis('equal')
            axis('off')
    show()
    
Logo

为开发者提供学习成长、分享交流、生态实践、资源工具等服务,帮助开发者快速成长。

更多推荐