Binary $k$-Center with Missing Entries:
Structure Leads to Tractability

Clustering is a fundamental problem in computer science with a wide range of applications [28, 30, 49, 25, 51], and has been thoroughly explored [39, 3, 31, 9, 52, 4]. In its most general formulation, given $n$ data points, the aim of a clustering algorithm is to partition these points into groups, called clusters, based on the similarity. The degree of similarity between points is modeled by a given distance function. Depending on the representation of the data points, several computationally different variants of the clustering problem arise. Commonly, the points are embedded in the $d$-dimensional Euclidean space ${ℝ}^d$ with the standard Euclidean distance or the distance given by $L_p$ norms, or in the space of binary strings equipped with the Hamming distance, or in a general metric space with an explicit given distance function.

Additionally, there exist different approaches to how the similarity is aggregated. For the purpose of this work, we focus on the classical center-based clustering objectives. In $k$-Center clustering, given a set of data points and the parameter $k$, the objective is to partition the points into $k$ clusters and identify for each cluster a point called the center, so that the maximum distance between the center and any point within its cluster is minimized. $k$-Median clustering is defined in the same way, except that the sum of distances between the data points and their respective centers is minimized, and $k$-Means clustering aims to minimize the sum of squared distances instead.

Unfortunately, virtually all versions of clustering are computationally hard, in the classical sense. $k$-Means in ${ℝ}^d$ is NP-hard even on the plane (dimension $d=2$) [42], and it is also NP-hard for $k=2$ clusters even when the vectors have binary entries [2, 18]. $k$-Median and $k$-Center are NP-hard in ${ℝ}^2$ as well [44]. Moreover, $k$-Center on binary strings under the Hamming distance is NP-hard even for $k=1$ cluster; that is, the problem of finding the binary string that minimizes the maximum distance to a given collection of strings is already NP-hard [19, 37]. The latter problem is well-studied in the literature under the name of Closest String, given its important applications in coding theory [35] to bioinformatics [50].

In order to circumvent the general hardness results and simultaneously increase the modeling power of the problem, we consider the following variant of clustering with missing entries. We assume that the data points are represented by vectors, where each entry is either an element of the original domain (e.g., in $ℝ$ or $\{0,1\}$), or the special element “?”, which corresponds to an unknown entry. For the clustering objective, the distance is computed normally between the known entries, but the distance to an unknown entry is always zero. In this way, we can define the problems $k$-Center with Missing Entries and $k$-Means with Missing Entries. For formal definitions, we refer to Section 2.

These problems have a wide range of applications. For an example in predictive analytics, consider the setting of the well-known Netflix Prize challenge¹¹1The problem description and the dataset is available at https://www.kaggle.com/netflix-inc/netflix-prize-data.. The input is a collection of user-movie ratings, and the task is to predict unknown ratings. The data can be represented in the matrix form, where the rows correspond to the users and the columns to the movies, and naturally most of the entries in this matrix would be unknown [46]. Clustering in this setting is then an important tool for grouping/labelling similar users or similar movies, based on the available data.

Clustering with missing entries is also closely related to string problems that arise in bioinformatics applications. In the fundamental genome phasing problem (known also as “haplotype assembly”), the input is a collection of reads, i.e., short subsequences of the two copies of the genome, and the task is to reconstruct both of the original sequences. Finding the best possible reconstruction in the presence of errors is then naturally modeled as an instance of $k$-Means with Missing Entries with $k=2$, where the data points correspond to the individual reads, using missing entries to mark the unknown parts of each read; the target centers represent the desired complete genomic sequences; and the clustering objective represents the total number of errors between the known reads and the desired sequences, which needs to be minimized. In fact, [47] use exactly this formalization of the phasing problem (under the name of “Weighted Minimum Error Correction”) as the algorithmic core of their WhatsHap phasing software. Note that in both examples above, the known entries lie in a finite, small domain. For the technical results in this work, we focus on vectors where the vectors have binary values, e.g., in $\{0,1\}$; the results however can be easily extended to the bounded domain setting.

Generally speaking, $k$-Center with Missing Entries and $k$-Means with Missing Entries cannot be easier than their fully-defined counterparts. While greatly increasing the modeling power of the problem, the introduction of missing entries poses also additional technical challenges. In particular, most of the methods developed for the standard, full-information versions of clustering are no longer applicable, since the space formed by vectors with missing entries is not necessarily metric: the distances may violate triangle inequality. On the positive side, one can observe that in practical applications the structure of the missing entries is not completely arbitrary. In particular, the known entries are often sparse – for example, in the above-mentioned Netflix Prize challenge, only about 1% of the user-movie pairs have a known rating [46]. Therefore, the “hard” cases coming from the standard fully-defined versions of $k$-Center/$k$-Means are quite far from the instances arising in applications of clustering with missing entries. This motivates the aim to identify tractable cases of $k$-Center with Missing Entries based on the structure of the missing entries, since the general hardness results for $k$-Center are not applicable in this setting. Formally, we use the framework of parameterized complexity in order to characterize such cases. We refer to standard textbooks on parameterized complexity for a more thorough introduction to the subject [11, 10].

