python函⼀维聚类_⼀维数组的聚类
需求:分析订单的价格分布
⽅案:按照100为梯度,分析不同价格区间的订单量
缺陷:现实⽣活中,定价存在⼀些⾃然的价格分隔,如果按照步距划分可能存在⼀些偏差,⽐如airbnb的价格筛选显⽰出的房价分布:
解决上述缺陷最好的⽅式是对价格进⾏聚类,出做合适的价格区间。
在学习聚类算法的过程中,学习到的聚类算法⼤部分都是针对n维的,针对⼀维数据的聚类⽅式较少,今天就来学习下如何给⼀维的数据进⾏聚类。
⽅案⼀:采⽤K-Means对⼀维数据聚类
Python代码如下:
from sklearn.cluster import KMeans
import numpy as np
x = np.random.random(10000)
y = x.reshape(-1,1)
km = KMeans()
km.fit(y)
核⼼的操作是y = x.reshape(-1,1),含义为将⼀维数据变成只有1列,⾏数不知道多少(-1代表根据剩下的维度计算出数组的另外⼀个shape属性值)。
⽅案⼆:采⽤⼀维聚类⽅法Jenks Natural Breaks
Jenks Natural Breaks(⾃然断点分类)。⼀般来说,分类的原则就是差不多的放在⼀起,分成若⼲类。统计上可以⽤⽅差来衡量,通过计算每类的⽅差,再计算这些⽅差之和,⽤⽅差和的⼤⼩来⽐较分类的好坏。因⽽需要计算各种分类的⽅差和,其值最⼩的就是最优的分类结果(但并不唯⼀)。这也是⾃然断点分类法的原理。另外,当你去看数据的分布时,可以⽐较明显的发现断裂之处,这些断裂之处和Jenks Natural Breaks⽅法算出来也是⼀致的。因⽽这种分类法很“⾃然”。
Jenks Natural Breaks和K Means在⼀维数据时,完全等价。它们的⽬标函数⼀样,但是算法的步骤不
完全相同。K Means是先设定好K 个初始随机点。⽽Jenks Breaks则是⽤遍历的⽅法,⼀个点⼀个点地移动,直到达到最⼩值。
Natural Breaks算法⼜有两种:
Jenks-Caspall algorithm(1971),是Jenks和Caspall发明的算法。原理就如前所述,实现的时候要将每种分类情况都计算⼀遍,到⽅差和最⼩的那⼀种,计算量极⼤。n个数分成k类,就要从n-1个数中k-1个组合,这个数⽬是很惊⼈的。数据量较⼤时,如果分类⼜多,以当时的计算机⽔平根本不能穷举各种可能性。
Fisher-Jenks algorithm(1977),Fisher(1958)发明了⼀种算法提⾼计算效率,不需要进⾏穷举。Jenks将这种⽅法引⼊到数据分类中。但后来者⼏乎只知道Jenks⽽不知Fisher了。
具体算法实现:
和K-Means⼀样,使⽤Jenks Natural Breaks需要先确定聚类数量K值。常见的⽅法是:GVF(The Goodness of Variance Fit)。GVF,翻译过来是“⽅差拟合优度”,公式如下:
其中,SDAM是the Sum of squared Deviations from the Array Mean,即原始数据的⽅差;SDCM是the Sum of
squared Deviations about Class Mean,即每⼀类⽅差的和。显然,SDAM是⼀个常数,⽽SDCM与分类数k有关。⼀定范围内,GVF 越⼤,分类效果越好。SDCM越⼩,GVF越⼤,越接近于1。⽽SDCM随k的增⼤⽽⼤,当k等于n时,SDMC=0,GVF=1。
GVF⽤于判定不同分类数的分类效果好坏。以k和GVF做图可得:
随着k的增⼤,GVF曲线变得越来越平缓。特别是在红线处(k=5),曲线变得基本平坦(之前起伏较⼤,之后起伏较⼩),k(5)也不是很⼤,所以可以分为5类。⼀般来说,GVF>0.7就可以接受了,当然越⾼越好,但⼀定要考虑k不能太⼤。显然,这是⼀个经验公式,但总⽐没有好吧。linspace函数python
代码⽰例:
from jenkspy import jenks_breaks
import numpy as np
def goodness_of_variance_fit(array, classes):
# get the break points
classes = jenks_breaks(array, classes)
# do the actual classification
classified = np.array([classify(i, classes) for i in array])
# max value of zones
maxz = max(classified)
# nested list of zone indices
zone_indices = [[idx for idx, val in enumerate(classified) if zone + 1 == val] for zone in range(maxz)]
# sum of squared deviations from array mean
sdam = np.sum((array - an()) ** 2)
# sorted polygon stats
array_sort = [np.array([array[index] for index in zone]) for zone in zone_indices]
# sum of squared deviations of class means
sdcm = sum([np.sum((classified - an()) ** 2) for classified in array_sort])
# goodness of variance fit
gvf = (sdam - sdcm) / sdam
return gvf
def classify(value, breaks):
for i in range(1, len(breaks)):
if value < breaks[i]:
return i
return len(breaks) - 1
if __name__ == '__main__':
gvf = 0.0
nclasses = 2
array = np.random.random(10000)
while gvf < .8:
gvf = goodness_of_variance_fit(array, nclasses)
print(nclasses, gvf)
nclasses += 1
参考链接:
⽅案三:核密度估计Kernel Density Estimation
所谓核密度估计,就是采⽤平滑的峰值函数(“核”)来拟合观察到的数据点,从⽽对真实的概率分布曲线进⾏模拟。核密度估计更多详细内容,可以参考先前的Mean Shift聚类中的相关说明。
使⽤⽰例:
import numpy as np
from scipy.signal import argrelextrema
import matplotlib.pyplot as plt
ighbors.kde import KernelDensity
a = np.array([10, 11, 9, 23, 21, 11, 45, 20, 11, 12]).reshape(-1, 1)
kde = KernelDensity(kernel='gaussian', bandwidth=3).fit(a)
s = np.linspace(0, 50)
e = kde.score_shape(-1, 1))
plt.plot(s, e)
plt.show()
mi, ma = argrelextrema(e, np.less)[0], argrelextrema(e, np.greater)[0]
print("Minima:", s[mi])
print("Maxima:", s[ma])
print(a[a < mi[0]], a[(a >= mi[0]) * (a <= mi[1])], a[a >= mi[1]])
plt.plot(s[:mi[0] + 1], e[:mi[0] + 1], 'r',
s[mi[0]:mi[1] + 1], e[mi[0]:mi[1] + 1], 'g',
s[mi[1]:], e[mi[1]:], 'b',
s[ma], e[ma], 'go',
s[mi], e[mi], 'ro')
plt.show()
输出内容:
Minima: [17.34693878 33.67346939]
Maxima: [10.20408163 21.42857143 44.89795918]
[10 11  9 11 11 12] [23 21 20] [45]
参考链接:

版权声明:本站内容均来自互联网,仅供演示用,请勿用于商业和其他非法用途。如果侵犯了您的权益请与我们联系QQ:729038198,我们将在24小时内删除。