python计算机视觉编程第六章——图像聚类

一、K-means聚类

K-means是一种将输入数据划分成k个簇的简单的聚类算法。K-means反复提炼初始评估的类中心,步骤如下:

  1. 以随机或猜测的方式初始化类中心u,,i=1…k;
  2. 将每个数据点归并到离它距离最近的类中心所属的类c;
  3. 对所有属于该类的数据点求平均,将平均值作为新的类中心;
  4. 重复步骤(2)和步骤(3)直到收敛。
     

K-means 试图使类内总方差最小:

V=\sum_{i=1}^{k}\sum_{x_{j\in c_{i}}}(x_{j}-\mu _{i})^{2}

1 SciPy聚类包 

为便于说明,我们先生成简单的二维数据: 

from scipy.cluster.vq import *
from pylab import *

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],'*')
ndx = where(code==1)[0]
plot(features[ndx,0],features[ndx,1],'r.')
plot(centroids[:,0],centroids[:,1],'go')
axis('off')
show()
运行结果:类中心标记为绿色大圆环,预测出的 类分别标记为蓝色星号和红色点。

2 图像聚类

我们利用之前计算过的前 40 个主 成分进行投影,用投影系数作为每幅图像的向量描述符。用 pickle 模块载入模型文 件,在主成分上对图像进行投影,然后用下面的方法聚类:
import imtools
import pickle
from scipy.cluster.vq import *

from PIL import Image
from pylab import *