In order to apply the existing machinery and to put the parameters we consider into perspective, we encode the arrangement of the missing entries into a graph. We say that the incidence graph of a given instance is the following bipartite graph: the vertices are the data points and the coordinates, and the edge between a point and a coordinate is present when the respective entry is known. When interpreting the input as a matrix, where the data points are the rows, replacing the known entries by “1” and missing entries by “0” results exactly in the biadjacency matrix of the incidence graph. We call this matrix the mask matrix of the instance, denoted by $𝑴$, and denote the incidence graph of the instance by $G_{𝑴}$; see Figure 1 for an illustration.

We mainly consider the following three fundamental sparsity parameters of the incidence graph $G_{𝑴}$:

• $\mathrm{■}$

  Vertex cover number $\mathrm{𝚟𝚌}(G_{𝑴})$: the smallest number of vertices that are necessary to cover all edges of the graph;

• $\mathrm{■}$

  Fracture number $\mathrm{𝚏𝚛}(G_{𝑴})$: the smallest number of vertices one needs to remove so that the connected components of the remaining graph are small, i.e., their size is bounded by the same number;

• $\mathrm{■}$

  Treewidth $\mathrm{𝚝𝚠}(G_{𝑴})$: the classical decomposition parameter measuring how “tree-like” the graph is.

See Figures 2 and 3 for a visual representation of instances with small vertex cover and fracture number, respectively. Treewidth is the most general parameter out of the three, and it encompasses a wide range of instances. For example, the main algorithmic ingredient in the WhatsHap genomic phasic software [47] is the FPT algorithm for $k$-Means with Missing Entries parameterized by the maximum number of known entries per column²²2Called coverage in their work., in the special case where $k=2$ and the known entries in each row form a continuous subinterval. In this setting, the treewidth of the incidence graph is never larger than this parameter. Vertex cover, on the other hand, is the most restrictive of the three; however, comparatively small vertex cover may be a feasible model in settings such as the Netflix Prize challenge, where most users would only have ratings for a relatively small collection of the most popular movies. The fracture number aims to generalize the setting of small vertex cover, to also allow for an arbitrary number of small local “information patches”, outside of the few “most popular” rows and columns. Note that fracture number is a strictly more general parameter than vertex cover number, since removing any vertex cover from the graph results in connected components of size one; that is, for any instance, $\mathrm{𝚏𝚛}(G_{𝑴})\le \mathrm{𝚟𝚌}(G_{𝑴})$. On the other hand, it holds that $\mathrm{𝚝𝚠}(G_{𝑴})\le 2\mathrm{𝚏𝚛}(G_{𝑴})$, see Section 6 for the details.

Previously, the perspective outlined above has been successfully applied to $k$-Means with Missing Entries. Ganian et al. [20] show that the problem admits an FPT algorithm parameterized by treewidth of the incidence graph in case of the bounded domain, as well as further FPT algorithms for real-valued vectors in more restrictive parametrization. However, similar questions for $k$-Center with Missing Entries remain widely open.

In this work, we aim to close this gap and investigate parameterized algorithms for the $k$-Center objective on classes of instances where the known entries are “sparse”, in the sense of the structural parameters above. Our motivation stems from the following. First, $k$-Center is a well-studied and widely applicable similarity objective, which in certain cases might be preferable over $k$-Means; specifically, whenever the cost of clustering is associated with each individual data point and has to be equally small, as opposed to minimizing the total, “social”, cost spread out over all data points. Second, $k$-Center is interesting from a theoretical perspective, being a natural optimization target that on the technical level behaves very differently from $k$-Means. Our findings, as described next, show that the $k$-Center objective is in fact more challenging than $k$-Means in this context, and we prove that $k$-Center with Missing Entries is as general as a wide class of integer linear programs.

