<commit_before>\documentclass[conference, 10pt, final]{IEEEtran}

\usepackage{cite}
\usepackage[cmex10]{amsmath}
\usepackage{amssymb}
\usepackage[pdftex]{graphicx}
\usepackage{algorithmic}
\usepackage{array}
\usepackage[caption=false,font=footnotesize]{subfig}
\usepackage{stfloats}
\usepackage{url}


\begin{document} 

\title{Improved Hierarchical Clustering Through Nesting}


\author{\IEEEauthorblockN{<NAME>, <NAME>, and <NAME>}
\IEEEauthorblockA{Operations Research and Information Engineering, Cornell University, Ithaca, NY.  \\
Email: \{ ajp336, cjw278, dal299 \} @<EMAIL>}}


\maketitle

\begin{abstract}
Hierarchical clustering organizes large data into smaller, more informative groups and offers a global view of the data structure by producing a hierarchy of clusters starting from individual data points and merging towards the cluster of all data. Typically, hierarchical clustering is done by greedily either merging or dividing clusters, called agglomerative and divisive clustering, respectively. Lin et al.~\cite{Lin} proposed a nesting framework that produces a hierarchal clustering from partitional clusterings for every possible number of clusters. This algorithm was introduced as an approximation algorithm for the hierarchical $k$-median problem and guarantees the resulting hierarchical layers are a limited distance from the original partitional clusters. However, this algorithm has not been previously placed in a clustering context. We test this algorithm using various methods for producing the partitional clusterings against the traditional hierarchical clustering methods and show both the approximation ratio and FScore measure improve significantly. We also propose a heuristic for this nesting framework which improves the run time and only suffers a slight decrease in approximation and FScore.
\end{abstract}


\section{Introduction and Notation}
As the data we have access to grows in variety and size and becomes cheaper to store, we need methods to efficiently learn from this data. 
Partitioning data (clustering) into similar groups (clusters) is one of the key tools to learn about a set of unlabeled data such as a set of customers who have visited a website or a set of handwritten documents.
Clustering in these situations helps us learn underlying structure and classifications of our data, and it allows large data sets to be presented in a simplified and compressed view \cite{Jain}.

There are two main branches of clustering: partitional (or flat) and hierarchical. In partitional clustering, we often specify a number $k$ and pairwise distances between the data points, representing a measure of dissimilarity such as Euclidean distance. The clustering algorithm then partitions the data into $k$ disjoint sets where the sets represent relatively similar data. In hierarchical clustering, on the other hand, we are interested in finding a nested tree of partitions. This offers the user a global view of the data and allows the user to choose a ``natural'' level of clustering \cite{ESL}. For example, a retailer could use a hierarchical clustering of customers to first find regions of customers in order to tailor the products offered in each region. The retailer could then also offer specialized discounts to subgroups in each region based on each cluster's shopping records. Examples of both types of clustering are given in Fig.~\ref{fig_part}-\ref{fig_hier}.

Clustering algorithms date back more than 50 years, and there have been many variations introduced since then \cite{Jain}. However, clustering algorithms are often formulated as NP-hard problems with non-convex objective functions, which are inherently hard to solve \cite{Jain}.
One common objective function for partitional clustering is to minimize the sum of distances from each point to the median point of its cluster, called $k$-median or $k$-medoids clustering. More specifically, suppose we are given a set $X$ and a distance (dissimilarity) metric $d(i,j)$ for all $i,j \in X$. The goal is to find a set of centers $S \subseteq X$ such that $|S| \leq k$ so as to minimize 
\[ f(S) = \sum_{ i\in X} \min_{j \in S} d(i,j). \]
The clustering is obtained from the set of centers $S$ by assigning each data point to its closest center.  An $\alpha$-approximation algorithm for this problem produces a set $S$ in polynomial time such that $f(S) \leq \alpha \cdot f(O)$ for all $O \subseteq X$ with $|O| = k$. 

\begin{figure*}[!]
\centerline{\subfloat[Partitional Clustering]{\includegraphics[width =1.75in]{partitional_clustering}
\label{fig_part}}
\hfil
\subfloat[Hierarchical Clustering]{\includegraphics[width=2.5in]{hier_clustering}
\label{fig_hier}}}
\caption{A partitional clustering for $k=3$, where the colors represent the different partitions, and a dendrogram representation of a hierarchical clustering of the data where edges represent merges of clusters and the clustering on the third level down (represented by the red line) corresponds to the partitional clustering on the left.}
\end{figure*}

This objective function can then be extended to the hierarchical $k$-median problem which represents an objective function relating to the hierarchical clustering problem. Suppose we represent a solution to the $k$-median problem as a pair $(S, a)$ where $S \subseteq X$ such that $|S| \leq k$ and $a$ is an assignment function that maps each $i \in X$ to a center in $S$. Then, for any $i \in S$, $a^{-1}(i)$ is the cluster centered at $i$. We define $f(S,a) := \sum_{i \in X} d(i, a(i))$. This is the same objective function as above except the assignment is given rather than implied.

Given solution pairs $(S_1, a_1)$ and $(S_2, a_2)$, we say that $(S_1, a_1)$ is \emph{nested} in $(S_2, a_2)$ (denoted $(S_1, a_1) \subseteq (S_2, a_2)$) if the following conditions hold:
\begin{enumerate}[\IEEEsetlabelwidth{3)}]
	\item $S_1 \subset S_2$, 
	\item $\forall  j \in X$, if $a_2(j) \in S_1$, then $a_1(j) = a_2(j)$, 
	\item $\forall j, k \in X$, if $a_2(j) = a_2(k)$, then $a_1(j) = a_1(k)$.
