标准 1 错误的差距统计

Question

我已经使用 Scikit Learn 命令实现了 Kmeans，并尝试了 Elbow 和 Silhoutte Coefficient 来找到最佳 K。我计划使用间隙统计来进一步验证我的结果。

def optimalK(data, nrefs=3, maxClusters=15):


gaps = np.zeros((len(range(1, maxClusters)),))
resultsdf = pd.DataFrame({'clusterCount':[], 'gap':[]})
for gap_index, k in enumerate(range(1, maxClusters)):

    # Holder for reference dispersion results
    refDisps = np.zeros(nrefs)

    for i in range(nrefs):

        # Create new random reference set
        randomReference = np.random.random_sample(size=data.shape)

        # Fit to it
        km = KMeans(k)
        km.fit(randomReference)

        refDisp = km.inertia_
        refDisps[i] = refDisp

    km = KMeans(k)
    km.fit(data)

    origDisp = km.inertia_

    # Calculate gap statistic
    gap = np.log(np.mean(refDisps)) - np.log(origDisp)

    # Assign this loop's gap statistic to gaps
    gaps[gap_index] = gap

    resultsdf = resultsdf.append({'clusterCount':k, 'gap':gap}, ignore_index=True)

return (gaps.argmax() + 1, resultsdf)  

但是，我的差距统计图正在增加，因此最佳聚类数始终是我的聚类范围的终点。假设我将集群范围定义为 1 到 10，那么最佳值为 10。

根据互联网网站和原始论文，解决方法是实现标准 1 错误，其中

GAP(K)> GAP(K+1)- S(K+1)

谁能向我解释如何在上面的代码中实现这一点？我不知道如何计算 S(k+1)，因为它涉及找到参考分布的标准差。

s(k+1) = sd(k+1)*square_root(1+(1/B))

B 是 Monte Carlo Samples 的副本数。我查看了不同的网站，但似乎他们没有实现标准 1 错误的差距统计。

Answer 1

def gap_stat(data,label):
    k = len(np.unique(label))
    n = data.shape[0]
    p = data.shape[1]
    D_r = []
    C_r = []
    for label_number in range(0,k):
        this_label_index = np.where(label==label_number)[0]
        temp_sum = 0
        pairwise_distance_matrix = 
    euclidean_distances(data[this_label_index],squared=True)
        D_r.append(np.sum(pairwise_distance_matrix))
        C_r.append(float(len(this_label_index)))

    W_r = np.sum(np.asarray(D_r)/(2*np.asarray(C_r)))
    gap_stats = np.log(float(p*n)/12)-(2/float(p))*np.log(k)- 
    np.log(W_r)
    return(gap_stats)