In this section, we introduce key definitions and notations used throughout the paper. For an integer $n$ we write $[n]$ to denote the set $\{1,\mathrm{\dots },n\}$. We use ${\mathbb{Z}}_{+}$ to denote the set of non-negative integers. For a vector of real numbers $𝒗\in {ℝ}^{\mathrm{\ell }}$ with length $\mathrm{\ell }$, its $i$-th entry is denoted by $𝒗[i]$ and $v_i$ interchangeably. By ${𝒗}^T$, we denote the transpose of the vector. Similarly, for a binary string $s\in {\{0,1\}}^{\mathrm{\ell }}$ of length $\mathrm{\ell }$, its $i$-th bit is denoted by $s[i]$, and $s_i$ and we write $s=[s_1,\mathrm{\dots },s_{\mathrm{\ell }}]$ to denote the complete string. For a set $I\subseteq [\mathrm{\ell }]$, the substring of $s$ restricted to the indices in $I$ is written as $s[I]\in {\{0,1\}}^{|I|}$. For two indices $i,j\in [\mathrm{\ell }]$, the substring from $s[i]$ throughout $s[j]$ (inclusive) is denoted by $s[i,j]\in {\{0,1\}}^{j-i+1}$. The complement of $s$ is represented by $\overline{s}$. The Hamming distance between two binary strings $s_1,s_2\in {\{0,1\}}^{\mathrm{\ell }}$ is denoted by $d_H:{\{0,1\}}^{\mathrm{\ell }}\times {\{0,1\}}^{\mathrm{\ell }}\to {\mathbb{Z}}_{+}$ with $d_H(s_1,s_2)=\left|\{i\in [\mathrm{\ell }]:s_1[i]\ne s_2[i]\}\right|$. For a graph $G$, the set of vertices and edges are denoted by $V(G)$ and $E(G)$ respectively.
Consider a matrix $𝑨$ with $n$ rows and $m$ coordinates. The sets of rows and coordinates of $𝑨$ are denoted by $R_{𝑨}=\{r_1,r_2,\mathrm{\dots },r_n\}$ and $C_{𝑨}=\{c_1,c_2,\mathrm{\dots },c_m\}$ respectively. For $i\in [n]$ and $j\in [m]$, the $i$-th row of $𝑨$ is denoted by $𝑨[i]$ and the $j$-th entry of $𝑨[i]$ is written as $𝑨[i][j]$. We may also use $a_i$ and $a_{ij}$ respectively to refer to $𝑨[i]$ and $𝑨[i][j]$, when clear from context. For matrices $𝑨$ and $𝑩$, with the same number of rows and columns, their entry-wise subtraction is written as $𝑨-𝑩$, meaning $(𝑨-𝑩)[i][j]=𝑨[i][j]-𝑩[i][j]$. Similarly, their entry-wise product is denoted by $𝑨\circ 𝑩$, where $(𝑨\circ 𝑩)[i][j]=𝑨[i][j]\cdot 𝑩[i][j]$. For two binary matrices $𝑨$ and $𝑩$, we define their row-wise Hamming distance as a column vector $d_H(𝑨,𝑩)$ of length $n$, where, $d_H(𝑨,𝑩)\left[i\right]=d_H(𝑨[i],𝑩[i])$. For two vectors $𝒂$ and $𝒃$ of the same length, $𝒂\le 𝒃$ if for all coordinate $i$, $𝒂[i]\le 𝒃[i]$.

Moreover, a missing entry is denoted by “?”. We define the Hamming distance between a known entry (0 or 1) and an unknown entry (?) as zero, i.e., $d_H(0,\mathrm{?})=d_H(1,\mathrm{?})=0$. We remark here that assuming the distance to a missing entry to be any other constant (e.g., 1 instead of 0) does not introduce additional complexity, as the number of missing entries in any given row is fixed by the input and their contribution to the cost of clustering does not depend on the solution. It is easy to observe that all the algorithms and techniques presented here would still apply in a straightforward manner to the problem with arbitrary cost of a missing entry, by adjusting the target radius of each row. For clarity of the exposition and consistency with existing literature on missing entries, we define the cost of a missing entry to be zero, as above.

For a binary matrix $𝑴\in {\{0,1\}}^{n\times m}$, its incidence graph is an undirected bipartite graph defined as $G_{𝑴}=(V_{𝑴},E_{𝑴})$, where $V_{𝑴}=(R_{𝑴}\cup C_{𝑴})$ and $E_{𝑴}=\{(u,v):u\in R_{𝑴},v\in C_{𝑴},𝑴[u][v]=1\}$. We use row (resp. coordinate) vertices to represent the vertices in $G_{𝑴}$ that correspond to the rows (coordinates) of $𝑴$. For a visual representation of the incidence graph, refer to Figure 1.