\end{enumerate}
A solution to the hierarchical $k$-median problem is a nested chain of solutions, and the hierarchical clustering it generates is found by looking at the clusters given by each $(S_i, a_i)$. An example is given in Fig.~\ref{fig_chain1}-\ref{fig_chain4}. However, we want this solution to stay close to the optimal $k$-median solutions when looking at any level of the hierarchy. An $\alpha$-approximation algorithm for this problem generates a nested chain of solutions in polynomial time such that for all $k = 1, \ldots, |X|$ there exists $(S_i, a_i)$ in the chain such that $|S_i| \leq k$ and $f(S_i, a_i) \leq \alpha \cdot f(O)$ for all $O \subseteq X$ with $|O| = k$.

\begin{figure*}[b]
\centering
\subfloat[$(S_1, a_1)$]{\includegraphics[width =1.75in]{chain1}
\label{fig_chain1}} \hspace{15mm}
\subfloat[$(S_2, a_2)$]{\includegraphics[width=1.75in]{chain2}
\label{fig_chain2}} \\
\subfloat[$(S_3, a_3)$]{\includegraphics[width=1.75in]{chain3}
\label{fig_chain3}}  \hspace{15mm}
\subfloat[$(S_4, a_4)$]{\includegraphics[width=1.75in]{chain4}
\label{fig_chain4}}
\caption{An example of a nested chain of solutions. $S_i$ is indicated by the colored nodes and $a_i$ by the edges. The clusters on that level are given by the connected components. Note that $S_{i-1}$ is a subset of $S_i$ and the assignments stay consistent. }
\label{fig_sim}
\end{figure*}

While some partitional clustering methods use the $k$-median objective function to inform their decisions, popular methods for hierarchical clustering do not and rely on greedy algorithms which have no guarantee on the solution produced with respect to the objective function but are typically fast and easy to implement. This paper looks at how using this objective function can vastly improve the hierarchical clusterings produced. 
In particular, we show how an approximation framework for this problem, which takes in partitional clusterings for all values of $k$, compares in practice to some of the most widely used clustering algorithms, both in the quality of the solution and in computation time. The algorithms we use to produce these partitional clusterings include both approximation algorithms for the $k$-median problem, which offer an overall approximation guarantee for the hierarchical $k$-median problem, and some heuristic algorithms that have been shown to perform well in practice. 

\section{Literature Review}

\subsection{Partitional Clustering}
Arguably the most common partitional clustering method is $k$-means, which chooses $k$ cluster centers points in space to minimize the sum of squared distances from each point to its closest cluster center \cite{Jain}.  $k$-means is usually solved using Lloyd's algorithm (also referred to simply as the $k$-means algorithm), which iteratively greedily updates the centers of each cluster (assigning points to the closest center) until the solution becomes stable \cite{Jain}. Lloyd's algorithm is only guaranteed to reach a locally optimum solution and can be very sensitive to outliers and the initial points chosen \cite{Kanungo}.

The $k$-median objective function removes some worst-case behavior of $k$-means by restricting the cluster centers to be data points.  There are two common algorithms for this problem. The first is a local search algorithm called partition around medoids (PAM), which starts with $k$ randomly chosen data points as centers, and looks to swap one of the centers with another data point to improve the objective value. This repeats until no more improving swaps exist \cite{ESL}. The generalization of this algorithm which considers swaps of size $\leq p$ is a $(3 + \frac{2}{p})$-approximation algorithm \cite{Arya}.
Another common heuristic (called $k$-medoids throughout this paper), which is more computationally efficient than PAM \cite{Park}, is a variation on Lloyd's $k$-means algorithm.  Starting with randomly chosen centers, we assign each data point to its closest center, and then update the center within each cluster \cite{Park}. The algorithm repeats this until the solution becomes stable.

The $k$-median problem is also a fundamental problem in approximation algorithms where it is often presented in the context of opening service centers or factories to serve clients. However, these algorithms are often not implemented for clustering despite their improved guarantees on the solution. The best known approximation algorithm for this problem has an approximation guarantee of 2.592 and is produced using linear programming \cite{Wu}. A simpler LP rounding algorithm by Charikar and Li has a slightly worse guarantee of 3.25 but is easier to implement in practice \cite{Charikar}. 

\subsection{Hierarchical Clustering}
The most common method for hierarchical clustering is agglomerative clustering. Agglomerative clustering starts with each data point in its own cluster, and proceeds to recursively merge the two closest clusters to form the next clustering level in the hierarchy, repeating until the entire data set is a single cluster.  There are several choices on how to measure distance in order to choose which clusters to merge. Single linkage (also called nearest neighbor) defines the distance between two clusters to be the minimum pairwise distance between the data points in each cluster.  Complete linkage clustering instead defines the distance to be the largest pairwise distance. A method between these two extremes is group average linkage which looks at the average pairwise distance between clusters. This strikes a balance between the closeness and compactness of the two techniques above and gives more reasonable results in practice \cite{ESL}. Lastly, there is Ward's method which chooses the two clusters to merge so as to minimize the sum of each cluster's variance \cite{Ward}.

Approximation algorithms for hierarchical clustering have also been studied, where we aim to identify a hierarchy of clusters so as to minimize the ratio of the objective value on every level in the hierarchy of size $k$ to the optimal objective value of the partitional clustering problem with $k$ clusters. Note that it may not always be possible to find a solution with approximation ratio 1. The hierarchical $k$-median problem was originally introduced by Plaxton who gave a 238.88 approximation algorithm \cite{Plaxton}. This was improved by Lin et al.~\cite{Lin} who gave an incremental approximation algorithm framework. In particular, their approach relies on a black box $k$-median solver and then uses these solutions for $k=1, \ldots, n$ to greedily construct a hierarchical clustering that increases the objective value by at most a factor of $20.71$ (or $10.03$ when randomized) at each level. Using the local $p$-swap $(3+\frac{2}{p})$-approximation algorithm mentioned above as the black-box $k$-median solver performed well in experiments on hierarchical $k$-median test problems \cite{Nagarajan}. While the authors tested the variants of this algorithm for different $k$-median approximation algorithms, it has not been tested against traditional clustering algorithms or on any labeled data sets. 

