Kmeans算法找相似商品、FP树找频繁项集
第一题:用Kmeans处理数据集
数据集下载地址
https://archive.ics.uci.edu/ml/datasets/Sales_Transactions_Dataset_Weekly
运行环境
python3.7、PyCharm 2018.2.4 (Community Edition)
思路
根据所给数据集及其说明可以看出数据集共有812行,其中第一行是每一列的说明,每行数据共有107个数据,第1列数据是商品的代码,第2-53列数据是此商品在一年52周中每周卖出的数量,紧接着2列数据是分别是此种商品在这52周中的销量的最小值和最大值,后面的52列数据是根据最小值和最大值对52周的销量进行正则化处理后得到的52个0-1之间的浮点数。
目标是用K聚类的方法将这些商品分类,找出那些具有相似特征的商品,在一年中的某个时间段内销量尽量相似。所给论文中用SPSS将这811种商品的52个周销量数据合并为26个每两周销量数据,然后将这811种商品分成200簇。
仿照论文的做法,这里用python及numpy等工具包将这811种商品用Kmeans算法分为200簇,采用的是二分Kmeans聚类算法,即从1个簇逐步地选择一个离散度最高的、熵值越大的簇进行二分,直至分成200个簇,最后将簇分类的结果保存到文件中,然后利用matlibplot和scipy包对这些簇中有效的、离散度最低的、熵值最低的簇进行可视化,分析结果。
源代码
Kmeans.py
#-*- coding: utf-8 -*-
#Author:Yinli
import pickle
import numpy as np
#从文件中读取数据
#返回原始周销量数据和正则化后的周销量数据
def loadData(filename):
#逐行读取文件
fr = open(filename, 'r', encoding='utf-8')
arrayOfLines = fr.readlines()
#读取商品数目并初始化数据集
m = len(arrayOfLines)-1
dataSet = np.zeros((m,52))
normalizedDataSet = np.zeros((m,52))
#逐行处理数据
for i in range(1,len(arrayOfLines)):
#去掉句尾的换行符并将数据转为列表
arrayOfLines[i] = arrayOfLines[i].rstrip('\n')
currentLine = arrayOfLines[i].split(',')
#分出原始周销量数据和正则化周销量数据并写入数据集
currentLine0 = currentLine[1:53]
currentLine1 = currentLine[55:]
dataSet[i-1, :] = list(map(int, currentLine0))
normalizedDataSet[i-1,:] = list(map(float,currentLine1))
#返回数据集
return dataSet, normalizedDataSet
#计算向量A和向量B的欧式距离并返回
def distEclud(vecA, vecB):
return np.sqrt(np.sum(np.power(vecA - vecB, 2)))
#在数据集中随机选择k个质心
def randCent(dataSet, k):
#得到特征数目,初始化质心矩阵
n = np.shape(dataSet)[1]
centroids = np.mat(np.zeros((k, n)))
#遍历每个特征
for j in range(n):
#计算每个特征的最小值和最大值及范围
minJ = min(dataSet[:, j])
rangeJ = float(max(dataSet[:, j]) - minJ)
#质心取范围中的一个随机值
centroids[:, j] = np.mat(minJ + rangeJ * np.random.rand(k, 1))
#返回质心矩阵
return centroids
#传入数据集和选定簇数k,将数据集分为k个簇
#距离计算函数默认取欧式距离,初始质心生成函数默认取随机生成
def kMeans(dataSet, k, distMeas=distEclud, createCent=randCent):
m = np.shape(dataSet)[0]
#构造m*2矩阵用来记录每个数据所属的质心以及离质心的距离
clusterAssment = np.mat(np.zeros((m, 2)))
#随机生成初始质心
centroids = createCent(dataSet, k)
#用标志记录簇是否改变,若没有改变则分类完成
clusterChanged = True
while clusterChanged:
#重置标志
clusterChanged = False
#遍历每组数据
for i in range(m):
#初始化簇序号和离质心的距离
minDist = np.inf
minIndex = -1
#遍历每个质心
for j in range(k):
#计算此数据离此质心的距离,若小于当前最小距离则更新
distJI = distMeas(centroids[j, :], dataSet[i, :])
if distJI < minDist:
minDist = distJI
minIndex = j
#若更新后的最近质心与更新前的不一样则标志记为True
if clusterAssment[i, 0] != minIndex: clusterChanged = True
#更新记录矩阵
clusterAssment[i, :] = minIndex, minDist ** 2
#每遍历完所有数据之后更新每个质心的坐标
for cent in range(k):
#得到所有属于此质心的数据
ptsInClust = dataSet[np.nonzero(clusterAssment[:, 0].A == cent)[0]]
#重新计算质心
centroids[cent, :] = np.mean(ptsInClust, axis=0)
#返回最终质心矩阵和记录矩阵
return centroids, clusterAssment
#二分Kmeans聚类算法,输入数据集和指定簇数,返回质心矩阵和记录矩阵
def biKmeans(dataSet, k, distMeans=distEclud):
#初始化记录矩阵
m = np.shape(dataSet)[0]
clusterAssment = np.mat(np.zeros((m, 2)))
#初始化第一个质心为所有数据的质心并放入质心矩阵
centroid0 = np.mean(dataSet, axis=0)[0].tolist()
centList = [centroid0]
#根据初始质心更新记录矩阵
for j in range(m):
clusterAssment[j, 1] = distMeans(np.mat(centroid0), dataSet[j, :]) ** 2
#当质心数量小于指定的k时循环
while (len(centList) < k):
#用来观察运行进度
print("第%d次循环" % len(centList))
#初始化最佳分类方案的总SSE
lowestSSE = np.inf
#遍历质心矩阵中的每个质心
for i in range(len(centList)):
#得到属于此质心的簇的数据
ptsInCurrCluster = dataSet[np.nonzero(clusterAssment[:, 0].A == i)[0],:]
#将这个簇分为两个簇
centroidMat, splitClustAss = kMeans(ptsInCurrCluster, 2, distMeans)
#计算分类后的总SSE是否比当前最小SSE小
sseSplit = np.sum(splitClustAss[:, 1])
sseNotSplit = np.sum(clusterAssment[np.nonzero(clusterAssment[:, 0].A != i)[0], 1])
#若小,则更新最佳的分类簇和最小的SSE
if (sseSplit + sseNotSplit) < lowestSSE:
bestCentToSplit = i
bestNewCents = centroidMat
bestClustAss = splitClustAss.copy()
lowestSSE = sseSplit + sseNotSplit
#将最佳分类簇分离的两个簇记为新的簇
bestClustAss[np.nonzero(bestClustAss[:, 0].A == 1)[0], 0] = len(centList)
bestClustAss[np.nonzero(bestClustAss[:, 0].A == 0)[0], 0] = bestCentToSplit
#追踪结果
print('the bestCentToSplit is: ', bestCentToSplit)
print('the len of bestClustAss is: ', len(bestClustAss))
#更新质心矩阵
centList[bestCentToSplit] = bestNewCents[0, :].tolist()[0]
centList.append(bestNewCents[1, :].tolist()[0])
#更新记录矩阵
clusterAssment[np.nonzero(clusterAssment[:, 0].A == bestCentToSplit)[0],:] = bestClustAss
#返回质心矩阵和记录矩阵
return np.mat(centList), clusterAssment
#用pickle方法将结果存在文件中
def storeResult(input, storeFilename):
with open(storeFilename, 'wb') as fw:
pickle.dump(input, fw)
#主函数
if __name__ == '__main__':
#从文件中读取数据
filename = r"D:\python_things\code\第5次作业更新\Kmeans数据集\Sales_Transactions_Dataset_Weekly.csv"
dataSet, normalizedDataSet = loadData(filename)
#用二分Kmeans算法对数据集进行分类,簇数定为200簇
centList, clusterAssment = biKmeans(dataSet, 200)
normalizedCentList, normalizedClusterAssment = biKmeans(normalizedDataSet, 200)
#将分类结果存入文件
normalizedStoreFilename = r"D:\python_things\code\第5次作业更新\normalizedClusterAssment.pkl"
storeFilename = r"D:\python_things\code\第5次作业更新\clusterAssment.pkl"
storeResult(clusterAssment, storeFilename)
storeResult(normalizedClusterAssment, normalizedStoreFilename)
plotData.py
#-*- coding: utf-8 -*-
#Author:Yinli
import pickle
import numpy as np
import matplotlib.pyplot as plt
from scipy import interpolate
from pylab import mpl
#从文件中读取数据
#返回原始周销量数据和正则化后的周销量数据
def loadData(filename):
#逐行读取文件
fr = open(filename, 'r', encoding='utf-8')
arrayOfLines = fr.readlines()
#读取商品数目并初始化数据集
m = len(arrayOfLines)-1
dataSet = np.zeros((m,52))
normalizedDataSet = np.zeros((m,52))
#逐行处理数据
for i in range(1,len(arrayOfLines)):
#去掉句尾的换行符并将数据转为列表
arrayOfLines[i] = arrayOfLines[i].rstrip('\n')
currentLine = arrayOfLines[i].split(',')
#分出原始周销量数据和正则化周销量数据并写入数据集
currentLine0 = currentLine[1:53]
currentLine1 = currentLine[55:]
dataSet[i-1, :] = list(map(int, currentLine0))
normalizedDataSet[i-1,:] = list(map(float,currentLine1))
#返回数据集
return dataSet, normalizedDataSet
#从文件中读取分类结果
def loadResult(storeFilename):
fr = open(storeFilename,'rb')
return pickle.load(fr)
#找到SSE最小的簇,对其中的数据进行插值处理并返回
def preprocessMinSSEData(dataSet, clusterAssment):
# 用SSEList列表记录每个簇的SSE
SSEList = []
for i in range(200):
example = dataSet[np.nonzero(clusterAssment[:, 0].A == i)[0]]
# 若簇的数据数目大于1则计算SSE
if (len(example) > 1):
SSEOfI = np.sum(clusterAssment[np.nonzero(clusterAssment[:, 0].A == i)[0], 1])
SSEList.append(SSEOfI)
# 否则SSE记为无限大
else:
SSEList.append(np.inf)
# 得到SSE最小的簇的索引,并筛选出这个簇的数据并打印其数据个数查看
minIndex = SSEList.index(min(SSEList))
example = dataSet[np.nonzero(clusterAssment[:, 0].A == minIndex)[0]]
numOfData = len(example)
print(len(example))
#x为原始横坐标,xnew为插值后的横坐标
x = np.arange(0, 52, 1)
xnew = np.arange(0, 51, 0.2)
#ynew为插值后的数据矩阵
ynew = np.zeros((numOfData, 255))
#对每组数据进行插值处理
for i in range(numOfData):
func = interpolate.interp1d(x, example[i], kind='cubic')
ynew[i] = func(xnew)
#返回处理后的数据
return xnew,ynew
#主函数
if __name__ == '__main__':
#从文件中读取数据集
filename = "D:\python_things\code\第5次作业更新\Kmeans数据集\Sales_Transactions_Dataset_Weekly.csv"
dataSet, normalizedDataSet = loadData(filename)
#从文件中读取分类结果
storeFilename = "D:\python_things\code\第5次作业更新\clusterAssment.pkl"
normalizedStoreFilename = r"D:\python_things\code\第5次作业更新\normalizedClusterAssment.pkl"
clusterAssment = loadResult(storeFilename)
normalizedClusterAssment = loadResult(normalizedStoreFilename)
#得到两个数据集中SSE最小的簇的插值后的数据
xnew, ynew = preprocessMinSSEData(dataSet, clusterAssment)
normalizedXnew, normalizedYnew = preprocessMinSSEData(normalizedDataSet, normalizedClusterAssment)
#设定字体,初始化画布
mpl.rcParams['font.sans-serif'] = ['SimHei']
plt.figure()
#画第一幅图
plt.subplot(2,1,1)
#有多少组数据画多少条曲线
for i in range(np.shape(ynew)[0]):
plt.plot(xnew,ynew[i])
plt.xlabel(u"周数")
plt.ylabel(u"每周卖出的商品数")
plt.title(u"上——SSE最小的簇的商品销量走势图 下——SSE最小的簇的商品销量走势图(正则化后)")
#画第二幅图
plt.subplot(2,1,2)
#有多少组数据画多少条曲线
for i in range(np.shape(normalizedYnew)[0]):
plt.plot(normalizedXnew, normalizedYnew[i])
plt.xlabel(u"周数")
plt.ylabel(u"每周卖出的商品比例")
plt.show()
运行结果
对原始周销量数据和正则化后的周销量数据都选取了SSE最小的簇中的数据进行了可视化,对这不同的两组数据都画出了它们的平滑过后的销量时间图如下:

