Scipyのcut_tree（）は、要求された数のクラスタを返さず、scipyとfastclusterで得られたリンケージ行列は一致しません

Question

scipy.cluster.hierarchyモジュール関数に関連してfastclusterパッケージを使用して凝集クラスタークラスタリング（AHC） Python 3にあり、私はcut_tree()機能の困惑している動作を見つけました。

私は問題なしでデータをクラスタ化し、を使用してmethod=wardを使用してZを取得します。次に、一定の数のクラスター（例えば、33）を得るために樹状図ツリーを切断したいと思って、これを適切にcut_tree(Z, n_clusters=33)を使用して行います。 （AHCは、ツリーのリーフにあるすべてのデータポイントを結ぶバイナリツリーを生成する決定論的な方法であることを覚えておいてください。このツリーをどのレベルで見ても、最後に必要なクラスタの数を見ることができます。 cut_tree（）は、 'n_cluster'整数ラベルのセットを0からn_clusters - 1に戻すことです。

私はこれを他の実験で何度もやっていますが、私はいつも私が要求するクラスタの数。問題は、この1つのデータセットでは、cut_tree()に33クラスタを依頼すると、私には32だけが与えられます。なぜこれが当てはまるのかわかりません。それはバグでしょうか？ cut_tree()のバグを認識していますか？私はこの動作をデバッグしようとし、scipyのlinkage()機能を使って同じクラスタリング実験を行った。得られたリンケージマトリックスをcut_tree()に入力すると、予想外の数のクラスタが出力として取得されませんでした。私はまた、2つの方法で出力されたリンケージ行列が等しくないことを確認しました。

私が使用している[dataset]は、それぞれが20次元の10680個のベクトルで構成されています。以下の実験を確認：あなたはデータセットが少なくとも正確なコピーで37個のベクトルが含まれて気づいたかもしれません

import numpy as np 
import fastcluster as fc 
import scipy.cluster.hierarchy as hac 
from scipy.spatial.distance import pdist 

### *Load dataset (10680 vectors, each with 20 dimensions)* 
X = np.load('dataset.npy') 

### *Hierarchical clustering using traditional scipy method* 
dists = pdist(X) 
Z_1 = hac.linkage(dists, method='ward') 

### *Hierarchical clustering using optimized fastcluster method* 
Z_2 = fc.linkage_vector(X, method='ward') 

### *Comparissons* 

## Are the linkage matrices equal? 
print("Z_1 == Z_2 ? ", np.allclose(Z_1, Z_2)) 

## Is scipy's cut_tree() returning the requested number of clusters when using Z_2? 
print("Req.\tGot\tequal?") 
for i in range(1,50): 
    cut = hac.cut_tree(Z_2, i) 
    uniq = len(np.unique(cut)) 
    print(i,"\t",uniq,"\t",i==uniq) 

## The same as before, but in condensed form. When requesting cut_tree() for clusters 
# in the range [1,50] does it return wrong results at some point? 
print("Any problem cutting Z_1 for n_clusters in [1,50]? ", not np.all([len(np.unique(
             hac.cut_tree(Z_1, i)))==i for i in range(1,50)])) 
print("Any problem cutting Z_2 for n_clusters in [1,50]? ", not np.all([len(np.unique(
             hac.cut_tree(Z_2, i)))==i for i in range(1,50)])) 

#Output: 
# 
#Z_1 == Z_2 ? False 
# 
#Req. Got  equal? 
#1  1  True 
#2  2  True 
#3  3  True 
#4  4  True 
#5  5  True 
#6  6  True 
#7  7  True 
#8  8  True 
#9  9  True 
#10  10  True 
#11  11  True 
#12  12  True 
#13  13  True 
#14  14  True 
#15  15  True 
#16  16  True 
#17  17  True 
#18  18  True 
#19  19  True 
#20  20  True 
#21  21  True 
#22  22  True 
#23  23  True 
#24  24  True 
#25  25  True 
#26  26  True 
#27  27  True 
#28  28  True 
#29  29  True 
#30  30  True 
#31  31  True 
#32  32  True 
#33  32  False 
#34  33  False 
#35  34  False 
#36  35  False 
#37  36  False 
#38  37  False 
#39  38  False 
#40  39  False 
#41  40  False 
#42  41  False 
#43  42  False 
#44  43  False 
#45  44  False 
#46  45  False 
#47  46  False 
#48  47  False 
#49  48  False 
# 
#Any problem cutting Z_1 for n_clusters in [1,50]? False 
#Any problem cutting Z_2 for n_clusters in [1,50]? True 

を、そしてすべてのコピーをカウントすることは、データセット内の少なくともコピーで、合計55個のベクトルがあります。

検査のために、私は2つのリンケージマトリックスの浅い深さレベルまでプロットすることに決めました。画像は下の図のようにZ_1、下にはZ_2です。括弧内の数字は、そのブランチの蛇腹を含む集団を示します。かっこのない数字はツリーのリーフです（数値はXマトリックスのベクトルのインデックスです）。 （プロットされたレベルでの）唯一の違いは、赤い四角でマークされた枝にあり、重なり合ったベクトルを含んでいるので0の距離で合体していることが分かります。

だから、私は再び、前のコードに示すように、クラスタリング手順を実行したが、少なくともコピーが55のベクターを含むデータのサブセットのみを使用してこの時間。したがって、fastclusterとscipyのダウンロードは、同じ結果を返すされていません

