EKM: An exact, polynomial-time algorithm for the $K$-medoids problem

TL;DR

This work addresses the exact K-medoids clustering problem by introducing EKM, an algorithm with worst-case time $O(N^{K+1})$ derived using the Bird-Meertens formalism and shortcut fusion to guarantee global optimality. EKM replaces traditional exponential-time approaches with a provably correct, recursive combinatorial generator that fuses evaluation and selection into the recursion, enabling polynomial-time worst-case analysis and parallelizable implementation. Empirical results on real-world and synthetic data demonstrate that EKM achieves the exact optimum and scales to datasets up to $N=5000$, often outperforming MIP-based solvers in wall-clock time. The work also discusses memory considerations, potential constraints handling via semiring lifting, and outlines directions for further parallel optimization and practical enhancements.

Abstract

The $K$-medoids problem is a challenging combinatorial clustering task, widely used in data analysis applications. While numerous algorithms have been proposed to solve this problem, none of these are able to obtain an exact (globally optimal) solution for the problem in polynomial time. In this paper, we present EKM: a novel algorithm for solving this problem exactly with worst-case $O\left(N^{K+1}\right)$ time complexity. EKM is developed according to recent advances in transformational programming and combinatorial generation, using formal program derivation steps. The derived algorithm is provably correct by construction. We demonstrate the effectiveness of our algorithm by comparing it against various approximate methods on numerous real-world datasets. We show that the wall-clock run time of our algorithm matches the worst-case time complexity analysis on synthetic datasets, clearly outperforming the exponential time complexity of benchmark branch-and-bound based MIP solvers. To our knowledge, this is the first, rigorously-proven polynomial time, practical algorithm for this ubiquitous problem.

Figures (2)

• Figure 1: A generic recursive generator of all possible sublists for list $\left[\boldsymbol{x}_{1},\boldsymbol{x}_{2},\boldsymbol{x}_{3},\boldsymbol{x}_{4}\right]$, in each recursive stage $n$, the sublists of size $k$ (combinations) are automatically grouped into list $\mathcal{S}_{k}^{n}$, for all $k\in\left\{ 0,1,\ldots n\right\}$.
• Figure 2: Log-log wall-clock run time (seconds) for our algorithm (EKM) tested on synthetic datasets (left panel). The run-time curves from left to right (corresponding to $K=2,3,4,5$ respectively), have slopes 3.005, 4.006, 5.018, and 5.995, an excellent match to the predicted worst-case run-time complexity of $O\left(N^{3}\right)$, $O\left(N^{4}\right)$, $O\left(N^{5}\right)$, and $O\left(N^{6}\right)$ respectively. Log-linear wall-clock run-time (seconds) comparing EKM algorithm against a classical MIP (BnB) solver (GLPK) on synthetic datasets with $K=3$ (right panel). On this log-linear scale, exponential run-time appears as a linear function of problem size $N$, whereas polynomial run-time is a logarithmic function of $N$.