With these definitions, we can now proceed to formalize the main problem in question. Note that we assume the input matrix $𝑨$ to only contain $\{0,1\}$ entries, as the missing entries are encoded by the matrix $𝑴$, and when $𝑴[i][j]=0$, its product with $𝑨[i][j]$ is always 0.

Note that in the context of a $k$-Center with Missing Entries instance, we use the term “point” interchangeably with “row” and “coordinate” interchangeably with “column”, since the points to be clustered are represented as rows of the matrices $𝑨$ and $𝑴$. We also say that an entry of $𝑨$ is a missing entry or “?” whenever the corresponding entry of $𝑴$ is zero.

Next, we define the related problems.

Finally, we recall the known algorithms for Closest String/ILP Feasibility with respect to the number of strings/number of constraints. Knop et al. [33] designed a general approach to solve ILP formulations of various problems, including Closest String and Non-uniform Closest String, which is also applicable to versions with missing entries.

For ILP Feasibility where the coefficients are bounded, Eisenbrand and Weismantel [16] show a fast dynamic programming algorithm that has a similar running time.

This section is dedicated to the reduction from ILP Feasibility to Closest String that preserves the number of rows up to a factor of $𝒪\left(\log {\parallel 𝑨\parallel }_{\mathrm{\infty }}\right)$, where $𝑨$ is the constraint coefficient matrix in the ILP Feasibility instance. We show the sequence of reductions that lead to the proof of Theorem 11, which is then used to establish Theorem 4. While complete proofs of Lemma 12 to Lemma 18 are provided in the full version, we state the results formally next.

We start with a lemma showing that a general ILP Feasibility instance, where the coefficients may be negative, and variables may have negative lower bounds, can be transformed, in linear time with respect to the number of variables, into an equivalent instance with non-negative variable domains, i.e., of form $0\le y_i\le u_i^{\prime }$, non-negative coefficient matrix and non-negative right-hand side, while mostly preserving the number of constraints, variables and the magnitude of the entries, by a constant factor.

By Lemma 12, we observe that a general ILP Feasibility instance is equivalent to an instance where everything is non-negative, with a bounded blow-up in the number of constraints, variables and the magnitude of values. Thus, we continue our chain of reductions from a non-negative ILP. In the next step, we reduce from arbitrary bounds on the variables to binary variables, as stated in the following lemma.

Next, we show that an instance of ILP Feasibility with non-negative coefficients and binary variables can be reduced to an instance with binary coefficients. The proof is similar in spirit to the proof of Lemma 8 in [34], however we are interested in ILP Feasibility instances where variables are only allowed to take values in $\{0,1\}$, which is the main technical difference.

Now we introduce a lemma that allows to reduce a $\{0,1\}$-coefficient matrix to a matrix with values only in $\{-1,1\}$.

It will be more convenient to work with ILP instances of form $𝑨𝒙\le 𝒃$, which motivates the following lemma.

Before reducing to Closest String itself, we first show a reduction to a slightly more general version of the problem, Non-uniform Closest String (Definition 7), where each string has its own upper bound on the distance to the closest string.

Finally, we show a reduction from Non-uniform Closest String to Closest String.

With the reductions above, we can conclude with the proof of Theorem 11 as outlined below.

In this section, we describe an FPT algorithm for $k$-Center with Missing Entries parameterized by the vertex cover number of the incidence graph, $G_{𝑴}$, derived from mask matrix $𝑴$. More formally, we prove the following: See 1

This section is dedicated to a fixed-parameter algorithm for the problem parameterized by the treewidth of the incidence graph $G_{𝑴}$, the number of clusters $k$ and the maximum permissible radius of the cluster, $d$. We restate the formal result next for convenience.

See 3 We briefly sketch the intuition of our approach. For a complete proof and discussion, refer to the full version. Informally, the fact that $G_{𝑴}$ has treewidth at most $t$ means that the graph can be constructed in a tree-like fashion, where at each point only a vertex subset of size at most $t$ is “active”. This small subset is called a bag, and in particular is a separator for the graph: there are no edges between the “past”, already constructed part of the graph and the “future”, not encountered yet part of the graph; all connections between these parts are via the bag itself. In terms of the $k$-Center with Missing Entries instance, this means that for each “past” row, all its entries that correspond to “future” columns are missing (as the respective entry of the mask matrix $𝑴$ has to be $0$), and the same holds for “past” columns and “future” rows.