Several algorithms have been built upon agglomerative clustering that we choose not to compare with our proposed algorithms. Divisive clustering, similar to agglomerative, starts with the data in one cluster and repeatedly splits a cluster on each level. This method has not been studied nearly as much as agglomerative clustering \cite{ESL} but has had relative success in classifying document datasets \cite{Zhao}.  Algorithms such as the CURE, ROCK, and Chameleon \cite{Guha_CURE, Guha_ROCK, Karypis} incorporate non-spherical cluster shapes and cluster density, but do not relate as much to the above approximation algorithm in their objective. We also choose not to compare our algorithms with Bayesian hierarchical methods since these methods require many more parameters and add additional prior information that our algorithms do not use. 

\section{Proposed Algorithms}

\subsection{Hierarchical Algorithm} 

Recall the hierarchical $k$-median problem: given a set $X$ and a distance metric $d(i,j)$ for all $i,j \in X$, we represent a solution to the $k$-median problem as a pair $(S, a)$ where $S \subseteq X$ such that $|S| \leq k$ and $a$ is an assignment function that maps each $i \in X$ to a center in $S$.
 A solution to the hierarchical $k$-median problem is a chain of nested solutions $(S_1, a_1) \subseteq (S_2, a_2) \subseteq \ldots \subseteq (S_l, a_l)$. However, we want this solution to stay close to the optimal $k$-median solutions when looking at any level of the hierarchy.  
Consider any $k$-median algorithm $\mathcal{A}$. After running $\mathcal{A}$ for all potential values of $k$, Lin et al.~\cite{Lin} introduce an augmentation step that nests these clusterings into the hierarchical framework and only increase costs by a bounded factor. 

First, set $\beta = 3+\sqrt{3}$, and $\beta_0 = 1$ (for the randomized algorithm we switch $\beta = 6.355$ and set $\beta_0 = \beta^X$, where $X$ is chosen from a uniform distribution on $[0,1]$). Let $V_i$ be the solutions generated from running $\mathcal{A}$ for $k=n-i+1$ for $i=1, \ldots, n$. The cost of these solutions is determined by assigning each point to the nearest cluster center and is scaled so that the minimum cost for $k<n$ is at least one. We then bucket these solutions by cost into buckets of the form $[0,0], (\beta_0, \beta_0 \beta], (\beta_0 \beta, \beta_0 \beta^2], \ldots $. 

From each non-empty bucket, pick the solution with the smallest number of cluster centers, and let these solutions be $\overline{V}_1, \overline{V}_2, \ldots \overline{V}_r$. We then inductively build our nested solutions: Set $S_1 = \overline{V}_1$ and $a_1$ to be the assignment that assigns each point to the closest point in $S_1$. Given $(S_{i-1}, a_{i-1})$ and $\overline{V}_i$, we construct $S_i$ by  adding the closest point $j \in S_{i-1}$ to $l$ for each $l \in \overline{V}_i$. Note that $|S_{i}| \leq |\overline{V}_i|$. 