#Z_1 == Z_2 ? False 
#Req. Got  equal? 
#1  1  True 
#2  2  True 
#3  2  False 
#4  3  False 
#5  4  False 
#6  5  False 
#7  6  False 
#8  7  False 
#9  8  False 
#10  9  False 
#11  10  False 
#12  11  False 
#13  12  False 
#14  13  False 
#15  14  False 
#16  15  False 
#17  16  False 
#18  17  False 
#19  18  False 
#20  20  True 
#21  21  True 
#22  22  True 
#23  23  True 
#24  24  True 
#25  25  True 
#26  26  True 
#27  27  True 
#28  28  True 
#29  29  True 
#30  30  True 
#31  31  True 
#32  32  True 
#33  33  True 
#34  34  True 
#35  35  True 
#36  36  True 
#37  37  True 
#38  38  True 
#39  39  True 
#40  40  True 
#41  41  True 
#42  42  True 
#43  43  True 
#44  44  True 
#45  45  True 
#46  46  True 
#47  47  True 
#48  48  True 
#49  49  True 
#Any problem cutting Z_1 for n_clusters in [1,50]? False 
#Any problem cutting Z_2 for n_clusters in [1,50]? True 

、それが唯一の原因重複の点にある場合、これは可能性：この時間に

uniqs, uniqs_indices, uniqs_count = np.unique(X, axis=0, return_index=True, return_counts=True) 
duplicate_rows_indices = list(set(range(len(X))) - set(uniqs_indices)) 
number_of_duplicate_rows = len(X)-len(uniqs) # 37 

all_duplicate_rows = set() 
for i in duplicate_rows_indices: 
    _rows = set(np.where(X == X[i])[0]) 
    for j in _rows: 
     all_duplicate_rows.add(j) 

rows_with_at_least_a_copy = list(all_duplicate_rows) 
number_of_rows_with_at_least_a_copy = len(rows_with_at_least_a_copy) # 55 

X_subset = X[rows_with_at_least_a_copy] 

と私の出力ました：私はX_subsetを得ましたそのクラスタリング状況のあいまいさのために受け入れられるべきである。しかし、問題はcut_tree()であり、linkage_vector()によって得られた連結マトリックスが与えられた場合、これらの場合に要求された数のクラスターを返さないことがある。どのようにこれを修正することができますか？使用

ライブラリのバージョン：scipyのダウンロード '0.19.1'、numpyの '1.13.3'、fastcluster '1.1.24'

編集：はまた、ここに掲載されます：https://github.com/scipy/scipy/issues/7977。

Answer 1

まず、私は3つの方法の出力の違いについて心配していませんscipy.cluster.hierarchy.linkage(),fastcluster.linkage()およびfastcluster.linkage_vector()。浮動小数点演算の制限によるクラスタ距離で

• 小さな差（代数的に等価な式が異なる結果をもたらす）：違いは、2つの理由のために生じ得ます。

• 算術的な相違点とは別に、アルゴリズムでは、異なる結びつきを解決することができます。この文脈において、「タイ」とBとの間の正確に同じ距離を有するクラスタ（A、B）、（C、D）の二対であり、CとD。 OPが書いたように、例えば、プロセスの開始時に距離0の多数の点のペア。アルゴリズムは任意の順序でクラスタ化されます。

第二に、scipy.cluster.hierarchy.cut_tree()は、クラスタの希望数を得ない理由の質問はここで本当に面白い問題です。代数的には、「ワード（Ward）」クラスタリング方法は、デンドログラムにおいていわゆる「反転」を生成することができない。 （1ステップのクラスタリング距離が次のステップのクラスタリング距離よりも大きい場合には反転が生じる）すなわち、階段状樹状図の第3列の距離のシーケンス（linkage()の出力）は、理想的には単調増加シーケンス「ワード」メソッドの場合しかしながら、浮動小数点の不正確さ、OPによって供給されたデータセットに、ステップ35においてクラスタリング距離は次のステップで2.2E-16 fastcluster.linkage_vector()によって、しかし、0.0として報告されている36

scipyのダウンロードのcut_tree()は残念ながらありませんこれらの反転をうまく処理してください。デンドログラム（現在の場合のように）が "深く"ダウンしても、誤って形成されたクラスタの影響を受けて、残りのマージ処理のアルゴリズムが混乱します。 （ちなみに、私は樹状図が間違ったレベルで「カット」されているだけでなく、クラスターが実際には間違っていると信じています）。

これは逆転としてさらに厄介ですこの場合のように数値的に不正確であるだけでなく、「重心」と「中央値」のクラスタリング方法は、完全な算術を使用しても非常に頻繁に生成されます。

最後に、不完全な解決策：次の行によって

for i, node in enumerate(nodes): 
    idx = node.pre_order() 
    this_group = last_group.copy() 
    # ------------- Replace this: 
    this_group[idx] = last_group[idx].min() 
    this_group[this_group > last_group[idx].max()] -= 1 
    # ------------- 
    if i + 1 in cols_idx: 
     groups[np.where(i + 1 == cols_idx)[0]] = this_group 
    last_group = this_group 

：問題を解決するために、scipyのダウンロードのcut_tree()関数内のループに2つの標線置き換える

unique_idx = np.unique(last_group[idx]) 
    this_group[idx] = unique_idx[0] 
    for ui in unique_idx[:0:-1]: 
     this_group[this_group > ui] -= 1 

を（見て文脈のためのSciPy source codeで。）

しかし、現在の実装が非常に複雑で、実行時の複雑さの観点から効率的にタスクを実行できるので、cut_tree()メソッドを最初から再実装することをお勧めします。