分级聚类算法（集体智慧编程）

分级聚类是通过连续不断地将最为相似的群组两两合并，来构造出一个群组的层级结构。其中的每个群组都是从单一元素开始的。如图所示：

元素的相似程序是通过它们的相对位置来体现的，距离越近越相似。两两合并，直到合并最后两个群组。

 

聚类是无监督学习的一个例子。与神经网络或决策树不同，无监督学习算法不是利用带有正确答案的样本数据进行“训练”。它们的目的是要在一组数据中找寻某种结构，而这些数据本身不是我们要找的答案。聚类算法的目标是采集数据，然后从中找出不同的群组。

 

下面我就利用分级聚类算法对iteye和csdn的博客进行分类。

 

1.准备数据。利用feedparser解析博客rss订阅源，利用jieba中文拆词获取关键词及其出现次数。

2.合并数据。利用皮尔逊相关度算法计算数据的相似度，两两合并相似度最小的数组，直到只剩下一个数级为止。

3.绘制树状图。利用PIL将合并后的数据做成图片，以此清晰地展现数据结构。

 

下面开始第一步：

安装feedparser.

feedparser是一个python库，可以用以分析RSS和Atom的订阅源。下载地址https://code.google.com/p/feedparser/downloads/list。解压后在feedparser目录下命令行输入python setup.py install。

下面去寻找一些iteye和csdn博客的订阅源。下面是我随便找的10个，写入到feedlist.txt中。

 

Python代码  

1. <span>http://blog.csdn.net/xxg2810/rss/list  
2. http://lobert.iteye.com/rss  
3. http://blog.csdn.net/lfsf802/rss/list  
4. http://blog.csdn.net/david_lv/rss/list  
5. http://blog.csdn.net/m13666368773/rss/list  
6. http://blog.csdn.net/shulianghan/rss/list  
7. http://blog.csdn.net/withiter/rss/list  
8. http://blog.csdn.net/lfsf802/rss/list  
9. http://blog.csdn.net/pan_tian/rss/list  
10. http://blog.csdn.net/beijiguangyong/rss/list</span>  

 

 

因为我们在天朝，所以用的当然是中文，因此需要一个分词库，这里我使用jieba(开始使用的是mmseg-cpp，但在实际运行中会出现卡死的现象)。

下载地址：https://github.com/fxsjy/jieba，解压后将文件夹直接放到python默认库地址即可。

准备工作完毕，开始mark了。

 

Python代码  

1. #encoding=utf-8           #注意编码  
2. import sys  
3. from imp import reload
   
4. reload(sys)   
5. import feedparser  
6. import re  
7. import jieba  
8.   
9. # Returns title and dictionary of word counts for an RSS feed  
10. def getwordcounts(url):  
11.   # Parse the feed  
12.   d=feedparser.parse(url)          
13.   wc={}  
14.   
15.   # Loop over all the entries  
16.   for e in d.entries:  
17.     if 'summary' in e: summary=e.summary  
18.     else: summary=e.description  
19.   
20.     # Extract a list of words  
21.     words=getwords(e.title+' '+summary)  
22.     for word in words:  
23.       wc.setdefault(word,0)  
24.       wc[word]+=1  
25.   return d.feed.title,wc  
26.   
27.   
28. def getwords(html):  
29.   # Remove all the HTML tags  
30.   txt=re.compile(r'<[^>]+>').sub('',html)  
31.   algor = jieba.cut(txt,cut_all=True)   
32.   
33.   # Split words by all non-alpha characters  
34.   # words=re.compile(r'[^A-Z^a-z]+').split(txt)  
35.   
36.   # Convert to lowercase  
37.   return [tok.lower() for tok in algor if tok!='']  
38.   # return [word.lower() for word in words if word!='']  
39.     
40.   
41. apcount={}  
42. wordcounts={}  
43. feedlist=[line for line in open('feedlist.txt')]  
44. for feedurl in feedlist:  
45.   try:  
46.     title,wc=getwordcounts(feedurl)  
47.     wordcounts[title]=wc  
48.     for word,count in wc.items():  
49.       apcount.setdefault(word,0)  
50.       if count>1:  
51.         apcount[word]+=1  
52.   except:  
53.     print 'Failed to parse feed %s' % feedurl  
54.   
55. wordlist=[]  
56. for w,bc in apcount.items():  
57.   frac=float(bc)/len(feedlist)  
58.   if frac>0.1 and frac<0.5:       //因为像“我”这样的单词几乎到处都是，而像“PYTHON”这样的单词则有可能只出现在个别博客中，所以通过只选择介于某个百分比范围内的单词，我们可以减少需要考查的单词总量。我们这里将1/10为下界，5/10为上界。  
59.     wordlist.append(w)  
60.   
61. out=open('blogdata1.txt','w')  
62. out.write('Blog')  
63. for word in wordlist: out.write('\t%s' % word.strip())  
64. out.write('\n')  
65. for blog,wc in wordcounts.items():  
66.   out.write(blog)  
67.   for word in wordlist:  
68.     if word in wc: out.write('\t%d' % wc[word])  
69.     else: out.write('\t0')  
70.   out.write('\n')  

 程序运行后数据保存到文件中，我这里跑了40s，如无意外，出现的是一个表格。其中的每一列对应一个单词，每一行对应一个博客。表格的数字就是单词出现在博客中的次数。