\begin{figure*}[!t]
\centering
\subfloat[Nested solution $(S_{i-1}, a_{i-1})$ ]{\includegraphics[width =1.75in]{Si-1}
\label{fig_Si-1}} \hspace{15mm}
\subfloat[Next set $\bar{V}_i$]{\includegraphics[width=1.75in]{Vi}
\label{fig_Vi}} \\
\subfloat[Next nested centers $S_i$]{\includegraphics[width=1.75in]{Si}
\label{fig_Si}}  \hspace{15mm}
\subfloat[Next nested solution $(S_i, a_i)$]{\includegraphics[width=1.75in]{Siai}
\label{fig_Siai}}
\caption{An example of nesting $\bar{V}_i$ with $(S_{i-1}, a_{i-1})$. Dark  nodes represent chosen centers and edges represent assignments to that center's cluster. First, node 2 chooses node from 1 from $S_{i-1}$ as its closest center and node 5 chooses node 4. Then, each cluster in $(S_{i-1}, a_{i-1})$ is assigned to the cheapest (lowest sum of distances) center.}
\label{fig_sim}
\end{figure*}

We then need to define an assignment $a_i$. For each $j \in X$ such that $a_{i-1}(j) \in S_{i}$, we set $a_i(j) = a_{i-1}(j)$. Otherwise, we set $a_i(j)$ to be the facility $i \in S_i$ such that assigning all points in $a_{i-1}^{-1}(j)$ (i.e. all points in $j$'s cluster in $S_{i-1}$) would have minimum cost. Using this assignment, $(a_i, S_i)$ is nested in $(a_{i-1}, S_{i-1})$. An example of this nesting is given in Fig.~\ref{fig_Si-1}-\ref{fig_Siai}.

We continue this nesting until $|S_i| = 1$, reaching the top level of the hierarchy where all data points are in the same cluster. After forming this hierarchical clustering, for any $k \in 1, \ldots, n$ \cite{Lin} guarantees that the cost of $S_i$ is at most $20.71$ times the cost of $V_k$ where $|S_{i-1}| \geq |V_k| \geq |S_{i}|$ (in the randomized case the expected cost of $S_i$ is at most $10.03$ times the cost of $V_k$).


\subsection{Partitional Clustering Algorithms}
This nesting framework allows us to choose the partitional clustering algorithm $\mathcal{A}$. As mentioned above, PAM and $k$-medoids are two potential algorithms. PAM has the added benefit of being a 5-approximation algorithm for the partitional $k$-median problem which yields to an overall approximation ratio of $20.71 \cdot 5 = 103.55$ for the hierarchical $k$-median problem ($10.03 \cdot 5 = 50.15$ randomized). We looked at increasing the size of the local swaps in PAM to improve the approximation ratio but found that the approximation ratio of constructed solutions for given values of $k$ did not significantly improve while the computation time increased greatly. 

We also consider three linear programming based algorithms, one of which has an approximation guarantee. We first present the linear programming formulation for the $k$-median problem. Suppose we have a set $X$. Let $s_i \in \{0,1\}$ be the variable representing whether or not we want to make $i$ a center for all $i \in X$. We also create a variable $a_{i,j} \in \{0, 1\}$ which represents whether or not we assign $i \in X$ to the cluster with center $j \in X$. Given this set and assignment, the $k$-median objective function is 
\begin{equation}
\sum_{i \in X} \sum_{j \in X} a_{i,j} d(i,j), 
\end{equation}
which we want to minimize. However, we are restricted by the fact that we can choose at most $k$ centers, every point must be assigned to a single cluster center, and that we can only assign a point to a center that is actually in $S$. This is represented by the following three equations, respectively:
\begin{eqnarray}
\sum_{i \in X} s_i   \leq k & ,\\
\sum_{j \in X} a_{i,j}  = 1 & \forall i \in X, \\
a_{i,j}   \leq s_j & \forall i, j \in X.
\end{eqnarray}
Therefore, solving the following integer program would directly find the optimal set of centers and assignment for the $k$-median problem.
\begin{eqnarray}
\text{Min } & \sum_{i \in X} \sum_{j \in X} a_{i,j} d(i,j) \\
\text{s.t. } & \sum_{i \in X} s_i   \leq k  \\
& \sum_{j \in X} a_{i,j}  = 1 \quad \forall i \in X \\
& a_{i,j}   \leq s_j \quad \forall i, j \in X \\
& s_i \in \{0,1\} \quad \forall i \in X \\
& a_{i,j} \in \{0,1\} \quad \forall i, j \in X
\end{eqnarray}
Relaxing the constraints so that $s_i \in [0,1]$ for all $i \in X$ and $a_{i,j} \in [0,1]$ for all $i,j \in X$ yields the linear programming (LP) relaxation for the $k$-median problem, which we can solve in polynomial time.

Interestingly, the solution to the linear programming is often integer valued \cite{Nagarajan}. Based on this fact, we propose two LP-based clustering algorithms: one that solves the LP relaxation and constructs $S$ by taking the points in $X$ with the $k$ highest $s_i$ values and one that solves the LP relaxation and constructs $S$ by choosing $i$ to be in $S$ with probability $s_i$, sampling for each $i \in X$ until $k$ centers are chosen. We refer to these algorithms as LP-greedy and LP-rand, respectively.
Charikar and Li \cite{Charikar} introduce a 3.25-approximation algorithm that uses this LP solution to smartly construct a set of centers. We refer to this algorithm as LP-Charikar. We choose to implement this LP algorithm rather than the LP-based approximation algorithm with best approximation ratio because the algorithm is simpler to implement (and we expect it to have a lower running time). 

The algorithm works by first choosing a subset of data points $C \subseteq X$ to represent the demand. For each data point $i$, let $d_{av}(i) = \sum_{j \in X : a_{i,j} > 0} s_j \cdot d(i,j)$. We start with $C = \emptyset$. At each step, we choose a point $j \in X$ to add to $C$, and all points $j'$ such that $d(j,j') \leq 4 d_{av}(j') $ are removed from $X$. This is continued until $X$ is empty. 

$C$ now represents data points that are relatively far apart. These points then choose a bundle of points $U_j$ that represent potential centers they would like to be assigned to. In particular, for each point $j \in C$, we look at all points $i \in X$ such that $0 < d(i,j) < 1.5 R_j$, where $R_j$ is half the distance from $j$ to its nearest point in $C$, and bundle these points into $U_j$ (if a point is in multiple bundles, choose one arbitrarily). Now these bundles of points have a sum of $s_i$ values of at least one half. For each $U_j$, we define the volume of $U_j$ to be $\mathrm{vol}(U_j) = \sum_{i \in U_j} s_i$. Lastly, these bundles are matched together greedily by distance. 

Now that the algorithm has grouped and matched together points so the volume is at least 1, the algorithm starts choosing points from these various groups through sampling. Here, opening $U_j$ refers to selecting a point $i \in U_j$ with probabilities proportional to the $s_i$ values and adding that point to $S$. 
\begin{enumerate}[\IEEEsetlabelwidth{3)}]
\item For each matched $U_j$ and $U_{j'}$, the algorithm opens $U_j$ with probability $1-\mathrm{vol}(U_{j'})$, opens $U_{j'}$ with probability $1- \mathrm{vol}(U_j)$, and opens both $U_{j}$ and $U_{j'}$ with probability $\mathrm{vol}(U_j) + \mathrm{vol}(U_{j'}) -1$. 
\item For any unmatched $U_j$, it is opened with probability $\mathrm{vol}(U_j)$.
\item For any $i$ not in any $U_j$ such that $s_i > 0$, add $i$ to $S$ with probability $s_i$. 
\end{enumerate}
This sampling process returns a set of centers of size $k$  in expectation and produces a set $S$ such that $\mathbb{E}[f(S)] \leq 3.25 \cdot f(O)$ for all $O \subseteq X$ such that $|O| = k$. For our purposes, we sample until we obtain a solution with exactly $k$ centers. 
<commit_msg>typos<commit_after>\documentclass[conference, 10pt, final]{IEEEtran}

\usepackage{cite}
\usepackage[cmex10]{amsmath}
\usepackage{amssymb}
\usepackage[pdftex]{graphicx}
\usepackage{algorithmic}
\usepackage{array}
\usepackage[caption=false,font=footnotesize]{subfig}
\usepackage{stfloats}
\usepackage{url}


\begin{document} 

\title{Improved Hierarchical Clustering Through Nesting}


\author{\IEEEauthorblockN{<NAME>, <NAME>, and <NAME>}
\IEEEauthorblockA{Operations Research and Information Engineering, Cornell University, Ithaca, NY.  \\
Email: \{ ajp336, cjw278, dal299 \} @corn<EMAIL>}}


\maketitle

\begin{abstract}
Hierarchical clustering organizes large data into smaller, more informative groups and offers a global view of the data structure by producing a hierarchy of clusters starting from individual data points and merging towards the cluster of all data. Typically, hierarchical clustering is done by greedily either merging or dividing clusters, called agglomerative and divisive clustering, respectively. Lin et al.~\cite{Lin} proposed a nesting framework that produces a hierarchal clustering from partitional clusterings for every possible number of clusters. This algorithm was introduced as an approximation algorithm for the hierarchical $k$-median problem and guarantees the resulting hierarchical layers are a limited distance from the original partitional clusters. However, this algorithm has not been previously placed in a clustering context. We test this algorithm using various methods for producing the partitional clusterings against the traditional hierarchical clustering methods and show both the approximation ratio and FScore measure improve significantly. We also propose a heuristic for this nesting framework which improves the run time and only suffers a slight decrease in approximation and FScore.
\end{abstract}


\section{Introduction and Notation}
As the data we have access to grows in variety and size and becomes cheaper to store, we need methods to efficiently learn from this data. 
Partitioning data (clustering) into similar groups (clusters) is one of the key tools to learn about a set of unlabeled data such as a set of customers who have visited a website or a set of handwritten documents.
Clustering in these situations helps us learn underlying structure and classifications of our data, and it allows large data sets to be presented in a simplified and compressed view \cite{Jain}.

There are two main branches of clustering: partitional (or flat) and hierarchical. In partitional clustering, we often specify a number $k$ and pairwise distances between the data points, representing a measure of dissimilarity such as Euclidean distance. The clustering algorithm then partitions the data into $k$ disjoint sets where the sets represent relatively similar data. In hierarchical clustering, on the other hand, we are interested in finding a nested tree of partitions. This offers the user a global view of the data and allows the user to choose a ``natural'' level of clustering \cite{ESL}. For example, a retailer could use a hierarchical clustering of customers to first find regions of customers in order to tailor the products offered in each region. The retailer could then also offer specialized discounts to subgroups in each region based on each cluster's shopping records. Examples of both types of clustering are given in Fig.~\ref{fig_part}-\ref{fig_hier}.

Clustering algorithms date back more than 50 years, and there have been many variations introduced since then \cite{Jain}. However, clustering algorithms are often formulated as NP-hard problems with non-convex objective functions, which are inherently hard to solve \cite{Jain}.
One common objective function for partitional clustering is to minimize the sum of distances from each point to the median point of its cluster, called $k$-median or $k$-medoids clustering. More specifically, suppose we are given a set $X$ and a distance (dissimilarity) metric $d(i,j)$ for all $i,j \in X$. The goal is to find a set of centers $S \subseteq X$ such that $|S| \leq k$ so as to minimize 
\[ f(S) = \sum_{ i\in X} \min_{j \in S} d(i,j). \]
The clustering is obtained from the set of centers $S$ by assigning each data point to its closest center.  An $\alpha$-approximation algorithm for this problem produces a set $S$ in polynomial time such that $f(S) \leq \alpha \cdot f(O)$ for all $O \subseteq X$ with $|O| = k$. 

\begin{figure*}[!]
\centerline{\subfloat[Partitional Clustering]{\includegraphics[width =1.75in]{partitional_clustering}
\label{fig_part}}
\hfil
\subfloat[Hierarchical Clustering]{\includegraphics[width=2.5in]{hier_clustering}
\label{fig_hier}}}
\caption{A partitional clustering for $k=3$, where the colors represent the different partitions, and a dendrogram representation of a hierarchical clustering of the data where edges represent merges of clusters and the clustering on the third level down (represented by the red line) corresponds to the partitional clustering on the left.}
\end{figure*}

This objective function can then be extended to the hierarchical $k$-median problem which represents an objective function relating to the hierarchical clustering problem. Suppose we represent a solution to the $k$-median problem as a pair $(S, a)$ where $S \subseteq X$ such that $|S| \leq k$ and $a$ is an assignment function that maps each $i \in X$ to a center in $S$. Then, for any $i \in S$, $a^{-1}(i)$ is the cluster centered at $i$. We define $f(S,a) := \sum_{i \in X} d(i, a(i))$. This is the same objective function as above except the assignment is given rather than implied.

Given solution pairs $(S_1, a_1)$ and $(S_2, a_2)$, we say that $(S_1, a_1)$ is \emph{nested} in $(S_2, a_2)$ (denoted $(S_1, a_1) \subseteq (S_2, a_2)$) if the following conditions hold:
\begin{enumerate}[\IEEEsetlabelwidth{3)}]
	\item $S_1 \subset S_2$, 
	\item $\forall  j \in X$, if $a_2(j) \in S_1$, then $a_1(j) = a_2(j)$, 
	\item $\forall j, k \in X$, if $a_2(j) = a_2(k)$, then $a_1(j) = a_1(k)$.
\end{enumerate}
A solution to the hierarchical $k$-median problem is a nested chain of solutions, and the hierarchical clustering it generates is found by looking at the clusters given by each $(S_i, a_i)$. An example is given in Fig.~\ref{fig_chain1}-\ref{fig_chain4}. However, we want this solution to stay close to the optimal $k$-median solutions when looking at any level of the hierarchy. An $\alpha$-approximation algorithm for this problem generates a nested chain of solutions in polynomial time such that for all $k = 1, \ldots, |X|$ there exists $(S_i, a_i)$ in the chain such that $|S_i| \leq k$ and $f(S_i, a_i) \leq \alpha \cdot f(O)$ for all $O \subseteq X$ with $|O| = k$.

\begin{figure*}[b]
\centering
\subfloat[$(S_1, a_1)$]{\includegraphics[width =1.75in]{chain1}
\label{fig_chain1}} \hspace{15mm}
\subfloat[$(S_2, a_2)$]{\includegraphics[width=1.75in]{chain2}
\label{fig_chain2}} \\
\subfloat[$(S_3, a_3)$]{\includegraphics[width=1.75in]{chain3}
\label{fig_chain3}}  \hspace{15mm}
\subfloat[$(S_4, a_4)$]{\includegraphics[width=1.75in]{chain4}
\label{fig_chain4}}
\caption{An example of a nested chain of solutions. $S_i$ is indicated by the colored nodes and $a_i$ by the edges. The clusters on that level are given by the connected components. Note that $S_{i-1}$ is a subset of $S_i$ and the assignments stay consistent. }
\label{fig_sim}
\end{figure*}

While some partitional clustering methods use the $k$-median objective function to inform their decisions, popular methods for hierarchical clustering do not and instead rely on greedy algorithms which have no guarantee on the solution produced with respect to the objective function but are typically fast and easy to implement. This paper looks at how using this objective function can vastly improve the hierarchical clusterings produced. 
In particular, we show how an approximation framework for this problem, which takes in partitional clusterings for all values of $k$, compares in practice to some of the most widely used clustering algorithms, both in the quality of the solution and in computation time. The algorithms we use to produce these partitional clusterings include both approximation algorithms for the $k$-median problem, which offer an overall approximation guarantee for the hierarchical $k$-median problem, and some heuristic algorithms that have been shown to perform well in practice. 

\section{Literature Review}

\subsection{Partitional Clustering}
Arguably the most common partitional clustering method is $k$-means, which chooses $k$ cluster centers (which are points in space, not necessarily points from the sample) to minimize the sum of squared distances from each point to its closest cluster center \cite{Jain}.  $k$-means is usually solved using Lloyd's algorithm (also referred to simply as the $k$-means algorithm), which iteratively greedily updates the centers of each cluster (assigning points to the closest center) until the solution becomes stable \cite{Jain}. Lloyd's algorithm is only guaranteed to reach a locally optimum solution and can be very sensitive to outliers and the initial points chosen \cite{Kanungo}.

The $k$-median objective function removes some worst-case behavior of $k$-means by restricting the cluster centers to be data points.  There are two common algorithms for this problem. The first is a local search algorithm called partition around medoids (PAM), which starts with $k$ randomly chosen data points as centers, and looks to swap one of the centers with another data point to improve the objective value. This repeats until no more improving swaps exist \cite{ESL}. The generalization of this algorithm which considers swaps of size $\leq p$ is a $(3 + \frac{2}{p})$-approximation algorithm \cite{Arya}.
Another common heuristic (called $k$-medoids throughout this paper), which is more computationally efficient than PAM \cite{Park}, is a variation on Lloyd's $k$-means algorithm.  Starting with randomly chosen centers, we assign each data point to its closest center, and then update the center within each cluster \cite{Park}. The algorithm repeats this until the solution becomes stable.

The $k$-median problem is also a fundamental problem in approximation algorithms where it is often presented in the context of opening service centers or factories to serve clients. However, these algorithms are often not implemented for clustering despite their improved guarantees on the solution. The best known approximation algorithm for this problem has an approximation guarantee of 2.592 and is produced using linear programming \cite{Wu}. A simpler LP rounding algorithm by Charikar and Li has a slightly worse guarantee of 3.25 but is easier to implement in practice \cite{Charikar}. 

\subsection{Hierarchical Clustering}
The most common method for hierarchical clustering is agglomerative clustering. Agglomerative clustering starts with each data point in its own cluster, and proceeds to recursively merge the two closest clusters to form the next clustering level in the hierarchy, repeating until the entire data set is a single cluster.  There are several choices on how to measure distance in order to choose which clusters to merge. Single linkage (also called nearest neighbor) defines the distance between two clusters to be the minimum pairwise distance between the data points in each cluster.  Complete linkage clustering instead defines the distance to be the largest pairwise distance. A method between these two extremes is group average linkage which looks at the average pairwise distance between clusters. This strikes a balance between the closeness and compactness of the two techniques above and gives more reasonable results in practice \cite{ESL}. Lastly, there is Ward's method which chooses the two clusters to merge so as to minimize the sum of each cluster's variance \cite{Ward}.

Approximation algorithms for hierarchical clustering have also been studied, where we aim to identify a hierarchy of clusters so as to minimize the ratio of the objective value on every level in the hierarchy of size $k$ to the optimal objective value of the partitional clustering problem with $k$ clusters. Note that it may not always be possible to find a solution with approximation ratio 1. The hierarchical $k$-median problem was originally introduced by Plaxton who gave a 238.88 approximation algorithm \cite{Plaxton}. This was improved by Lin et al.~\cite{Lin} who gave an incremental approximation algorithm framework. In particular, their approach relies on a black box $k$-median solver and then uses these solutions for $k=1, \ldots, n$ to greedily construct a hierarchical clustering that increases the objective value by at most a factor of $20.71$ (or $10.03$ when randomized) at each level. Using the local $p$-swap $(3+\frac{2}{p})$-approximation algorithm mentioned above as the black-box $k$-median solver performed well in experiments on hierarchical $k$-median test problems \cite{Nagarajan}. While the authors tested the variants of this algorithm for different $k$-median approximation algorithms, it has not been tested against traditional clustering algorithms or on any labeled data sets. 

Several algorithms have been built upon agglomerative clustering that we choose not to compare with our proposed algorithms. Divisive clustering, similar to agglomerative, starts with the data in one cluster and repeatedly splits a cluster on each level. This method has not been studied nearly as much as agglomerative clustering \cite{ESL} but has had relative success in classifying document datasets \cite{Zhao}.  Algorithms such as the CURE, ROCK, and Chameleon \cite{Guha_CURE, Guha_ROCK, Karypis} incorporate non-spherical cluster shapes and cluster density, but do not relate as much to the above approximation algorithm in their objective. We also choose not to compare our algorithms with Bayesian hierarchical methods since these methods require many more parameters and add additional prior information that our algorithms do not use. 

\section{Proposed Algorithms}

\subsection{Hierarchical Algorithm} 

Recall the hierarchical $k$-median problem: given a set $X$ and a distance metric $d(i,j)$ for all $i,j \in X$, we represent a solution to the $k$-median problem as a pair $(S, a)$ where $S \subseteq X$ such that $|S| \leq k$ and $a$ is an assignment function that maps each $i \in X$ to a center in $S$.
 A solution to the hierarchical $k$-median problem is a chain of nested solutions $(S_1, a_1) \subseteq (S_2, a_2) \subseteq \ldots \subseteq (S_l, a_l)$. However, we want this solution to stay close to the optimal $k$-median solutions when looking at any level of the hierarchy.  
Consider any $k$-median algorithm $\mathcal{A}$. After running $\mathcal{A}$ for all potential values of $k$, Lin et al.~\cite{Lin} introduce an augmentation step that nests these clusterings into the hierarchical framework and only increase costs by a bounded factor. 

First, set $\beta = 3+\sqrt{3}$, and $\beta_0 = 1$ (for the randomized algorithm we switch $\beta = 6.355$ and set $\beta_0 = \beta^X$, where $X$ is chosen from a uniform distribution on $[0,1]$). Let $V_i$ be the solutions generated from running $\mathcal{A}$ for $k=n-i+1$ for $i=1, \ldots, n$. The cost of these solutions is determined by assigning each point to the nearest cluster center and is scaled so that the minimum cost for $k<n$ is at least one. We then bucket these solutions by cost into buckets of the form $[0,0], (\beta_0, \beta_0 \beta], (\beta_0 \beta, \beta_0 \beta^2], \ldots $. 

From each non-empty bucket, pick the solution with the smallest number of cluster centers, and let these solutions be $\overline{V}_1, \overline{V}_2, \ldots \overline{V}_r$. We then inductively build our nested solutions: Set $S_1 = \overline{V}_1$ and $a_1$ to be the assignment that assigns each point to the closest point in $S_1$. Given $(S_{i-1}, a_{i-1})$ and $\overline{V}_i$, we construct $S_i$ by  adding the closest point $j \in S_{i-1}$ to $l$ for each $l \in \overline{V}_i$. Note that $|S_{i}| \leq |\overline{V}_i|$. 

\begin{figure*}[!t]
\centering
\subfloat[Nested solution $(S_{i-1}, a_{i-1})$ ]{\includegraphics[width =1.75in]{Si-1}
\label{fig_Si-1}} \hspace{15mm}
\subfloat[Next set $\bar{V}_i$]{\includegraphics[width=1.75in]{Vi}
\label{fig_Vi}} \\
\subfloat[Next nested centers $S_i$]{\includegraphics[width=1.75in]{Si}
\label{fig_Si}}  \hspace{15mm}
\subfloat[Next nested solution $(S_i, a_i)$]{\includegraphics[width=1.75in]{Siai}
\label{fig_Siai}}
\caption{An example of nesting $\bar{V}_i$ with $(S_{i-1}, a_{i-1})$. Dark  nodes represent chosen centers and edges represent assignments to that center's cluster. First, node 2 chooses node from 1 from $S_{i-1}$ as its closest center and node 5 chooses node 4. Then, each cluster in $(S_{i-1}, a_{i-1})$ is assigned to the cheapest (lowest sum of distances) center.}
\label{fig_sim}
\end{figure*}

We then need to define an assignment $a_i$. For each $j \in X$ such that $a_{i-1}(j) \in S_{i}$, we set $a_i(j) = a_{i-1}(j)$. Otherwise, we set $a_i(j)$ to be the facility $i \in S_i$ such that assigning all points in $a_{i-1}^{-1}(j)$ (i.e. all points in $j$'s cluster in $S_{i-1}$) would have minimum cost. Using this assignment, $(a_i, S_i)$ is nested in $(a_{i-1}, S_{i-1})$. An example of this nesting is given in Fig.~\ref{fig_Si-1}-\ref{fig_Siai}.

We continue this nesting until $|S_i| = 1$, reaching the top level of the hierarchy where all data points are in the same cluster. After forming this hierarchical clustering, for any $k \in 1, \ldots, n$ \cite{Lin} guarantees that the cost of $S_i$ is at most $20.71$ times the cost of $V_k$ where $|S_{i-1}| \geq |V_k| \geq |S_{i}|$ (in the randomized case the expected cost of $S_i$ is at most $10.03$ times the cost of $V_k$).


\subsection{Partitional Clustering Algorithms}
This nesting framework allows us to choose the partitional clustering algorithm $\mathcal{A}$. As mentioned above, PAM and $k$-medoids are two potential algorithms. PAM has the added benefit of being a 5-approximation algorithm for the partitional $k$-median problem which yields an overall approximation ratio of $20.71 \cdot 5 = 103.55$ for the hierarchical $k$-median problem ($10.03 \cdot 5 = 50.15$ randomized). We looked at increasing the size of the local swaps in PAM to improve the approximation ratio but found that the approximation ratio of constructed solutions for given values of $k$ did not significantly improve while the computation time increased greatly. 

We also consider three linear programming based algorithms, one of which has an approximation guarantee. We first present the linear programming formulation for the $k$-median problem. Suppose we have a set $X$. Let $s_i \in \{0,1\}$ be the variable representing whether or not we want to make $i$ a center for all $i \in X$. We also create a variable $a_{i,j} \in \{0, 1\}$ which represents whether or not we assign $i \in X$ to the cluster with center $j \in X$. Given this set and assignment, the $k$-median objective function is 
\begin{equation}
\sum_{i \in X} \sum_{j \in X} a_{i,j} d(i,j), 
\end{equation}
which we want to minimize. However, we are restricted by the fact that we can choose at most $k$ centers, every point must be assigned to a single cluster center, and that we can only assign a point to a center that is actually in $S$. This is represented by the following three equations, respectively:
\begin{eqnarray}
\sum_{i \in X} s_i   \leq k & ,\\
\sum_{j \in X} a_{i,j}  = 1 & \forall i \in X, \\
a_{i,j}   \leq s_j & \forall i, j \in X.
\end{eqnarray}
Therefore, solving the following integer program would directly find the optimal set of centers and assignment for the $k$-median problem.
\begin{eqnarray}
\text{Min } & \sum_{i \in X} \sum_{j \in X} a_{i,j} d(i,j) \\
\text{s.t. } & \sum_{i \in X} s_i   \leq k  \\
& \sum_{j \in X} a_{i,j}  = 1 \quad \forall i \in X \\
& a_{i,j}   \leq s_j \quad \forall i, j \in X \\
& s_i \in \{0,1\} \quad \forall i \in X \\
& a_{i,j} \in \{0,1\} \quad \forall i, j \in X
\end{eqnarray}
Relaxing the constraints so that $s_i \in [0,1]$ for all $i \in X$ and $a_{i,j} \in [0,1]$ for all $i,j \in X$ yields the linear programming (LP) relaxation for the $k$-median problem, which we can solve in polynomial time.

Interestingly, the solution to the linear programming is often integer valued \cite{Nagarajan}. Based on this fact, we propose two LP-based clustering algorithms: one that solves the LP relaxation and constructs $S$ by taking the points in $X$ with the $k$ highest $s_i$ values and one that solves the LP relaxation and constructs $S$ by choosing $i$ to be in $S$ with probability $s_i$, sampling for each $i \in X$ until $k$ centers are chosen. We refer to these algorithms as LP-greedy and LP-rand, respectively.
Charikar and Li \cite{Charikar} introduce a 3.25-approximation algorithm that uses this LP solution to smartly construct a set of centers. We refer to this algorithm as LP-Charikar. We choose to implement this LP algorithm rather than the LP-based approximation algorithm with best approximation ratio because the algorithm is simpler to implement (and we expect it to have a lower running time). 

The algorithm works by first choosing a subset of data points $C \subseteq X$ to represent the demand. For each data point $i$, let $d_{av}(i) = \sum_{j \in X : a_{i,j} > 0} s_j \cdot d(i,j)$. We start with $C = \emptyset$. At each step, we choose a point $j \in X$ to add to $C$, and all points $j'$ such that $d(j,j') \leq 4 d_{av}(j') $ are removed from $X$. This is continued until $X$ is empty. 

$C$ now represents data points that are relatively far apart. These points then choose a bundle of points $U_j$ that represent potential centers they would like to be assigned to. In particular, for each point $j \in C$, we look at all points $i \in X$ such that $0 < d(i,j) < 1.5 R_j$, where $R_j$ is half the distance from $j$ to its nearest point in $C$, and bundle these points into $U_j$ (if a point is in multiple bundles, choose one arbitrarily). Now these bundles of points have a sum of $s_i$ values of at least one half. For each $U_j$, we define the volume of $U_j$ to be $\mathrm{vol}(U_j) = \sum_{i \in U_j} s_i$. Lastly, these bundles are matched together greedily by distance. 

Now that the algorithm has grouped and matched together points so that the volume is at least 1, the algorithm starts choosing points from these various groups through sampling. Here, opening $U_j$ refers to selecting a point $i \in U_j$ with probabilities proportional to the $s_i$ values and adding that point to $S$. 
\begin{enumerate}[\IEEEsetlabelwidth{3)}]
\item For each matched $U_j$ and $U_{j'}$, the algorithm opens $U_j$ with probability $1-\mathrm{vol}(U_{j'})$, opens $U_{j'}$ with probability $1- \mathrm{vol}(U_j)$, and opens both $U_{j}$ and $U_{j'}$ with probability $\mathrm{vol}(U_j) + \mathrm{vol}(U_{j'}) -1$. 
\item For any unmatched $U_j$, it is opened with probability $\mathrm{vol}(U_j)$.
\item For any $i$ not in any $U_j$ such that $s_i > 0$, add $i$ to $S$ with probability $s_i$. 
\end{enumerate}
This sampling process returns a set of centers of size $k$  in expectation and produces a set $S$ such that $\mathbb{E}[f(S)] \leq 3.25 \cdot f(O)$ for all $O \subseteq X$ such that $|O| = k$. For our purposes, we sample until we obtain a solution with exactly $k$ centers. 

{'fraction_non_alphanumeric': 0.08710021321961621, 'fraction_numerical': 0.010191897654584222, 'mean_word_length': 4.63985088985089, 'pattern_counts': {'":': 0, '<': 17, '<?xml version=': 0, '>': 15, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'dir': 'git-commits-cleaned', 'id': '1957693', 'n_tokens_mistral': 13348, 'n_tokens_neox': 12547, 'n_words': 6859}