Our algorithm performs dynamic programming over this decomposition, supporting a collection of records for the current bag. For the rows of the current bag, we store the partition of the rows into clusters, and additionally for each row the distance to its center vector among the already encountered columns. For each cluster and for each column of the bag, we store the value of the cluster center in this column. These three characteristics (called together a fragment) act as a “trace” of a potential solution on the current bag. In our DP table, we store whether there exists a partial solution for each choice of the fragment. Because of the separation property explained above, only knowing the fragments is sufficient for computing the records for every possible update on the bag. The claimed running time follows from upper-bounding the number of possible fragments for a bag of size at most $t$.

In this section, we present a fixed-parameter algorithm for the $k$-Center with Missing Entries, where the parameter is the fracture number of the incidence graph $G_{𝑴}$. To provide context, we first define the key concepts of fracture modulator and fracture number. Following this, we exploit the relevant results from Theorem 1 and Theorem 3, which lead to the proof of Theorem 2. See 2

For a graph $G=(V,E)$, a fracture modulator is a subset of vertices $F\subseteq V$ such that, after removing the vertices in $F$, each remaining connected component contains at most $|F|$ vertices; the size of the smallest fracture modulator is denoted by $\mathrm{𝚏𝚛}(G)$. Consider the incidence graph $G_{𝑴}$ corresponding to the mask matrix $𝑴$, and let $F$ represent the fracture modulator of $G_{𝑴}$ with the smallest cardinality. We note that $F$ can be computed in time $𝒪\left({(\mathrm{𝚏𝚛}(G)+1)}^{\mathrm{𝚏𝚛}(G)}nm\right)$ by the algorithm of [12]. Let $R_F=(R_{𝑴}\cap F)$ and $C_F=(C_{𝑴}\cap F)$ represent the row and column vertices of $G_{𝑴}$ that belong to the fracture modulator $F$, respectively. Additionally, define $\overline{R_F}=R_{𝑴}\setminus R_F$ and $\overline{C_F}=C_{𝑴}\setminus C_F$. We refer to the rows in $R_F$ as long rows while the rows corresponding to $\overline{R_F}$ are called short rows. See Figure 3 for an illustration of the structure of the instance. Now, the algorithm considers two cases.

We have investigated the algorithmic complexity of $k$-Center with Missing Entries in the setting where the missing entries are sparse and exhibit certain graph-theoretic structure. We have shown that the problem is FPT when parameterized by $\mathrm{𝚟𝚌}+k$, $\mathrm{𝚏𝚛}+k$ and $\mathrm{𝚝𝚠}+k+d$, where $\mathrm{𝚟𝚌}$ is the vertex cover number, $\mathrm{𝚏𝚛}$ is the fracture number, and $\mathrm{𝚝𝚠}$ is the treewidth of the incidence graph. In fact, it is not hard to get rid of $k$ in the parameter; for example, in the enumeration of partial center assignments in the algorithm of Theorem 1 it can be assumed that $k\le 2^{\mathrm{𝚟𝚌}}$, as multiple centers that have the same partial assignment are redundant. However, this would increase the running time as a function of the parameter to doubly exponential, therefore we state the upper bounds with explicit dependence on $k$. It is, on the other hand, an interesting open question, whether $d$ in the parameter is necessary for the algorithm parameterized by $\mathrm{𝚝𝚠}$. As shown by [20], this is not necessary for $k$-Means with Missing Entries, since the problem admits an FPT algorithm when parameterized by treewidth alone. Yet, it does not seem that the dynamic programming approaches used in their work and in our work, can be improved to avoid the factor of $d^{𝒪\left(\mathrm{𝚝𝚠}\right)}$ in the case of $k$-Center with Missing Entries. Therefore, it is natural to ask whether it can be shown that $k$-Center with Missing Entries is W[1]-hard in this parameterization.

Another intriguing open question is the tightness of our algorithm in the parameterization by the vertex cover number (and fracture number). On the one hand, the running time we show is $2^{{\mathrm{𝚟𝚌}}^{2+o(1)}}\cdot \mathrm{poly}(nm)$ (for values of $k$ in $𝒪\left(vc\right)$), exceeding the “natural” single-exponential in $\mathrm{𝚟𝚌}$ time, and we are not aware of a matching lower bound that is based on standard complexity assumptions. On the other hand, we show that improving this running time improves also the best-known running time for ILP Feasibility when parameterized by the number of constraints and Closest String parameterized by the number of strings. The latter are major open questions; [48] show also that several other open problems are equivalent to these. On the positive side, our algorithm in fact reduces $k$-Center with Missing Entries to just $2^{\mathrm{𝚟𝚌}}$ instances of ILP Feasibility. Since practical ILP solvers are quite efficient, coupling our reduction with an ILP solver is likely to result in the running time that is much more efficient than prescribed by the upper bound of Theorem 1.