# 获取 selected-fontimages 文件下图像文件名,并保存在列表中
imlist = imtools.get_imlist('E:/data/a_selected_thumbs')
imnbr = len(imlist)
# 载入模型文件
with open('../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()

 运行结果:

3在主成分上可视图像 

改变投影为:
projected = array([dot(V[[0,2]],immatrix[i]-immean) for i in range(imnbr)])
PIL 中的 ImageDraw 模块进行可视化:
# 高和宽
h,w = 1200,1200
# 创建一幅白色背景图
img = Image.new('RGB',(w,h),(255,255,255))
draw = ImageDraw.Draw(img)
# 绘制坐标轴
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 = abs(projected).max(0)
scaled = floor(array([ (p / scale) * (w/2-20,h/2-20) + (w/2,h/2) for p in projected]))
# 粘贴每幅图像的缩略图到白色背景图片
for i in range(imnbr):
 nodeim = Image.open(imlist[i])
 nodeim.thumbnail((25,25))
 ns = nodeim.size
 img.paste(nodeim,(int(scaled[i][0]-ns[0]//2),int(scaled[i][1]-ns[1]//2),int(scaled[i][0]+ns[0]//2+1),int(scaled[i][1]+ns[1]//2+1)))
img.save('pca_font.jpg')
在成对主成分上投影的字体图像。左图用的是第一个和第二个主成分,右图用的是
第二个和第三个主成分:

 

4 像素聚类

下面的代码示例载入一幅图像,用一个步长为 steps 的方形网格在图像中滑动,每
滑一次对网格中图像区域像素求平均值,将其作为新生成的低分辨率图像对应位置
处的像素值,并用 K-means 进行聚类:
from PIL import Image, ImageDraw
from pylab import *
from scipy.cluster.vq import *
# from scipy.misc import imresize

steps = 50 # 图像被划分成 steps×steps 的区域
im = array(Image.open('E:/data/empire.jpg'))
dx = int(im.shape[0] / steps)
dy = int(im.shape[1] / steps)

# 计算每个区域的颜色特征
features = []
for x in range(steps):
 for y in range(steps):
    R = mean(im[x*dx:(x+1)*dx,y*dy:(y+1)*dy,0])
    G = mean(im[x*dx:(x+1)*dx,y*dy:(y+1)*dy,1])
    B = mean(im[x*dx:(x+1)*dx,y*dy:(y+1)*dy,2])
    features.append([R,G,B])
features = array(features,'f') # 变为数组
# 聚类
centroids,variance = kmeans(features,3)
code,distance = vq(features,centroids)
# 用聚类标记创建图像
codeim = code.reshape(steps,steps)
#codeim = imresize(codeim,im.shape[:2],interp='nearest')
# im = imresize(im, (h, int(w * aspect_ratio)), interp='bicubic')
codeim = np.array(Image.fromarray(codeim).resize(im.shape[:2]))

figure()
imshow(codeim)
show()

 

二、层次聚类

层次聚类(或凝聚式聚类)是另一种简单但有效的聚类算法,其思想是基于样本间 成对距离建立一个简相似性树。该算法首先将特征向量距离最近的两个样本归并为 一组,并在树中创建一个“平均”节点,将这两个距离最近的样本作为该“平均” 节点下的子节点;然后在剩下的包含任意平均节点的样本中寻找下一个最近的对, 重复进行前面的操作。在每一个节点处保存了两个子节点之间的距离。遍历整个树, 通过设定的阈值,遍历过程可以在比阈值大的节点位置终止,从而提取出聚类簇。
from itertools import combinations
from pylab import *


class ClusterNode(object):
    def __init__(self, vec, left, right, distance=0.0, count=1):
        self.left = left
        self.right = right
        self.vec = vec
        self.distance = distance
        self.count = count  # 只用于加权平均

    def extract_clusters(self, dist):
        """ 从层次聚类树中提取距离小于 dist 的子树簇群列表 """
        if self.distance < dist:
            return [self]
        return self.left.extract_clusters(dist) + self.right.extract_clusters(dist)

    def get_cluster_elements(self):
        """ 在聚类子树中返回元素的 id """
        return self.left.get_cluster_elements() + self.right.get_cluster_elements()

    def get_height(self):
        """ 返回节点的高度,高度是各分支的和 """
        return self.left.get_height() + self.right.get_height()

    def get_depth(self):
        """ 返回节点的深度,深度是每个子节点取最大再加上它的自身距离 """
        return max(self.left.get_depth(), self.right.get_depth()) + self.distance


class ClusterLeafNode(object):
    def __init__(self, vec, id):
        self.vec = vec
        self.id = id

    def extract_clusters(self, dist):
        return [self]

    def get_cluster_elements(self):
        return [self.id]

    def get_height(self):
        return 1

    def get_depth(self):
        return 0

    def L2dist(v1, v2):
        return sqrt(sum((v1 - v2) ** 2))

    def L1dist(v1, v2):
        return sum(abs(v1 - v2))

    def hcluster(features, distfcn=L2dist):

        """ 用层次聚类对行特征进行聚类 """

        # 用于保存计算出的距离
        distances = {}
        # 每行初始化为一个簇
        node = [ClusterLeafNode(array(f), id=i) for i, f in enumerate(features)]
        while len(node) > 1:
            closest = float('Inf')
            # 遍历每对,寻找最小距离
            for ni, nj in combinations(node, 2):
                if (ni, nj) not in distances:
                    distances[ni, nj] = distfcn(ni.vec, nj.vec)
                d = distances[ni, nj]
                if d < closest:
                    closest = d
                    lowestpair = (ni, nj)
            ni, nj = lowestpair
            # 对两个簇求平均
            new_vec = (ni.vec + nj.vec) / 2.0
            # 创建新的节点
            new_node = ClusterNode(new_vec, left=ni, right=nj, distance=closest)
            node.remove(ni)
            node.remove(nj)
            node.append(new_node)
        return node[0]

运行结果:

number of clusters 2
[168, 192, 128, 145, 143, 149, 116, 144, 130, 199, 118, 123, 103, 181, 165, 193, 117, 159, 190, 135, 174, 109, 108, 184, 137, 171, 115, 177, 189, 197, 146, 166, 161, 183, 141, 152, 179, 107, 195, 125, 188, 111, 127, 140, 191, 113, 100, 121, 136, 185, 194, 178, 101, 163, 169, 150, 126, 182, 142, 158, 186, 106, 151, 147, 112, 167, 132, 164, 134, 162, 120, 157, 170, 104, 129, 155, 156, 176, 133, 122, 153, 160, 198, 131, 124, 148, 105, 180, 139, 154, 102, 196, 175, 119, 114, 173, 138, 110, 187]
[52, 94, 33, 39, 90, 23, 28, 49, 69, 91, 26, 30, 71, 89, 66, 3, 27, 76, 87, 93, 96, 172, 61, 59, 84, 72, 46, 43, 50, 15, 42, 55, 78, 16, 14, 79, 88, 63, 29, 56, 95, 8, 44, 57, 10, 45, 2, 7, 1, 19, 62, 0, 99, 17, 81, 18, 80, 20, 74, 6, 98, 12, 13, 75, 21, 67, 9, 25, 54, 85, 41, 51, 64, 31, 73, 22, 32, 92, 58, 83, 35, 77, 86, 4, 70, 40, 65, 37, 82, 97, 48, 5, 47, 11, 24, 36, 38, 68, 34, 53, 60]

图像聚类

def draw_dendrogram(node, imlist, filename='clusters.jpg'):
    """ 绘制聚类树状图,并保存到文件中 """
    # 高和宽
    rows = node.get_height() * 20
    cols = 1200
    # 距离缩放因子,以便适应图像宽度
    s = float(cols - 150) / node.get_depth()
    # 创建图像,并绘制对象
    im = Image.new('RGB', (cols, rows), (255, 255, 255))
    draw = ImageDraw.Draw(im)
    # 初始化树开始的线条
    draw.line((0, rows / 2, 20, rows / 2), fill=(0, 0, 0))
    # 递归地画出节点
    node.draw(draw, 20, (rows / 2), s, imlist, im)
    im.save(filename)
    im.show()

#画树状图时对于每个节点用了 draw() 方法

    def draw(self, draw, x, y, s, imlist, im):
        """ 用图像缩略图递归地画出叶节点 """
        h1 = int(self.left.get_height() * 20 / 2)
        h2 = int(self.right.get_height() * 20 / 2)
        top = y - (h1 + h2)
        bottom = y + (h1 + h2)
        # 子节点垂直线
        draw.line((x, top + h1, x, bottom - h2), fill=(0, 0, 0))
        # 水平线
        ll = self.distance * s
        draw.line((x, top + h1, x + ll, top + h1), fill=(0, 0, 0))
        draw.line((x, bottom - h2, x + ll, bottom - h2), fill=(0, 0, 0))
        # 递归地画左边和右边的子节点
        self.left.draw(draw, x + ll, top + h1, s, imlist, im)
        self.right.draw(draw, x + ll, bottom - h2, s, imlist, im)

#在画实际图像缩略图时,叶节点有自己的方法

    def draw(self, draw, x, y, s, imlist, im):
        nodeim = Image.open(imlist[self.id])
        nodeim.thumbnail([20, 20])
        ns = nodeim.size
        im.paste(nodeim, [int(x), int(y - ns[1] // 2), int(x + ns[0]), int(y + ns[1] - ns[1] // 2)])

100 幅日落图像进行层次聚类的示例聚类簇,阈值集合设定为树中最大节点距离的 23%

三、谱聚类

谱聚类方法是一种有趣的聚类算法。

对于n个元素(如n幅图像),相似矩阵(或亲和矩阵,有时也称距离矩阵)是一个n×n的矩阵,矩阵每个元素表示两两之间的相似性分数。谱聚类是由相似性矩阵构建谱矩阵而得名的。对该谱矩阵进行特征分解得到的特征向量可以用于降维,然后聚类。

下面说明谱聚类的过程。给定一个n×n的相似矩阵S,S_{ij}为相似性分数,我们可以创建一个矩阵,称为拉普拉斯矩阵:

\mathbf{L=I-D^{-1/2}SD^{-1/2}}

拉普拉斯矩阵中的\mathbf{D^{-1/2}}为:


import pickle
from scipy.cluster.vq import *

from PIL import Image
from pylab import *


# 获取 selected-fontimages 文件下图像文件名,并保存在列表中
from py_cvp import imtools


imlist = imtools.get_imlist('E:/data/a_selected_thumbs')
imnbr = len(imlist)
# 载入模型文件
with open('../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)])
projected = array([dot(V[[1,2]],immatrix[i]-immean) for i in range(imnbr)])
# 进行 k-means 聚类
projected = whiten(projected)
centroids,distortion = kmeans(projected,4)
code,distance = vq(projected,centroids)



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()
 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()

 

 

 

 

 

  • 10
    点赞
  • 8
    收藏
    觉得还不错? 一键收藏
  • 0
    评论
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值