由于宽度原因，因此只截取部分数据。

到此，数据准备完毕。

 

 

第二步：mark

 

Python代码  

1. def readfile(filename):  
2.   lines=[line for line in file(filename)]  
3.     
4.   # First line is the column titles  
5.   colnames=lines[0].strip().split('\t')[1:]  
6.   rownames=[]  
7.   data=[]  
8.   for line in lines[1:]:  
9.     p=line.strip().split('\t')  
10.     # First column in each row is the rowname  
11.     rownames.append(p[0])  
12.     # The data for this row is the remainder of the row  
13.     data.append([float(x) for x in p[1:]])  
14.   return rownames,colnames,data  
15.   
16.   
17. from math import sqrt  
18.   
19. def pearson(v1,v2):  
20.   # Simple sums  
21.   sum1=sum(v1)  
22.   sum2=sum(v2)  
23.     
24.   # Sums of the squares  
25.   sum1Sq=sum([pow(v,2) for v in v1])  
26.   sum2Sq=sum([pow(v,2) for v in v2])    
27.     
28.   # Sum of the products  
29.   pSum=sum([v1[i]*v2[i] for i in range(len(v1))])  
30.     
31.   # Calculate r (Pearson score)  
32.   num=pSum-(sum1*sum2/len(v1))  
33.   den=sqrt((sum1Sq-pow(sum1,2)/len(v1))*(sum2Sq-pow(sum2,2)/len(v1)))  
34.   if den==0: return 0  
35.   
36.   return 1.0-num/den  
37.   
38. class bicluster:  
39.   def __init__(self,vec,left=None,right=None,distance=0.0,id=None):  
40.     self.left=left  
41.     self.right=right  
42.     self.vec=vec  
43.     self.id=id  
44.     self.distance=distance  
45.   
46. def hcluster(rows,distance=pearson):  
47.   distances={}  
48.   currentclustid=-1  
49.   
50.   # Clusters are initially just the rows  
51.   clust=[bicluster(rows[i],id=i) for i in range(len(rows))]  
52.   
53.   while len(clust)>1:  
54.     lowestpair=(0,1)  
55.     closest=distance(clust[0].vec,clust[1].vec)  
56.   
57.     # loop through every pair looking for the smallest distance  
58.     for i in range(len(clust)):  
59.       for j in range(i+1,len(clust)):  
60.         # distances is the cache of distance calculations  
61.         if (clust[i].id,clust[j].id) not in distances:   
62.           distances[(clust[i].id,clust[j].id)]=distance(clust[i].vec,clust[j].vec)  
63.   
64.         d=distances[(clust[i].id,clust[j].id)]  
65.   
66.         if d<closest:  
67.           closest=d  
68.           lowestpair=(i,j)  
69.   
70.     # 计算两个聚类的平均值  
71.     mergevec=[  
72.     (clust[lowestpair[0]].vec[i]+clust[lowestpair[1]].vec[i])/2.0   
73.     for i in range(len(clust[0].vec))]  
74.   
75.     # create the new cluster  
76.     newcluster=bicluster(mergevec,left=clust[lowestpair[0]],  
77.                          right=clust[lowestpair[1]],  
78.                          distance=closest,id=currentclustid)  
79.   
80.     # 不在原始集合中的聚类，其id为负数  
81.     currentclustid-=1  
82.     del clust[lowestpair[1]]  
83.     del clust[lowestpair[0]]  
84.     clust.append(newcluster)  
85.   
86.   return clust[0]  

 

 

