数据挖掘实验指导书
聚类分析
实验目的
1. 理解聚类算法的基本原理。
2. 根据实际应用，灵活运用聚类算法，比较分析不同聚类方法的结果。
实验数据
现有若干视频数据，每个视频由图像序列It(t=1..N)构成，每帧图像和下一帧（或上一帧）之间的运动可以通过光流估计算法得到运动图像Mt（运动图像中有的是错误的）,通过聚类算法从图像中分割出目标，并将其和真实分割中目标（白色区域）比较，计算全部视频数据的平均Recall，Precision，F1度量。
图像I1 图像I2 运动M1 真实分割
实验要求
- 运用的聚类算法中，至少有一种聚类算法是自己独立编程实现的，不能直接调用第三方库。其他聚类算法可以调用第三方库。
- 比较不同的聚类算法，不同的属性选择，聚类结果差别并解释原因。
报告内容：
一·实验内容：
通过聚类算法从图像中分割出目标。
二·设计思想：主要算法的基本思想。
- 对运动图像依次进行canny边缘检测，核填充，对目标打上Box标签。
- 分别对运动图像的Box位置，非Box位置进行聚类（Box区域聚为5个簇，非Box区域聚为1个簇），并返回各个簇质心像素。
- 得出Box区域的5个簇质心像素与非Box区域聚的1个簇质心像素之间的距离并从小到大进行排列，选取中间三个像素的平均像素作为目标像素，非Box区域质心像素之作为环境像素，分别计算Box区域的5个簇质心与目标像素，环境像素之间的距离，拿Box区域中的一个簇质心举例：令其到目标像素，环境像素的距离为a,b，若a/b<1,则表明该簇属于目标物体，否则属于背景环境。
- 新建空白灰度图像并向其添加灰度值，目标物体添加256，背景环境添加0。
三·主要算法的程序实现：
#对运动图像依次进行canny边缘检测，核填充，对目标打上Box标签。
edge_img = cv2.Canny(img, 0.5, 100,apertureSize=3)
# img = Image.fromarray(edge_img)
kernel = cv2.getStructuringElement(cv2.MORPH_ELLIPSE, (95, 70))
closed = cv2.morphologyEx(edge_img, cv2.MORPH_CLOSE, kernel)
# closed = cv2.erode(closed, None, iterations=16)
# closed = cv2.dilate(closed, None, iterations=16)
cnts, _ = cv2.findContours(closed.copy(),
cv2.RETR_LIST,
cv2.CHAIN_APPROX_SIMPLE)
c = sorted(cnts, key=cv2.contourArea, reverse=True)[0]
# compute the rotated bounding box of the largest contour
rect = cv2.minAreaRect(c)
box = np.int0(cv2.boxPoints(rect))
draw_img = cv2.drawContours(original_img.copy(), [box], -1, (0, 0, 255), 3)
return draw_img,crop_img,y1,y3,x1,x3 #容纳Box的正矩形在原图中的位置;左下，右上坐标。
#分别对运动图像的Box位置，非Box位置进行聚类（Box区域聚为5个簇，非Box区域聚为1个簇），判断簇的类别（目标还是环境）并分别向其中填充灰度值。
km = KMeans(init="k-means++",n_clusters=5) # 聚类中心个数为5
km2 = KMeans(init="k-means++", n_clusters=1)
f1 = km.fit(imgData)
f2 = km2.fit(imgData2)
f1center = f1.cluster_centers_
f2center = f2.cluster_centers_
cendis = []
for i in f1center:
cendis.append(distance(i,f2center[0]))
#获取数组中的N个最大值的索引
arr = np.array(cendis)
idx1 = arr.argsort()[2:4]
c = [0,0,0]
for i in idx1:
for j in range(len(f1center[i])):
c[j] += f1center[i][j]
mdis = [x/len(idx1) for x in c]
# maxdis,idx1 = max(cendis),cendis.index(max(cendis))
F = []
E = []
for i in range(len(f1center)):
a = distance(f1center[i],mdis)
b = distance(f1center[i],f2center[0])
if a/b<=1:
F.append(i)
else:
E.append(i)
# print(F)
label = km.predict(imgData)
print(f1center)
print(km.predict(f1center))
label = label.reshape([x2-x1+1,y2-y1+1]) # 这里用的numpy里面的一个新数组，row,col代表行与列
# print(label)
# 创建一张新的灰度图保存聚类后结果
pic_new = Image.new("L", (row, col))
# 根据所属类别向图片中添加灰度值
for i in range(row):
for j in range(col):
if i>=x1 and i<=x2 and j>=y1 and j<=y2:
if label[i-x1][j-y1] in F:
pic_new.putpixel((i, j), 256)
else:
pic_new.putpixel((i, j), 0)
else:
pic_new.putpixel((i, j), 0)
return pic_new
注：标黄的代码参数可调节：
km = KMeans(init="k-means++",n_clusters=5) # 聚类中心个数为5
km2 = KMeans(init="k-means++", n_clusters=1)
km1可再适当增大，将目标聚成更多的簇，效果还会好一些，但是会增加运行时间。
对于km2，也可再适当增大，假如增加到2，需要对每个目标簇中心计算与环境簇中心（2个）的最小距离并排序，得出来目标像素值（求平均）之后，再计算每个目标簇中心与目标像素值，每个环境簇中心的最小距离，然后判断每个目标簇的类别。环境簇不能太多，会增加目标簇判别的随机性。这是基于目标簇不一定是检测目标，但环境簇一定是背景环境的假设，该假设成立的关键在于Box的划分效果。在该实验中，Box的划分效果还是很不错的，在实际的测试中用了每种物体的前50张图像，共350张图像，Box划分失误的图像大概有4-5张。每类50张之后的Box划分会增加失误率，因为图像边缘不明显，canny效果不好，这时候可以再在真实图像上想办法，对真实图像再进行一些处理使得Box标记准确。
edge_img = cv2.Canny(img, 0.5, 100,apertureSize=3)
Canny对梯度大小的设置可根据具体图像轮廓明显度和背景嘈杂度进行更改，
图像轮廓越明显，背景嘈杂度越大，其梯度下界越高，反之越低。
kernel = cv2.getStructuringElement(cv2.MORPH_ELLIPSE, (95, 70))
kernel是用来填充canny的填充核，设置的太大会导致边缘模糊，太小又会导致填充不完全。在breakdance-flare图像部分表现效果比较好一点（因为其填充形状是细长形的）。效果好坏关键在于设置这个填充核的大小。识别物体越近似球形，面积越大，效果越差，反之越好。
idx1 = arr.argsort()[1:4]
这个是对目标簇中心与环境簇中心的距离进行从小到大排序并返回对应的索引值。所以调节时左边参数（上面是1），右边参数（上面是4）越小，平均像素值（目标像素值）越接近环境像素值，两边参数越大，则越远离环境像素值，换句话说就是越接近目标像素（这句话不严谨，但可以这样理解记忆，方便后面调节）。所以在调节时，若检测物体效果残缺，则可以调小两边参数，若检测物体效果包含环境部分，则把两边参数调小。这个调节我称之为粗调。
a/b<=1
这个调节我称之为细调，增大参数，即增加了目标簇成为目标的可能性，减小参数，即减小了目标簇成为目标的可能性。
在做本次实验时，我用的参数是
idx1 = arr.argsort()[1:4]
a/b<=1.2
其中有些图像进行了微调。
理论上说，对于一些效果不好的图像，可以调整参数进行调节，总能得出一个比较好的效果。
四·调试报告：调试过程中遇到的主要问题是如何解决的；对设计和编码的回顾讨论和分析以及对主要算法的改进设想。
1.调试报告：
在代码实现过程中，一开始是使用普通的聚类，聚成两类，但效果不好，于是便尝试聚成3，4，5.....类，因为不知道哪些簇属于物体，所以在向空白灰度图写入像素时，出现非常不好的效果，甚至不如分成两类，于是想如果能识别出来哪块是物体，哪块是环境，那么分别对环境，物体进行聚类应该会得到一个比较好的效果，但怎么样才能知道哪一块是物体，哪一块是环境呢？我当时没有想法，所以放弃了聚类，转向了边缘检测以及图像形态学进行分割图像，第一次是接触了一些算子和canny，然后通过边缘检测实现图像分割的硬币例子知道了分水岭注水，但它是应用在封闭的图形，结果效果不好，于是用其他的填充方法：cv2.getStructuringElement()，cv2.morphologyEx()，效果比分水岭强一些，但效果依然不理想，填充核kernel的调试并不能取得令人满意的效果，后来发现了Box标记：cv2.boxPoints(rect)，感觉一下子有想法了，因为这样能够知道大致了解哪块是物体，哪块是环境，所以有了后来的思路。
2.代码改进
1.由于时间有限，所以我只对部分效果不好的图像进行了重新调整参数（主要用到了粗调和细调），有时间可以对所有不好的图像进行重新调参。
2.环境簇，物体簇都可以聚成更多簇。
3.对于Box标记不准的问题，一般这种问题出现在canny梯度下界过高或者运动的图像质量太差（肉眼几乎分辨不出来物体与环境边界的哪种）。可降低梯度下界，设置Box的面积临界值，若小于临界值则对整张图像进行聚类，或者以Box标记的一小部分为中心设置一个面积合适的矩形作为Box,或者通过对原图再进行一些处理，用canny重新进行检测边缘。我认为最好的办法还是得引入深度学习，训练相应的卷积神经网络直接检测出物体，或者检测出物体的Box，然后再进行处理。
"""
coding:utf-8
@author: Li Sentan
@time:2021.12.12
@file:clustersys02.py
"""
import matplotlib.pyplot as plt
from scipy import ndimage as ndi
from skimage import morphology,color,filters
import cv2
from PIL import Image
import numpy as np
from sklearn.cluster import KMeans
import math
np.set_printoptions(threshold=np.inf) # np.inf表示正无穷
def distance(a,b): #计算两个数组像素之间的距离
distance = 0
for i in range(len(a)):
distance += (a[i]-b[i])**2
dis = math.sqrt(distance)
return dis
def loadData(filePath,y1,y2,x1,x2):
f = open(filePath,'rb')#以二进制打开
img = Image.open(f)#以列表形式返回图片像素值
data = []
data1 = []
m,n = img.size#获得图片大小
print(m,n)
for i in range(m):#将每一个像素点RGB颜色处理到0-1
for j in range(n):#范围内并存进data
x, y, z = img.getpixel((i, j))
if i >= x1 and i <= x2 and j >= y1 and j <= y2:
# data.append([x,y,z])
data.append([x/256.0,y/256.0,z/256.0])
else:
data1.append([x/256.0,y/256.0,z/256.0])
f.close()
return np.mat(data),np.mat(data1),m,n#以矩阵形式返回data，以及图片大小
def clustersys(imgpath,km,km2,y1,y2,x1,x2):
imgData,imgData2, row, col = loadData(imgpath,y1,y2,x1,x2) # 加载数据
# img = cv2.imread(imgpath,0)
# 聚类获得每一个像素所属的类别
f1 = km.fit(imgData)
f2 = km2.fit(imgData2)
f1center = f1.cluster_centers_
f2center = f2.cluster_centers_
cendis = []
for i in f1center:
cendis.append(distance(i,f2center[0]))
## 获取numn数组中的N个最大值的索引
arr = np.array(cendis)
idx1 = arr.argsort()[2:4]
c = [0,0,0]
for i in idx1:
for j in range(len(f1center[i])):
c[j] += f1center[i][j]
mdis = [x/len(idx1) for x in c]
# maxdis,idx1 = max(cendis),cendis.index(max(cendis))
F = []
E = []
for i in range(len(f1center)):
a = distance(f1center[i],mdis)
b = distance(f1center[i],f2center[0])
if a/b<=1:
F.append(i)
else:
E.append(i)
# print(F)
label = km.predict(imgData)
print(f1center)
print(km.predict(f1center))
label = label.reshape([x2-x1+1,y2-y1+1]) # 这里用的numpy里面的一个新数组，row,col代表行与列
# print(label)
# 创建一张新的灰度图保存聚类后结果
pic_new = Image.new("L", (row, col))
# 根据所属类别向图片中添加灰度值
for i in range(row):
for j in range(col):
if i>=x1 and i<=x2 and j>=y1 and j<=y2:
if label[i-x1][j-y1] in F:
pic_new.putpixel((i, j), 256)
else:
pic_new.putpixel((i, j), 0)
else:
pic_new.putpixel((i, j), 0)
return pic_new
def drawcnts_and_cut(original_img, box):
# 因为这个函数有极强的破坏性，所有需要在img.copy()上画
# draw a bounding box arounded the detected barcode and display the image
draw_img = cv2.drawContours(original_img.copy(), [box], -1, (0, 0, 255), 3)
Xs = [i[0] for i in box]
Ys = [i[1] for i in box]
x1 = min(Xs)
x2 = max(Xs)
y1 = min(Ys)
y2 = max(Ys)
hight = y2 - y1
width = x2 - x1
y3 = y1+hight
x3 = x1+width
if y1<0:
y1 = 0
if x1<0:
x1 = 0
if y3>=480:
y3 = 439
if x3 >= 854:
x3 = 853
crop_img = original_img[y1:y1 + hight, x1:x1 + width]
return draw_img,crop_img,y1,y3,x1,x3
def cv2_canny(img):
'''
由于Canny只能处理灰度图，所以将读取的图像转成灰度图。
用高斯平滑处理原图像降噪。
调用Canny函数，指定最大和最小阈值，其中apertureSize默认为3。
:param img:
:return:
'''
img = cv2.imread(img, 0)
img = cv2.GaussianBlur(img, (3, 3), 2)
edge_img = cv2.Canny(img, 0.5, 100,apertureSize=3)
# img = Image.fromarray(edge_img)
kernel = cv2.getStructuringElement(cv2.MORPH_ELLIPSE, (95, 70))
closed = cv2.morphologyEx(edge_img, cv2.MORPH_CLOSE, kernel)
# closed = cv2.erode(closed, None, iterations=16)
# closed = cv2.dilate(closed, None, iterations=16)
cnts, _ = cv2.findContours(closed.copy(),
cv2.RETR_LIST,
cv2.CHAIN_APPROX_SIMPLE)
c = sorted(cnts, key=cv2.contourArea, reverse=True)[0]
# compute the rotated bounding box of the largest contour
rect = cv2.minAreaRect(c)
box = np.int0(cv2.boxPoints(rect))
return img,edge_img,closed,box
if __name__ == '__main__':
km = KMeans(init="k-means++",n_clusters=5) # 聚类中心个数为3
km2 = KMeans(init="k-means++", n_clusters=1)
for w in range(17,18):
print("The {}th picture:".format(w))
imgpath = 'cluster data/flowsports/bear/000{}.jpg'.format("%02d" % w)
imgpath2 = 'cluster data/truepictures/bear/000{}.jpg'.format("%02d" % w)
img,edge_img,closed,box = cv2_canny(imgpath)
original_img = cv2.imread(imgpath2)
img3,img4,y1, y2, x1, x2 = drawcnts_and_cut(original_img,box)
print(y1,y2,x1,x2)
pic_new = clustersys(imgpath,km,km2,y1,y2,x1,x2)
fig, axes = plt.subplots(nrows=3, ncols=2, figsize=(6, 6))
axes = axes.ravel()
ax0, ax1, ax2, ax3,ax4,ax5 = axes
ax0.imshow(img, cmap=plt.cm.gray, interpolation='nearest')
ax0.set_title("Original")
ax1.imshow(edge_img, cmap=plt.cm.gray, interpolation='nearest')
ax1.set_title("canny")
ax2.imshow(closed, cmap=plt.cm.gray, interpolation='nearest')
ax2.set_title("closed")
ax3.imshow(img3, cmap=plt.cm.gray, interpolation='nearest')
ax3.set_title("box")
ax4.imshow(img4, cmap=plt.cm.gray, interpolation='nearest')
ax4.set_title("cut")
ax5.imshow(pic_new, cmap=plt.cm.gray, interpolation='nearest')
ax5.set_title("new")
for ax in axes:
ax.axis('off')
fig.tight_layout()
plt.savefig("F:/Desktop/clusters/bear/canny-{}.jpg".format("%02d" % w))
plt.show()
939
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
