Selecting the number of clusters k
Here is an example using exclusively PyCVI in order to guess the number of clusters in a dataset. The preprocessing steps and the clustering steps can be integrated into the PyCVI pipeline by providing sklearn-like classes of clustering models (e.g. KMeans) and data preprocessor (e.g. StandardScaler).
In this example, we use time-series data and non-time-series data. In addition we use classes from scikit-learn, scikit-learn extra and aeon in order to illustrate the compatibility of PyCVI with sklearn-like libraries.
Here we assume that we are in real conditions, which means that we don’t have access to the true labels (except that we plot the true data for illustrating purpose). We then don’t use the features included in the pycvi.vi module.
If you wish to run the example scripts on your own computer, please first follow the instructions detailed in Running example scripts on your computer.
import numpy as np
import time
from sklearn.cluster import AgglomerativeClustering, KMeans
from sklearn.preprocessing import StandardScaler
from sklearn_extra.cluster import KMedoids
from aeon.clustering import TimeSeriesKMeans
from pycvi.cluster import generate_all_clusterings, get_clustering
from pycvi.cvi import CVIs
from pycvi.compute_scores import compute_all_scores
from pycvi.vi import variation_information
from pycvi.datasets.benchmark import load_data
from pycvi.exceptions import SelectionError
from pycvi_examples_utils import plot_true_best, plot_hist_selected, plot_only_selected
def pipeline(
X: np.ndarray,
y: np.ndarray,
model_class,
model_kw: dict,
k_max: int = 25,
scaler = StandardScaler(),
DTW:bool = False,
fig_title: str = "",
fig_name: str = "",
) -> None:
"""
This function gives an example of typical pipeline using PyCVI.
In this example we assume that we are in real conditions, which
means that we don't have access to the true labels (except for the
final figure). We then don't use the features included in the
:mod:`pycvi.vi` module. In this function we:
- Standardize the data
- Generate all clusterings for a given range of number of clusters.
- For each CVI available in PyCVI, compute their values for all
generated clusterings
- For each CVI, select the best clustering according to its CVI
values
- Create a summary plot containing the true clustering, the
clustering assuming the correct number of clusters, for each
potential number of clusters :math:`k`, the number of CVI that
selected :math:`k` and the selected clusterings of each CVI.
"""
print(f'\n ***** {fig_title} ***** \n')
k_range = range(k_max)
# ------------------------------------------------------------------
# ------------------ Define true clustering -----------------------
# ------------------------------------------------------------------
# From the label for each datapoint to a list of
# datapoints for each cluster.
# true clusters: List[List[int]]
true_clusters = get_clustering(y)
k_true = len(true_clusters)
# ------------------------------------------------------------------
# ------------------ Generate clusterings -------------------------
# ------------------------------------------------------------------
t_start = time.time()
clusterings = generate_all_clusterings(
X,
model_class,
model_kw=model_kw,
n_clusters_range = k_range,
DTW = DTW,
scaler=scaler,
)
t_end = time.time()
dt = t_end - t_start
print(f"Clusterings generated in: {dt:.2f}s")
# ------------------------------------------------------------------
# ------------ Compute CVI values and select k ---------------------
# ------------------------------------------------------------------
summary = {}
for cvi_class in CVIs:
# Instantiate a CVI model
cvi = cvi_class()
t_start = time.time()
print(f" ====== {cvi} ====== ")
# Compute CVI values for all clusterings
scores = compute_all_scores(
cvi,
X,
clusterings,
DTW=DTW,
scaler=StandardScaler(),
)
t_end = time.time()
dt = t_end - t_start
# Print CVI values and selected k
for k in clusterings:
print(k, scores[k])
print('Code executed in %.2f s' %(dt))
# Select k
try:
k_selected = cvi.select(scores)
# If no k could be selected with this CVI don't do anything...
except SelectionError as e:
k_selected = None
# Otherwise update summary info related to the selected clusterings
else:
if k_selected not in summary:
summary[k_selected] = {}
summary[k_selected]["clustering"] = clusterings[k_selected]
summary[k_selected]["#CVI"] = summary[k_selected].pop("#CVI", 0) + 1
ax_title = f'k={k_selected}, #CVI={summary[k_selected]["#CVI"]}'
summary[k_selected]["ax_title"] = ax_title
finally:
print(f"Selected k: {k_selected} | True k: {k_true}")
# ------------------------------------------------------------------
# --------------------- Plot true clustering -----------------------
# ------------------------------------------------------------------
best_clusters = clusterings[k_true]
# -- Plot true clusters & clusters when assuming k_true clusters ---
vi = variation_information(true_clusters, clusterings[k_true])
fig = plot_true_best(
X, y, best_clusters, n_plots=len(summary) + 2, VI_best=vi
)
# ------------------------------------------------------------------
# ----------------------- Summary plot -----------------------------
# ------------------------------------------------------------------
fig_hist = plot_hist_selected(summary)
fig_hist.savefig(fig_name + "-histogram.png")
fig = plot_only_selected(X, summary, fig)
fig.suptitle(fig_title)
fig.savefig(fig_name + ".png")
# ======================================================================
# PyCVI on non time-series data
# ======================================================================
# ------------- KMeans ------------------------
X, y = load_data("zelnik1", "barton")
DTW = False
k_max = 10
model_class = KMeans
model_kw = {}
scaler = StandardScaler()
fig_title = "Non time-series data with KMeans"
fig_name = "Barton_data_KMeans"
pipeline(X, y, model_class, model_kw, k_max, scaler, DTW, fig_title, fig_name)
# --------- AgglomerativeClustering ----------
X, y = load_data("zelnik1", "barton")
DTW = False
k_max = 10
model_class = AgglomerativeClustering
# sklearn kwargs for AgglomerativeClustering
model_kw = {"linkage" : "single"}
scaler = StandardScaler()
fig_title = "Non time-series data with AgglomerativeClustering-Single"
fig_name = "Barton_data_AgglomerativeClustering_Single"
pipeline(X, y, model_class, model_kw, k_max, scaler, DTW, fig_title, fig_name)
# ======================================================================
# PyCVI on time series data
# ======================================================================
X, y = load_data("Trace", "UCR")
# ==========================
# PyCVI using DTW
# ==========================
DTW = True
model_class = TimeSeriesKMeans
# aeon kwargs for TimeSeriesKMeans
model_kw = {
"distance" : "dtw",
"distance_params" : {"window": 0.2},
}
scaler = StandardScaler()
fig_title = "Time-series data using DTW with TimeSeriesKMeans"
fig_name = "UCR_data_DTW_TimeSeriesKMeans"
pipeline(X, y, model_class, model_kw, k_max, scaler, DTW, fig_title, fig_name)
# ==========================
# PyCVI not using DTW
# ==========================
DTW = False
model_class = KMeans
model_kw = {}
scaler = StandardScaler()
fig_title = "Time-series data without DTW with KMeans"
***** Non time-series data with KMeans *****
Clusterings generated in: 0.23s
====== Hartigan_monotonous ======
0 231.31309486007345
1 119.53017448266189
2 117.61566423046895
3 75.96832482873897
4 95.0656001459015
5 57.30656920145313
6 51.687707164140214
7 59.10992898474898
8 62.139148832813106
9 None
Code executed in 0.01 s
Selected k: 1 | True k: 3
====== CalinskiHarabasz_original ======
0 None
1 0.0
2 119.53017448266186
3 142.03940200451063
4 143.99871868455529
5 166.202708270155
6 169.8882540473872
7 174.67976151161565
8 187.9661510034123
9 206.79303061701037
Code executed in 0.00 s
Selected k: 9 | True k: 3
====== GapStatistic_original ======
0 None
1 0.006792153575967852
2 -0.15580672581154076
3 -0.3206999727877289
4 -0.514634591757889
5 -0.43655796457858465
6 -0.43133194253892704
7 -0.4440182169815916
8 -0.4464424904973012
9 -0.36545923888575427
Code executed in 1.48 s
Selected k: 1 | True k: 3
====== Silhouette ======
0 None
1 None
2 0.2942832123441149
3 0.31865535671016215
4 0.30974666695744546
5 0.3758686461703666
6 0.40405671642495117
7 0.4088473534363824
8 0.42222482515851223
9 0.4565121748825993
Code executed in 0.40 s
Selected k: 9 | True k: 3
====== ScoreFunction ======
0 None
1 0.2674380756639073
2 0.1551856671096168
3 0.08708406007354763
4 0.04494313141469164
5 0.026431295225813223
6 0.019498243943131932
7 0.01490000324244145
8 0.01017110630012763
9 0.010637054270650426
Code executed in 0.01 s
Selected k: 1 | True k: 3
====== MaulikBandyopadhyay_absolute ======
0 None
1 0.0
2 0.9214600886394652
3 1.0082609641450428
4 1.2705442508674658
5 1.5944981818915518
6 1.8107219181634564
7 1.7428007278275894
8 1.600835581724516
9 1.6759021856870129
Code executed in 0.01 s
Selected k: 6 | True k: 3
====== SD ======
0 None
1 None
2 421342.90086491493
3 307736.47143837134
4 249804.42413208555
5 246241.62504429856
6 180150.512201299
7 172320.84357779022
8 161626.49256517325
9 120547.3174032998
Code executed in 0.06 s
Selected k: 9 | True k: 3
====== SDbw ======
0 None
1 None
2 1.8032473069048642
3 1.7969053311140333
4 1.300067878472333
5 0.8701205967905576
6 0.6296604573287411
7 0.49424608481289706
8 0.4381264613272081
9 0.34855988186953746
Code executed in 0.02 s
Selected k: 9 | True k: 3
====== Dunn ======
0 None
1 None
2 0.015708595341348486
3 0.0209067934147451
4 0.021605533072104735
5 0.0256934220523547
6 0.03066494732712681
7 0.011161718323872688
8 0.02240543249776124
9 0.031620837931671024
Code executed in 0.01 s
Selected k: 9 | True k: 3
====== XB ======
0 None
1 None
2 0.5636434569575249
3 0.3626913977341436
4 0.37083302077190367
5 0.2616447893620177
6 0.22242970911602786
7 0.32616085580193266
8 0.4167901373745973
9 0.3395171720212871
Code executed in 0.01 s
Selected k: 6 | True k: 3
====== XB_star ======
0 None
1 None
2 0.5870203811596103
3 0.3896403429652914
4 0.4421355679560954
5 0.40272258313838305
6 0.3024200706004079
7 0.4701424556905313
8 0.7846567580586185
9 0.5580896588638549
Code executed in 0.01 s
Selected k: 6 | True k: 3
====== DB ======
0 None
1 None
2 inf
3 inf
4 inf
5 inf
6 inf
7 inf
8 inf
9 inf
Code executed in 0.01 s
Selected k: None | True k: 3
====== Inertia_sum ======
0 521.512709434739
1 349.0214135130419
2 305.7195598105686
3 266.3545714359798
4 243.23239039370878
5 209.33097190122675
6 193.642076374862
7 175.03286012458017
8 158.3277907300725
9 145.75774903454402
Code executed in 0.00 s
Selected k: 1 | True k: 3
====== Diameter_max ======
0 6.206391095165446
1 4.737154354950527
2 4.658670749537977
3 3.8030177966354515
4 3.2470070563718183
5 3.0945238538036546
6 2.5928271318607923
7 2.462705545799648
8 2.595285468655654
9 2.2992948404264295
Code executed in 0.00 s
Selected k: 1 | True k: 3
***** Non time-series data with AgglomerativeClustering-Single *****
Clusterings generated in: 0.02s
====== Hartigan_monotonous ======
0 229.88306277879883
1 0.017211966847463733
2 0.011381644034122118
3 34.58290182504769
4 35.42587757523952
5 1.2950951116888347
6 60.33785949064848
7 114.59147910894787
8 1.8539186101937388
9 None
Code executed in 0.01 s
Selected k: 1 | True k: 3
====== CalinskiHarabasz_original ======
0 None
1 0.0
2 0.01721196684745756
3 0.0142681588644381
4 11.538225250696614
5 18.52000323839848
6 15.089892803326242
7 25.17787001256915
8 46.34652110326673
9 40.903946263374344
Code executed in 0.01 s
Selected k: 8 | True k: 3
====== GapStatistic_original ======
0 None
1 0.005679770808947104
2 -0.49967338072825473
3 -0.974896598975767
4 -1.2931520658810287
5 -1.370677078478475
6 -1.5551364962805523
7 -1.5661162048009238
8 -1.3999772484305049
9 -1.5063524913352189
Code executed in 1.75 s
Selected k: 1 | True k: 3
====== Silhouette ======
0 None
1 None
2 0.19510829001256025
3 0.12276877018507622
4 0.2666718402395447
5 0.31348555456323224
6 0.36558953430010227
7 0.3731558417803549
8 0.4027680246588693
9 0.4166374121873863
Code executed in 0.32 s
Selected k: 9 | True k: 3
====== ScoreFunction ======
0 None
1 0.2674380756639073
2 0.05526989678867322
3 0.03930535158636095
4 0.02794482468974735
5 0.022399044305932536
6 0.02434609671721777
7 0.01898781804203109
8 0.015903727593373773
9 0.016224624922212327
Code executed in 0.01 s
Selected k: 1 | True k: 3
====== MaulikBandyopadhyay_absolute ======
0 None
1 0.0
2 0.0001624885689488273
3 0.00011507955476191883
4 0.5826084938588393
5 0.5140899521509295
6 0.37313655054752143
7 0.43569610764369654
8 0.6956504829553489
9 0.5612754701918737
Code executed in 0.00 s
Selected k: 8 | True k: 3
====== SD ======
0 None
1 None
2 783980.3997208492
3 549586.6578302318
4 432027.73873390944
5 341231.9036186975
6 275681.20322029595
7 240429.69366993732
8 195142.12738698898
9 169718.94353532573
Code executed in 0.05 s
Selected k: 9 | True k: 3
====== SDbw ======
0 None
1 None
2 2.350951786083218
3 2.1462945100388757
4 2.7385019142539297
5 inf
6 inf
7 inf
8 inf
9 inf
Code executed in 0.01 s
Selected k: 3 | True k: 3
====== Dunn ======
0 None
1 None
2 0.22248894960242496
3 0.09930825342423276
4 0.07694370083375403
5 0.06829859591450839
6 0.06533380025087833
7 0.06397229612167823
8 0.060820685569083446
9 0.06081222025210591
Code executed in 0.00 s
Selected k: 2 | True k: 3
====== XB ======
0 None
1 None
2 3821.6310955892613
3 5646.566350558779
4 3300.5590447734985
5 2945.622749207921
6 2932.6600403937264
7 2430.4419997127675
8 1743.770908279943
9 1732.6941019304875
Code executed in 0.00 s
Selected k: 9 | True k: 3
====== XB_star ======
0 None
1 None
2 9446.229005417188
3 13957.603039842357
4 8989.974686370095
5 8641.081615403018
6 8641.081615403018
7 7751.014168430341
8 5854.809633100446
9 5854.809633100446
Code executed in 0.00 s
Selected k: 8 | True k: 3
====== DB ======
0 None
1 None
2 inf
3 inf
4 inf
5 inf
6 inf
7 inf
8 inf
9 inf
Code executed in 0.01 s
Selected k: None | True k: 3
====== Inertia_sum ======
0 532.0548916172369
1 349.0214135130419
2 349.0095479285938
3 349.0536409924574
4 325.43055668310694
5 303.72575042981157
6 300.60978707111775
7 272.065674157918
8 233.85530390170882
9 231.56485922473567
Code executed in 0.00 s
Selected k: 1 | True k: 3
====== Diameter_max ======
0 6.180876027930203
1 4.737154354950527
2 4.737154354950527
3 4.737154354950527
4 4.737154354950527
5 4.737154354950527
6 4.737154354950527
7 4.737154354950527
8 4.665215171139693
9 4.665215171139693
Code executed in 0.00 s
Selected k: 1 | True k: 3
***** Time-series data using DTW with TimeSeriesKMeans *****
Clusterings generated in: 10374.47s
====== Hartigan_monotonous ======
0 -56.77317302493791
1 1754.8645616800438
2 76.24579045110623
3 60.777332178599764
4 1.3609113235168597
5 175.00215889526612
6 34.87676809542157
7 10.407030564341193
8 7.616910942373827
9 None
Code executed in 33.30 s
Selected k: 2 | True k: 4
====== CalinskiHarabasz_original ======
0 None
1 0.0
2 1754.085447852665
3 1433.539675966072
4 1234.1460305847866
5 1315.4680324836097
6 1594.3418842965932
7 1492.1951418605763
8 2816.6989164255124
9 2544.0901695781668
Code executed in 36.69 s
Selected k: 8 | True k: 4
====== GapStatistic_original ======
0 None
1 -0.8817177006460835
2 2.0566606194001658
3 2.636178054341171
4 3.1160429279552018
5 3.121160451261151
6 4.165457949134533
7 4.471067105499183
8 4.571000739007487
9 4.642315267340773
Code executed in 164.42 s
Selected k: 4 | True k: 4
====== Silhouette ======
0 None
1 None
2 0.970434561006521
3 0.8527254054088229
4 0.772910182768726
5 0.6711379707049308
6 0.7072500915064953
7 0.605201899730207
8 0.7085115341237531
9 0.694960801190566
Code executed in 22.38 s
Selected k: 2 | True k: 4
====== ScoreFunction ======
0 None
1 0.0
2 1.0
3 1.0
4 1.0
5 1.0
6 1.0
7 1.0
8 1.0
9 1.0
Code executed in 38.07 s
Selected k: 2 | True k: 4
====== MaulikBandyopadhyay_absolute ======
0 None
1 0.0
2 13570219.549174054
3 17417652.449314255
4 46382693.55690764
5 27149877.110046364
6 253876820.053277
7 404368207.58905816
8 422468110.5107191
9 416660003.3039842
Code executed in 24.38 s
Selected k: 8 | True k: 4
====== SD ======
0 None
1 None
2 1.2933347837618612
3 0.8723177381652063
4 0.9913874511437437
5 2.1480767082751564
6 0.8039022581611647
7 1.4093900478184977
8 1.2930901937226948
9 2.545568946545605
Code executed in 78.60 s
Selected k: 6 | True k: 4
====== SDbw ======
0 None
1 None
2 0.047171971054305324
3 0.028330098513020513
4 0.10369451212216116
5 inf
6 inf
7 0.06909185190063981
8 inf
9 inf
Code executed in 60.47 s
Selected k: 3 | True k: 4
====== Dunn ======
0 None
1 None
2 0.8305380457449938
3 0.004488337269329974
4 0.0006448396453246943
5 0.004343978650946443
6 0.015622169774173958
7 0.015389723232957044
8 0.02449231792844653
9 0.021862745625592615
Code executed in 11.94 s
Selected k: 2 | True k: 4
====== XB ======
0 None
1 None
2 0.013293449506085436
3 0.12378343790297791
4 0.4640240146988178
5 0.6923161663334289
6 0.2176150373185443
7 0.44324421603345265
8 0.19277183877019782
9 0.4204044080444315
Code executed in 22.40 s
Selected k: 2 | True k: 4
====== XB_star ======
0 None
1 None
2 0.024136782030059964
3 0.319786471287527
4 1.4427515208743282
5 1.6243286027215955
6 0.6571422606352396
7 1.6243286027215955
8 0.38969393633100585
9 0.889768850921729
Code executed in 25.27 s
Selected k: 2 | True k: 4
====== DB ======
0 None
1 None
2 inf
3 inf
4 inf
5 inf
6 inf
7 inf
8 inf
9 inf
Code executed in 26.91 s
Selected k: None | True k: 4
====== Inertia_sum ======
0 7619.754120376843
1 11346.458540247657
2 507.8773337544358
3 298.28457577945136
4 247.84334363584722
5 162.79113450833808
6 116.23199827162976
7 104.22438807761891
8 47.47200812717924
9 45.342736278713616
Code executed in 40.22 s
Selected k: 2 | True k: 4
====== Diameter_max ======
0 145.34014637603363
1 746.9674811658102
2 272.39277521031187
3 272.39277521031187
4 272.39277521031187
5 42.9663835289349
6 11.243614881589929
7 12.127901842859023
8 6.454483761851528
9 7.301726657697256
Code executed in 17.88 s
Selected k: 5 | True k: 4
***** Time-series data without DTW with KMeans *****
Clusterings generated in: 3.01s
====== Hartigan_monotonous ======
0 15.8094041815658
1 180.2820962759362
2 53.846852033119134
3 22.394820456074406
4 25.202760282907153
5 21.024527659437258
6 15.545274133443503
7 10.063484555853018
8 11.422390923636229
9 None
Code executed in 0.06 s
Selected k: 2 | True k: 4
====== CalinskiHarabasz_original ======
0 None
1 0.0
2 180.2820962759362
3 166.17448068089337
4 142.93732434587815
5 140.81535014558358
6 140.85515053229975
7 138.3283859855348
8 131.68560991187803
9 129.8496200914157
Code executed in 0.03 s
Selected k: 2 | True k: 4
====== GapStatistic_original ======
0 None
1 0.15179053114040642
2 1.1770774621210798
3 1.6057907369955204
4 1.7986187638345736
5 2.02284127008378
6 2.210532306790167
7 2.3475005692662902
8 2.43906866332496
9 2.543070698952932
Code executed in 23.61 s
Selected k: None | True k: 4
====== Silhouette ======
0 None
1 None
2 0.5851724720291331
3 0.5053323519853271
4 0.4602712659960852
5 0.45701901067637046
6 0.4165628388106481
7 0.42845540182332315
8 0.4114924154798477
9 0.4027101265250719
Code executed in 0.31 s
Selected k: 2 | True k: 4
====== ScoreFunction ======
0 None
1 6.445205555927203e-07
2 1.0
3 1.0
4 1.0
5 1.0
6 0.980164633078687
7 0.0013957461556722839
8 2.1928616544708568e-05
9 1.1352474449388694e-07
Code executed in 0.02 s
Selected k: 2 | True k: 4
====== MaulikBandyopadhyay_absolute ======
0 None
1 0.0
2 425.1118921667764
3 349.1043610495471
4 293.3275902208881
5 238.83141407717827
6 209.6927519929607
7 179.30535090398828
8 169.59916219362452
9 150.86196150274426
Code executed in 0.03 s
Selected k: 2 | True k: 4
====== SD ======
0 None
1 None
2 3.0428152641301534
3 3.18555999614733
4 2.3873900047715813
5 2.2756740410384326
6 2.4693582300323067
7 2.574320154251908
8 2.2589722360137667
9 2.3146390543079653
Code executed in 0.36 s
Selected k: 8 | True k: 4
====== SDbw ======
0 None
1 None
2 inf
3 inf
4 inf
5 inf
6 inf
7 inf
8 inf
9 inf
Code executed in 0.40 s
Selected k: None | True k: 4
====== Dunn ======
0 None
1 None
2 0.6502295327816004
3 0.4183330532684955
4 0.05453741885973551
5 0.05453741885973551
6 0.06321374469697073
7 0.07400800378311924
8 0.06494750239192468
9 0.06961462292209639
Code executed in 0.02 s
Selected k: 2 | True k: 4
====== XB ======
0 None
1 None
2 0.133838927267939
3 0.2070968818131891
4 0.5675473082510523
5 0.4475521107815224
6 0.36502976295588885
7 0.3122660522513699
8 0.4660908666110322
9 0.4135987883245314
Code executed in 0.01 s
Selected k: 2 | True k: 4
====== XB_star ======
0 None
1 None
2 0.23315333559033175
3 0.3264007671262032
4 1.105363220173563
5 1.105363220173563
6 0.8913374973888161
7 0.7152928875287857
8 1.1857212478302734
9 1.1857212478302734
Code executed in 0.01 s
Selected k: 2 | True k: 4
====== DB ======
0 None
1 None
2 inf
3 inf
4 inf
5 inf
6 inf
7 inf
8 inf
9 inf
Code executed in 0.02 s
Selected k: None | True k: 4
====== Inertia_sum ======
0 1546.7954242293042
1 1411.2211210804048
2 766.9830051666322
3 633.3155315412099
4 525.9256968224321
5 463.6934090080416
6 429.1543473123384
7 406.28201524121323
8 366.4045266679531
9 347.462941844373
Code executed in 0.01 s
Selected k: 2 | True k: 4
====== Diameter_max ======
0 26.59865134240156
1 29.63943217185524
2 26.085661705315456
3 21.104625863172615
4 21.104625863172615
5 21.104625863172615
6 18.207936044532552
7 15.55226140068361
8 15.55226140068361
9 14.509602898389424
Code executed in 0.02 s
Selected k: 3 | True k: 4