结果分析
从两张图中可以看到,对于两个数据集,簇的分类结果中SSE最小的簇中的数据组数都是两组,即都是两种商品的关联性最大,而且从图中可以看到两组数据的曲线重合度很高,趋势相似,可以认为此分类是成功的。
第二题:用FP-Growth算法处理数据
数据来源
http://fimi.ua.ac.be/data/retail.dat
运行环境
python3.7、PyCharm 2018.2.4 (Community Edition)
思路
从所给数据及其说明文档可以看出此数据集是从购物数据中收集而来的,每行数据都是一条购物记录,每行数据中的元素包括了此次购物所买的东西,不考虑数量,只考虑种类,即每行数据都记录了此次购物购买了哪几种商品,每种商品用序号标识,共有88162条记录。现在要求找出哪些商品是顾客经常一起买的,一起买的频率是多大,即找出这些商品中的频繁项集,考虑采用FP-Growth算法进行分析,先构建此数据集的FP树,然后从这棵FP树中搜寻其频繁项集,最后返回打印频繁项集,另外还将FP树进行可视化。
代码包括两个文件,一个文件FP-Growth.py用来读取数据文件,用数据集构建FP树并保存,然后找出其频繁项集打印输出,另一个文件plotTree.py用来读取保存的FP树,然后对其进行可视化(采用graphviz包进行可视化)。
源代码
FP-Growth.py
#-*- coding: utf-8 -*-
#Author:Yinli
import pickle
def loadData(filename):
fr = open(filename, 'r', encoding='utf-8')
arrayOfLines = fr.readlines()
dataSet = []
for i in range(len(arrayOfLines)):
arrayOfLines[i] = arrayOfLines[i].rstrip('\n')
currentList = arrayOfLines[i].split(',')
dataSet.append(currentList)
return dataSet
#定义树节点数据结构
class treeNode:
#构造函数
def __init__(self, nameValue, numOccur, parentNode):
self.name = nameValue#节点名字
self.count = numOccur#节点的计数值
self.nodeLink = None#用于指向下一个名字相同的节点
self.parent = parentNode#指向父节点
self.children = {}#子节点用字典存储
#增加节点计数值
def inc(self, numOccur):
self.count += numOccur
#打印输出此节点为根的树
def disp(self, ind=1):
print(' ' *ind, self.name, ' ', self.count)
for child in self.children.values():
child.disp(ind +1)
#用传入的数据集构造一个字典,键为每行数据的不可变集合,值为1,即数量
def createInitSet(dataSet):
retDict = {}
for trans in dataSet:
retDict[frozenset(trans)] = 1
return retDict
def createTree(dataSet, minSup=1):
#初始化头指针表
headerTable = {}
#遍历数据集中每个单个数据,统计频度
for trans in dataSet:
for item in trans:
headerTable[item] = headerTable.get(item, 0) + dataSet[trans]
#遍历头指针表,删除频度小于指定值的元素
for k in list(headerTable.keys()):
if headerTable[k] < minSup:
del(headerTable[k])
#统计频繁单元素,构造集合
freqItemSet = set(headerTable.keys())
#如果没有频繁单元素返回None
if len(freqItemSet) == 0:
return None,None
#字典值中添加每个元素的指针,初始化为None
for k in headerTable:
headerTable[k] = [headerTable[k], None]
#根节点设为空
retTree = treeNode('Null Set', 1, None)
#遍历每组数据
for tranSet, count in dataSet.items():
#初始化记录每组数据的字典
localD = {}
#遍历数据中的每个元素,若频繁则记录其频度
for item in tranSet:
if item in freqItemSet:
localD[item] = headerTable[item][0]
#若这组数据中有频繁元素
if len(localD)>0:
#将数据中的元素按频度从大到小排序
orderedItems = [v[0] for v in sorted(localD.items(),key=lambda p:p[1], reverse=True)]
#更新节点
updateTree(orderedItems, retTree, headerTable, count)
return retTree, headerTable
#更新节点
#items是处理过后的只含频繁元素的按频度排序的一组数据
#inTree是要更新的节点,count是此组数据的频度
def updateTree(items, inTree, headerTable, count):
#若第一个元素在子节点中,增加频度,否则生成子节点
if items[0] in inTree.children:
inTree.children[items[0]].inc(count)
#否则生成子节点并更新头指针表
else:
inTree.children[items[0]] = treeNode(items[0], count, inTree)
#若当前元素在头指针表中的指针指向空,则令其指向当前子节点
if headerTable[items[0]][1] == None:
headerTable[items[0]][1] = inTree.children[items[0]]
#否则遍历当前元素头指针,将当前节点接在最末尾
else:
updateHeader(headerTable[items[0]][1], inTree.children[items[0]])
#若还有频繁元素,递归更新子节点
if len(items) > 1:
updateTree(items[1::], inTree.children[items[0]], headerTable, count)
#更新头指针表
def updateHeader(nodeToTest, targetNode):
#指针一路向下,直到最末尾
while(nodeToTest.nodeLink != None):
nodeToTest = nodeToTest.nodeLink
#将目标节点接在最末尾
nodeToTest.nodeLink = targetNode
#从leafNode节点一路向上开始找它的前缀路径
def ascendTree(leafNode, prefixPath):
if leafNode.parent != None:
prefixPath.append(leafNode.name)
ascendTree(leafNode.parent, prefixPath)
#找treeNode及其相似节点的前缀路径集合即条件模式基并存入返回
def findPrefixPath(basePat, treeNode):
#初始化前缀路径集合字典
condPats = {}
while treeNode != None:
#初始化前缀路径
prefixPath = []
#向上回溯记录前缀路径
ascendTree(treeNode, prefixPath)
#若前缀路径存在,存入字典,路径为键,频数为值
if len(prefixPath) > 1:
condPats[frozenset(prefixPath[1:])] = treeNode.count
#继续遍历下一个名字相同的相似节点
treeNode = treeNode.nodeLink
return condPats
#挖掘FP树中的频繁项集,minSup为阈值
#prefix为当前的前缀,freqItemList用来记录频繁项集
def mineTree(inTree, headerTable, minSup, prefix, freqItemList):
#首先将头指针表中的频繁元素从小到大排序
bigL = [v[0] for v in sorted(headerTable.items(), key=lambda p:p[0])]
#遍历每个频繁元素
for basePat in bigL:
#newFreqSet记录加了当前节点的前缀
newFreqSet = prefix.copy()
newFreqSet.add(basePat)
#将当前前缀加入到频繁项集列表中
freqItemList.append(newFreqSet)
#求当前元素的前缀路径集及条件模式基
condPattBases = findPrefixPath(basePat, headerTable[basePat][1])
#对此条件模式基构造FP树
myCondTree, myHead = createTree(condPattBases, minSup)
#若树不为空,即此条件模式基中含有频繁元素
#则递归地挖掘此FP树中的频繁项集
if(myHead != None):
mineTree(myCondTree, myHead, minSup, newFreqSet, freqItemList)
#用pickle方法将结果存在文件中
def storeResult(input, storeFilename):
with open(storeFilename, 'wb') as fw:
pickle.dump(input, fw)
#主函数
if __name__ == '__main__':
#从文件中导入数据
filename = r"D:\python_things\code\第5次作业更新\FP树数据集\retail.csv"
dataSet = loadData(filename)
#将数据集初始化为字典
dataSet = createInitSet(dataSet)
#设定阈值
minSup = 10000
#对数据集构造FP树并打印查看
myTree, myHeaderTable = createTree(dataSet, minSup)
myTree.disp()
treeFilename = r"D:\python_things\code\第5次作业更新\storedTree.pkl"
storeResult(myTree, treeFilename)
#初始化频繁项集列表,挖掘FP树的频繁项集并打印查看
freqItems = []
mineTree(myTree, myHeaderTable, 5000 ,set([]), freqItems)
print("频繁项集:")
for i in range(len(freqItems)):
if(len(freqItems[i]) > 1):
print(freqItems[i])
plotTree.py
#-*- coding: utf-8 -*-
#Author:Yinli
from graphviz import Digraph
import pickle
#定义树节点数据结构
class treeNode:
#构造函数
def __init__(self, nameValue, numOccur, parentNode):
self.name = nameValue#节点名字
self.count = numOccur#节点的计数值
self.nodeLink = None#用于指向下一个名字相同的节点
self.parent = parentNode#指向父节点
self.children = {}#子节点用字典存储
#增加节点计数值
def inc(self, numOccur):
self.count += numOccur
#打印输出此节点为根的树
def disp(self, ind=1):
print(' ' *ind, self.name, ' ', self.count)
for child in self.children.values():
child.disp(ind +1)
#从文件中读取分类结果
def loadResult(storeFilename):
fr = open(storeFilename,'rb')
return pickle.load(fr)
#传入结点,将以此结点为根的FP树画出来
def plotTree(inTree, dot):
#结点名字
currentNodeName = 'No.' + inTree.name + ' : ' + str(inTree.count)
#结点标签,用于标识每个唯一的结点
currentNodeLabel = inTree.name + '_' + str(inTree.count)
#将当前结点画出来
dot.node(currentNodeLabel, currentNodeName)
#如果此结点有孩子
if (inTree.children != {}):
#则对于每个孩子递归调用plotTree
for key in inTree.children.keys():
#返回此孩子结点的唯一标签
childrenNodeLabel = plotTree(inTree.children[key], dot)
#将当前结点和此孩子结点连接起来
dot.edge(currentNodeLabel,childrenNodeLabel)
#返回当前结点的唯一标签
return currentNodeLabel
#主函数
if __name__ == '__main__':
#读取FP树的根节点
treeFilename = r"D:\python_things\code\第5次作业更新\storedTree.pkl"
myTree = loadResult(treeFilename)
#将FP树画出来并保存
dot = Digraph(comment='tree')
plotTree(myTree, dot)
dot.render('tree.gv', view=True)
运行结果
将构建FP树的最小阈值设定为10000,即在这8w+条数据中频度超过10000才算频繁项,寻找频繁项集的最小阈值也设定为5000,得到的FP树和频繁项集如下所示:
FP树

频繁项集

结果分析
从结果来看,一起购买次数超过5000次的频繁项集有以上10组,这为卖场的营销策略提供了支持,可以认为分析工作是有效的。