第三步;借助树状图，我们可以更加理解聚类。

安装PIL，由于我安装了python(x,y)，里面已经包含了PIL，因此无需安装。

 

Python代码  

1. from PIL import Image,ImageDraw,ImageFont  
2. #计算聚类总体的高度  
3. def getheight(clust):  
4.   # Is this an endpoint? Then the height is just 1  
5.   if clust.left==None and clust.right==None: return 1  
6.   
7.   # Otherwise the height is the same of the heights of  
8.   # each branch  
9.   return getheight(clust.left)+getheight(clust.right)  
10. #计算误差  
11. def getdepth(clust):  
12.   # The distance of an endpoint is 0.0  
13.   if clust.left==None and clust.right==None: return 0  
14.   
15.   #一个枝节点的距离等于左右两侧分支中距离较大者加上该枝节点自身的距离  
16.   return max(getdepth(clust.left),getdepth(clust.right))+clust.distance  
17.   
18. #画图  
19. def drawdendrogram(clust,labels,jpeg='clusters.jpg'):  
20.   # height and width  
21.   h=getheight(clust)*30  
22.   w=1200  
23.   depth=getdepth(clust)  
24.   
25.   # width is fixed, so scale distances accordingly  
26.   scaling=float(w-350)/depth  
27.   
28.   # Create a new image with a white background  
29.   img=Image.new('RGB',(w,h),(255,255,255))  
30.   draw=ImageDraw.Draw(img)  
31.   
32.   draw.line((0,h/2,10,h/2),fill=(255,0,0))      
33.   
34.   # Draw the first node  
35.   drawnode(draw,clust,10,(h/2),scaling,labels)  
36.   img.save(jpeg,'JPEG')  
37.   
38. def drawnode(draw,clust,x,y,scaling,labels):  
39.   if clust.id<0:  
40.     h1=getheight(clust.left)*20  
41.     h2=getheight(clust.right)*20  
42.     top=y-(h1+h2)/2  
43.     bottom=y+(h1+h2)/2  
44.     # Line length  
45.     ll=clust.distance*scaling  
46.     # Vertical line from this cluster to children      
47.     draw.line((x,top+h1/2,x,bottom-h2/2),fill=(255,0,0))      
48.       
49.     # Horizontal line to left item  
50.     draw.line((x,top+h1/2,x+ll,top+h1/2),fill=(255,0,0))      
51.   
52.     # Horizontal line to right item  
53.     draw.line((x,bottom-h2/2,x+ll,bottom-h2/2),fill=(255,0,0))          
54.   
55.     # Call the function to draw the left and right nodes      
56.     drawnode(draw,clust.left,x+ll,top+h1/2,scaling,labels)  
57.     drawnode(draw,clust.right,x+ll,bottom-h2/2,scaling,labels)  
58.   else:     
59.     font = ImageFont.truetype('simsun.ttc',24)  
60.     # If this is an endpoint, draw the item label  
61.     draw.text((x+5,y-7),unicode(labels[clust.id],'utf-8'),(0,0,0),font=font)  

 

 

最后我们写个test跑下

 

Python代码  

1. #encoding=utf-8  
2. import clusters #里面包含了第二，三步的代码  
3. blognames,words,data = clusters.readfile('blogdata1.txt')  
4. clust = clusters.hcluster(data)  
5. clusters.printclust(clust,labels=blognames)  
6. clusters.drawdendrogram(clust,blognames,jpeg='blogclust.jpg')  

 

 

生成图片：