# 粘贴每幅图像的缩略图到白色背景图片
for i in range(imnbr):
nodeim = Image.open(imlist[i])
nodeim.thumbnail((25,25))
ns = nodeim.size
img.paste(nodeim,(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.save('pca_font.jpg')
图 6-3：在成对主成分上投影的字体图像。左图用的是第一个和第二个主成分，右图用的是第二个和第三个主成分
这里，我们用到了整数或 floor 向下取整除法运算符 //，通过移去小数点后面的部分，可以返回各个缩略图在白色背景中对应的整数坐标位置。
这类图像说明这些字体图像在 40 维里的分布情况，对于选择一个好的描述子很有帮助。可以很清楚地看到，二维投影后相似的字体图像距离较近。
6.1.4 像素聚类
在结束本节之前，我们来看一个对单幅图像中的像素而非全部图像进行聚类的例子。将图像区域或像素合并成有意义的部分称为图像分割，它是第 9 章的主题。除了在一些简单的图像上，单纯在像素水平上应用 K-means 得出的结果往往是毫无意义的。要产生有意义的结果，往往需要更复杂的类模型而非平均像素色彩或空间一致性。现在，我们仅会在 RGB 三通道的像素值上运用 K-means 进行聚类，分割问题的处理方法会在之后谈到（9.2 节）给予关注及给出细节部分。下面的代码示例载入一幅图像，用一个步长为 steps 的方形网格在图像中滑动，每滑一次对网格中图像区域像素求平均值，将其作为新生成的低分辨率图像对应位置处的像素值，并用 K-means 进行聚类：
from scipy.cluster.vq import *
from scipy.misc import imresize
steps = 50 # 图像被划分成 steps×steps 的区域
im = array(Image.open('empire.jpg'))
dx = im.shape[0] / steps
dy = 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')
figure()
imshow(codeim)
show()
K-means 的输入是一个有 steps×steps 行的数组，数组的每一行有 3 列，各列分别为区域块 R、G、B 三个通道的像素平均值。为可视化最后的结果 , 我们用 SciPy 的imresize() 函数在原图像坐标中显示这幅 steps×steps 的图像。参数 interp 指定插值方法；我们在这里采用最近邻插值法，以便在类间进行变换时不需要引入新的像素值。图 6-4 显示了用 50×50 和 100×100 窗口对两幅相对比较简单的示例图像进行像素聚类后的结果。注意，K-means 标签的次序是任意的（在这里的标签指最终结果中图像的颜色）。正如你所看到的，尽管利用窗口对它进行了下采样，但结果仍然是有噪声的。如果图像某些区域没有空间一致性，则很难将它们分开，如下方图中小男孩和草坪的图。空间一致性和更好的分割方法以及其他的图像分割算法会在后面讨论。现在，让我们继续来看下一个基本的聚类算法。
图 6-4：基于颜色像素值用 K-means 对像素进行聚类的结果。左边是原始图像；中间是用 k=3和 50×50 大小的窗口进行聚类的结果；右边是用 k=3 和 100×100 大小的窗口进行聚类的结果
6.2 层次聚类
层次聚类（或凝聚式聚类）是另一种简单但有效的聚类算法，其思想是基于样本间成对距离建立一个简相似性树。该算法首先将特征向量距离最近的两个样本归并为一组，并在树中创建一个“平均”节点，将这两个距离最近的样本作为该“平均”节点下的子节点；然后在剩下的包含任意平均节点的样本中寻找下一个最近的对，重复进行前面的操作。在每一个节点处保存了两个子节点之间的距离。遍历整个树，通过设定的阈值，遍历过程可以在比阈值大的节点位置终止，从而提取出聚类簇。层次聚类有若干优点。例如，利用树结构可以可视化数据间的关系，并显示这些簇是如何关联的。在树中，一个好的特征向量可以给出一个很好的分离结果。另外一个优点是，对于给定的不同的阈值，可以直接利用原来的树，而不需要重新计算。不足之处是，对于实际需要的聚类簇，我们需要选择一个合适的阈值。让我们看看层次聚类算法怎样在代码中体现 1 。创建文件 hcluster.py，将下面代码添加进去：
from itertools import combinations
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]
我们为树节点创建了两个类，即 ClusterNode 和 ClusterLeafNode，这两个类将用于创建聚类树，其中函数 hcluster() 用于创建树。首先创建一个包含叶节点的列表，然后根据选择的距离度量方式将距离最近的对归并到一起，返回的终节点即为树的根。对于一个行为特征向量的矩阵，运行 hcluster() 会创建和返回聚类树。距离度量的选择依赖于实际的特征向量，这里我们利用欧式距离 L2（同时提供了L1 距离度量函数），不过你可以创建任意距离度量函数，并将它作为参数传递给hcluster()。对于每个子树，计算其所有节点特征向量的平均值，作为新的特征向量来表示该子树，并将每个子树视为一个对象。当然，还有其他将哪两个节点合并在一起的方案，比如在两个子树中使用对象间距离最小的单向锁，及在两个子树中用对象间距离最大的完全锁。选择不同的锁会生成不同类型的聚类树。为了从树中提取聚类簇，需要从顶部遍历树直至一个距离小于设定阈值的节点终止，这通过递归很容易做到。ClusterNode 的 extract_clusters() 方法用于处理该过程，如果节点间距离小于阈值，则用一个列表返回节点，否则调用子节点（叶节点通常返回它们自身）。调用该函数会返回一个包含聚类簇的子树列表。对于每一个子聚类簇，为了得到包含对象 id 的叶节点，需要遍历每个子树，并用方法 get_cluster_elements() 返回一个包含叶节点的列表。下面，我们在一个简单的例子中观察该聚类过程。首先创建一些二维数据点（和之前 K-means 一样）：
class1 = 1.5 * randn(100,2)
class2 = randn(100,2) + array([5,5])
features = vstack((class1,class2))
对这些数据点进行聚类，设定阈值（这里的阈值设定为 5），从列表中提取这些聚类簇，并在控制台打印出来：
import hcluster
tree = hcluster.hcluster(features)
clusters = tree.extract_clusters(5)
print 'number of clusters', len(clusters)
for c in clusters:
print c.get_cluster_elements()
打印结果应该与下面类似：
number of clusters 2
[184, 187, 196, 137, 174, 102, 147, 145, 185, 109, 166, 152, 173, 180, 128, 163, 141,
178, 151, 158, 108,182, 112, 199, 100, 119, 132, 195, 105, 159, 140, 171, 191, 164, 130,
149, 150, 157, 176, 135, 123, 131,118, 170, 143, 125, 127, 139, 179, 126, 160, 162, 114,
122, 103, 146, 115, 120, 142, 111, 154, 116, 129,136, 144, 167, 106, 107, 198, 186, 153,
156, 134, 101, 110, 133, 189, 168, 183, 148, 165, 172, 188, 138,192, 104, 124, 113, 194,
190, 161, 175, 121, 197, 177, 193, 169, 117, 155]
[56, 4, 47, 18, 51, 95, 29, 91, 23, 80, 83, 3, 54, 68, 69, 5, 21, 1, 44, 57, 17, 90, 30,
22, 63, 41, 7, 14, 59, 96, 20, 26, 71, 88, 86, 40, 27, 38, 50, 55, 67, 8, 28, 79, 64,
66, 94, 33, 53, 70, 31, 81, 9, 75, 15, 32, 89, 6, 11, 48, 58, 2, 39, 61, 45, 65, 82, 93,
97, 52, 62, 16, 43, 84, 24, 19, 74, 36, 37, 60, 87, 92, 181, 99, 10, 49, 12, 76, 98, 46,
72, 34, 35, 13, 73, 78, 25, 42, 77, 85]
理想情况下，你应该得到两一个聚类簇，但是在实际数据中，你可能会得到三类或更多这主要依赖于实际生成的二维数据。在这个对二维数据聚类的简单例子中，一个类中的值应该小于 100，另外一个应该大于等于 100。
图像聚类
我们来看一个基于图像颜色信息对图像进行聚类的例子。文件 sunsets.zip 中包含 100 张图像，这些图像是用“sunset”和“sunsets”关键字在 Flickr 下载下来的。在这个例子中，我们用颜色直方图作为每幅图像的特征向量。虽然这样处理有些简单粗糙，但仍然能够很好地说明分层聚类的过程。在包含这些日落图像的文件夹中运行下面的代码：
import os
import hcluster
# 创建图像列表
path = 'flickr-sunsets/'
imlist = [os.path.join(path,f) for f in os.listdir(path) if f.endswith('.jpg')]
# 提取特征向量，每个颜色通道量化成 8 个小区间
features = zeros([len(imlist), 512])
for i,f in enumerate(imlist):
im = array(Image.open(f))
# 多维直方图
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)
我们将 R、G、B 三个颜色通道作为特征向量，将其传递到 NumPy 的 histogramdd()中，该函数能够计算多维直方图（本例中是三维）。我们在每个颜色通道中使用 8 个小区间进行量化，将三个通道量化后的小区间拉成一行后便可用 512（8×8×8）维的特征向量描述每幅图像。为避免图像尺寸不一致，我们用“normed=True”归一化直方图，并将每个颜色通道范围设置为 0...255。将 reshape() 第一个参数设置为 -1会自动确定正确的尺寸，故可以创建一个输入数组来计算以 RGB 颜色值为行向量的直方图。为了可视化聚类树，我们可以画出树状图。树状图是一种显示树布局的图表。在判定给出的描述子向量好坏，以及在特征场合考虑什么是相似的时候，树状图可以提供有用的信息。将下面的代码添加到 hcluster.py 中：
from PIL import Image,ImageDraw
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() 方法，将该方法添加到 ClusterNode 类中：
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)
在画实际图像缩略图时，叶节点有自己的方法，将该方法添加到 ClusterLeafNode类中：
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)])
树状图的高和子部分由距离决定，这些都需要调整，以适应所选择的图像分辨率。随着坐标向下传递到下一级，会递归绘制出这些节点，上述代码用 20×20 像素绘制叶节点的缩略图，使用 get_height() 和 get_depth() 这两个辅助函数可以获得树的高和宽。
树状图可以通过下面的代码绘制，并保存在 sunset.pdf 中：
hcluster.draw_dendrogram(tree,imlist,filename='sunset.pdf')
图 6-5 展示了日落图像聚类后的树状图。可以看到，树中颜色相似的图像距离较近。
图 6-6 中为三个示例簇。可以通过下面的代码提取该例子中的簇：
# 设置一些（任意的）阈值以可视化聚类簇
clusters = tree.extract_clusters(0.23*tree.distance)
# 绘制聚类簇中元素超过 3 个的那些图像
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()
图 6-5：用 100 幅日落图像进行层次聚类，将 RGB 空间的 512 个小区间直方图作为每幅图像的特征向量。树中挨的相近的图像具有相似的颜色分布
图 6-6：用 100 幅日落图像进行层次聚类的示例聚类簇，阈值集合设定为树中最大节点距离的 23%
作为最后一个例子，我们对前面的字体图像创建一个树状图：
tree = hcluster.hcluster(projected)
hcluster.draw_dendrogram(tree,imlist,filename='fonts.jpg')
其中，projected 和 imlist 是 6.1 节 K-means 例子中的变量。图 6-7 显示了对字体图像进行层次聚类后的树状图。
6.3 谱聚类
谱聚类方法是一种有趣的聚类算法，与前面 K-means 和层次聚类方法截然不同。对于 n 个元素（如 n 幅图像），相似矩阵（或亲和矩阵，有时也称距离矩阵）是一个n×n 的矩阵，矩阵每个元素表示两两之间的相似性分数。谱聚类是由相似性矩阵构建谱矩阵而得名的。对该谱矩阵进行特征分解得到的特征向量可以用于降维，然后聚类。
谱聚类的优点之一是仅需输入相似性矩阵，并且可以采用你所想到的任意度量方式构建该相似性矩阵。像 K-means 和层次聚类需要计算那些特征向量求平均；为了计算平均值，会将特征或描述子限制为向量。而对于谱方法，特征向量就没有类别限制，只要有一个“距离”或“相似性”的概念即可。
6-7：用 66 幅选定字体图像，使用 40 个主成分作为特征量，用层次聚类方法进行聚类
计算 L 的特征向量，并使用 k 个最大特征值对应的 k 个特征向量，构建出一个特征向量集（记住，我们可能并没有以任何东西来开始），从而可以找到聚类簇。创建一个矩阵，该矩阵的各列是由之前求出的 k 个特征向量构成，每一行可以看做一个新的特征向量，长度为 k。这些新的特征向量可以用诸如 K-means 方法进行聚类，生成最终的聚类簇。本质上，谱聚类算法是将原始空间中的数据转换成更容易聚类的新特征向量。在某些情况下，不会首先使用聚类算法。
讲解了足够的理论，我们来看看真实的例子中谱聚类算法的代码。我们再次使用1.3.6 节 K-means 例子中的字体图像。
from scipy.cluster.vq import *
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(path+imlist[ind[i]])
subplot(4,10,i+1)
imshow(array(im))
axis('equal')
axis('off')
show()
在本例中，我们用两两间的欧式距离创建矩阵 S，并对 k 个特征向量（eignvector）用常规的 K-means 进行聚类（在该例中，k=5）。注意，矩阵 V 包含的是对特征值进行排序后的特征向量。最后，绘制出这些聚类簇。图 6-8 显示了运行后的聚类簇；需要记住的是，在 K-means 阶段，每次运行的结果可能不同。
图 6-8：用拉普拉斯矩阵的特征向量对字体图像进行谱聚类
我们也可以在没有任何特征向量或没有严格定义相似性的例子中尝试该算法。2.3 节中带有地理标签的 Panoramio 图像是基于它们之间有多少匹配的局部描述符连接起来的。2.3.2 节的矩阵是一个用分数表示的相似性矩阵，其中分数等于匹配的特征数（没有归一化）。由于 imlist 列表包含了图像文件名，并已用 NumPy 的 savetxt() 将相似性矩阵保存到了文件中，所以我们只需修改上面代码的前面几行：
n = len(imlist)
# 载入相似矩阵并重新格式化
S = loadtxt('panoramio_matches.txt')
S = 1 / (S + 1e-6)
这里对分数进行转换，使得相似图像的分数值较小，这样我们就不需要修改上面的代码。我们添加了一个很小的数以防止与 0 相除，后面的代码不需要修改。在该例中选择 k 有些技巧。很多人会认为这里只有两类（即白宫的两侧），以及其他一些垃圾图像。用 k=2 可以得到类似图 6-9 的结果，其中一个聚类簇是包含很多白宫一侧的图像，另一个聚类簇是白宫另一侧的图像和其他所有垃圾图像。将 k 设定为一个较大的值，比如 k=10，则有些聚类簇可能只包含一幅图像（很可能是垃圾图像），另一些是真实的聚类簇。图 6-10 给出了上面示例代码运行的结果。在该例中，仅有两个真实的聚类簇，每个聚类簇包含白宫一个侧面的图像。
图 6-9：用 k=2、局部特征匹配数作为相似性分数对白宫地理图像进行谱聚类的结果
标签：特征向量,self,聚类,im,图像,节点 From： https://blog.csdn.net/qq_45632142/article/details/139919596