# Unsupervised Manifold Linearizing and Clustering

Tianjiao Ding1 Shengbang Tong2 Kwan Ho Ryan Chan1 Xili Dai3 Yi Ma2 Benjamin D.
Haeffele1

###### Abstract

Clustering data lying close to a union of low-dimensional manifolds, with each
manifold as a cluster, is a fundamental problem in machine learning. When the
manifolds are assumed to be linear subspaces, many methods succeed using low-
rank and sparse priors, which have been studied extensively over the past two
decades. Unfortunately, most real-world datasets can not be well approximated
by linear subspaces. On the other hand, several works have proposed to
identify the manifolds by learning a feature map such that the data
transformed by the map lie in a union of linear subspaces, even though the
original data are from non-linear manifolds. However, most works either assume
knowledge of the membership of samples to clusters, or are shown to learn
trivial representations. In this paper, we propose to simultaneously perform
clustering and learn a union-of-subspace representation via Maximal Coding
Rate Reduction. Experiments on synthetic and realistic datasets show that the
proposed method achieves clustering accuracy comparable with state-of-the-art
alternatives, while being more scalable and learning geometrically meaningful
representations.

11footnotetext: Mathematical Institute for Data Science, Johns Hopkins
University, USA 2Department of Electrical Engineering and Computer Sciences,
University of California, Berkeley, USA 3The Hong Kong University of Science
and Technology (Guangzhou), PRC

(a)

(b)

(c)

(d)

Figure 1: (a) Input data ${\bm{X}}$ of two manifolds each containing $100$
points. (b) Features ${\bm{Z}}_{\bm{\theta}}$ at random initialization. (c)
${\bm{Z}}_{\bm{\theta}}$ after self-supervised initialization. (d)
${\bm{Z}}_{\bm{\theta}}$ after MLC (4) training. Details are in the Appendix.

##  1 Introduction

###  1.1 Motivation and Contributions

Clustering is a fundamental problem in machine learning, allowing one to group
data into clusters based on assumptions about the geometry of clusters. For
example, when data are concentrated around distinct centroids, classical
k-means clustering [26, 14, 17, 29] and its variants [5, 2, 3] are able to
find the cluster centroids and assign membership to each data point. More
generally444This includes k-means-based methods, since a centroid is a
$0$-dimensional affine subspace., subspace clustering methods [11, 13, 27, 25,
16, 50] are designed to cluster data that lie close to a union of different
low-dimensional linear (or affine) subspaces, where each subspace defines a
cluster. Overall, those methods often enjoy theoretical guarantees of correct
clustering [40, 27, 13, 41, 47, 48, 51, 50, 44, 22, 49, 37] and find
applications in various problems such as image clustering, face recognition,
motion segmentation, and recently in popular Transformer architectures in deep
learning [38].

Despite the wide range of applications and theoretical guarantees, subspace
clustering methods rely on a crucial assumption that each cluster can be well
approximated by a linear/affine subspace, which is often not valid for many
real-world datasets. For instance, even in a dataset as simple as MNIST hand-
written digits, images of a single digit do not lie close to a low-dimensional
linear subspace, thus directly applying subspace clustering will fail.
Instead, it is more natural to assume the clusters are from non-linear low-
dimensional manifolds (one manifold per cluster), and attempt to learn or
design a non-linear embedding of the data so that the transformed data lies
close to distinct linear subspaces, with points from one manifold mapped to
the same subspace. For example, [24] shows that a subspace clustering method
can achieve $99\%$ clustering accuracy on MNIST images after embedding the
data with the scattering transform [6].

Beyond the above example, numerous other subspace clustering methods have
explored hand-designing an appropriate feature embedding (or kernel) such as
polynomial or exponential mappings [34]. However, these embeddings assume
specific families of manifolds, thus they need to be hand-crafted for various
tasks and datasets using domain knowledge, which makes their application
challenging for complicated data such as natural images. On the other hand,
[12] proposes to cluster data based on treating a local neighborhood of the
manifold approximately as a linear subspace. However, for this to succeed
sufficient sampling density is required, which implies a prohibitive number of
samples when the manifolds are of moderate dimension or are highly curved.
Further, for a new sample unseen at training time one needs to run the
algorithm with all samples to embed it or assign a membership to it, which is
expensive computationally. More recently, numerous works propose to learn an
appropriate linear embedding of the data via deep networks and then perform
subspace clustering in the feature space [35, 18, 1, 54, 20]. Unfortunately,
it has been shown that many of these formulations are provably ill-posed and
learn trivial representations555In this paper, we use ‘representation’ and
‘feature’ interchangeably to mean the image of data under a (learned)
transformation., with much of the claimed benefit coming from ad-hoc post-
processing rather than the method itself [15]. This leads to the following
question:

######  Question 1.

For data approximately supported on an underlying union of manifolds, can we
learn a transformation of the data, so that the transformed data lie in
distinct linear subspaces to be easily clustered?

Meanwhile, learning a representation from multi-modal data has been a topic of
its own interest in machine learning. An ideal property of the learned
representation often pursued is between-cluster discrimination, i.e., features
from different clusters should be well separated. Further, an important yet
often ignored property of the learned representation is that it maintains
within-cluster diversity. This is desirable as it allows distances of samples
within a cluster to be preserved under the learned transformation, which could
facilitate downstream tasks such as denoising, hierarchical clustering and
semantic interpretation. In the supervised setting, training with the cross-
entropy (CE) classification objective fails to achieve the second property, as
it has been shown empirically [32] and theoretically [43, 56] that the
representation learned by CE has the property that features from one cluster
tend to collapse to a single point. On the other hand, recent work has
proposed the principle of Maximal Coding Rate Reduction (MCR2) [53] as one of
the few methods that are able to achieve the two ideal properties by learning
a representation where features from each cluster are expected to lie close to
a low-dimensional subspace (within-cluster diverse), and the subspaces from
different clusters are orthogonal to each other (between-cluster
discriminative). However, for MCR2 to learn such orthogonal subspaces each
corresponding to one cluster, one needs the annotation of which sample belong
to which cluster. Such annotation might be expensive or impossible to acquire
for large-scale datasets. This motivates another question of interest:

######  Question 2.

Can we learn a union-of-orthogonal-subspace representation of data coming from
an underlying union of manifolds without access to the labels?

This paper gives positive answers to the two interrelated questions by making
the following contributions.

  1. 1.

We propose to simultaneously cluster the data and learn a union-of-orthogonal-
subspace representation via MCR2, when data is assumed to lie close to a union
of manifolds. This is achieved by formulation (4), which optimizes over both
the representation and a doubly stochastic membership formulation inspired by
the state-of-the-art subspace clustering result [24].

  2. 2.

Since the membership has as many entries as the square of the batch size of
the input data, we give a parameterization of the membership (Figure 2).
Further, as problem (4) is highly non-convex, we give a meta-algorithm
(Algorithm 1) on how to initialize the variables and to optimize it.

  3. 3.

We conduct experiments on simulation and CIFAR10 to demonstrate some desirable
properties of the proposed method. We further experiment on datasets with
larger number of clusters and imbalanced clusters such as CIFAR100-20,
CIFAR100-100, and TinyImageNet200, and show that the proposed method achieves
state-of-the-art performance.

###  1.2 Additional Related Work

Beyond the above, we make connections to a few important works that are
related to this paper.

Deep Clustering and Representation Learning. Recently, there is an interesting
line of research in representation learning and clustering that takes
advantage of pseudo-labelling and semi/self-supervised learning [7, 45, 33,
30]. Specifically, one first identifies a subset of samples (often termed
reliable samples) based on geometric or statistical criteria in the learned
representation and cluster prediction, and then uses the predicted labels for
those reliable samples as if they are ground-truth labels to refine the
representation and cluster prediction of other samples. Despite the promising
clustering performance, the representation learned by these methods are not
constrained to be both between-cluster discriminative and within-cluster
diverse. In contrast, the proposed method learns a representation with these
two ideal properties (see Figure 4) and also achieves state-of-the-art
clustering performance (see Tables 2 and 4).

Neural Manifold Clustering and Embedding (NMCE). A recent preprint [23] also
proposes a solution to the same problem we study, i.e., clustering the data
and learning an union-of-orthogonal-subspace representation. In particular,
[23] proposes to model the point-to-cluster membership and optimize MCR2 [53]
over both the representation and the membership. In this paper, we adopt a
similar formulation, but we propose to model the point-to-point affinity using
a doubly stochastic matrix, inspired by the state-of-the-art subspace
clustering methods (§2.2). Aside from having different conceptual formulations
and algorithms, our formulation is much more stable with respect to
initialization and is naturally suitable for hierarchical clustering. We
detail these distinctions in §2.2. Experiments (Table 2) further show that the
proposed method (MLC) achieves higher accuracy than [23] (NMCE) on large scale
realistic datasets.

##  2 Problem Formulation

We start by defining the problem that we study. Suppose
${\bm{X}}=[{\bm{x}}_{1},\ldots,{\bm{x}}_{n}]\in{\mathbb{R}}^{D\times n}$ is a
dataset of $n$ samples drawn from a union of $k$ underlying manifolds
$\bigcup_{j=1}^{k}{\mathcal{M}}_{j}$ and ${\bm{y}}\in{\mathbb{R}}^{n}$ their
memberships to the manifolds, i.e., ${\bm{x}}_{i}\in{\mathcal{M}}_{y(i)}$.

######  Problem 1 (Unsupervised Manifold Linearizing and Clustering).

Given the dataset ${\bm{X}}$, can we simultaneously 1) cluster the samples,
i.e., estimate ${\bm{y}}$, and 2) learn a linear representation for manifolds,
i.e., find a transformation $f:{\mathbb{R}}^{D}\to{\mathbb{R}}^{d}$, such that
the image of each manifold $f({\mathcal{M}}_{i})$ is a low-dimensional linear
subspace of ${\mathbb{R}}^{d}$, and the subspaces satisfy desired properties
(§1), i.e., they are between-cluster discriminative and within-cluster
diverse?

Here we base our approach on the principle of Maximal Coding Rate Reduction
(MCR2) which is designed to learn ideal representations in the supervised
case, i.e., when the membership ${\bm{y}}$ is given (§2.1). Then we discuss
the challenges of simultaneously clustering and learning representation
(Problem 1), and propose our MCR2 clustering objective to solve Problem 1 with
those challenges in mind (§2.2). We further give an algorithm to optimize the
proposed objective (§2.3).

###  2.1 Supervised Manifold Linearizing via MCR2

In the case when the labels ${\bm{y}}$ are given as supervision, MCR2 [53]
aims to address part 2) of Problem 1. Let
$f_{\bm{\theta}}:{\mathbb{R}}^{D}\to{\mathbb{S}}^{d-1}$ be a featurizer
parameterized by a neural network, which in turn gives an embedding
${\bm{Z}}_{\bm{\theta}}:=[{\bm{z}}_{1},\dots,{\bm{z}}_{n}]\in{\mathbb{R}}^{d\times
n}$ of data with
${\bm{z}}_{i}:=f_{\bm{\theta}}({\bm{x}}_{i})\in{\mathbb{S}}^{d-1}$. MCR2 aims
to learn an ideal representation by optimizing

$\displaystyle\max_{\bm{\theta}}\quad
R({\bm{Z}}_{\bm{\theta}};\,\epsilon)-R_{c}({\bm{Z}}_{\bm{\theta}},\bm{\Pi};\,\epsilon)\quad\text{s.t.}\quad{\bm{Z}}_{\bm{\theta}}\in{\mathcal{S}}$
(1) $\displaystyle\text{where
}R({\bm{Z}}_{\bm{\theta}};\,\epsilon):=\log\det\left({\bm{I}}+\frac{d}{n\epsilon^{2}}{\bm{Z}}_{\bm{\theta}}{\bm{Z}}_{\bm{\theta}}^{\top}\right),$
$\displaystyle\text{and }R_{c}({\bm{Z}}_{\bm{\theta}},\bm{\Pi};\,\epsilon):=$
$\displaystyle\sum_{j=1}^{k}\frac{\langle\bm{\Pi}_{j},{\bm{1}}\rangle}{n}\log\det\left({\bm{I}}+\frac{d}{\langle\bm{\Pi}_{j},{\bm{1}}\rangle\epsilon^{2}}{\bm{Z}}_{\bm{\theta}}\operatorname{Diag}(\bm{\Pi}_{j}){\bm{Z}}_{\bm{\theta}}^{\top}\right).$

Here ${\mathcal{S}}$ is the set of matrices whose columns all have unit
$\ell_{2}$ norm666This can be easily achieved by having the last layer of the
neural network $f_{\bm{\theta}}$ be a normalization layer.,
$\bm{\Pi}\in{\mathbb{R}}^{n\times k}$ is a given membership matrix such that
$\Pi_{ij}=1$ if $j=y(i)$ and $\Pi_{ij}=0$ otherwise, $\epsilon>0$ is a
prescribed precision parameter, $\bm{\Pi}_{j}\in{\mathbb{R}}^{n}$ denotes the
$j^{\text{th}}$ column of $\bm{\Pi}$, ${\bm{1}}$ is a vector of all ones so
that $\langle\bm{\Pi}_{j},{\bm{1}}\rangle$ is the number of points in cluster
$j$, and finally for ${\bm{v}}\in{\mathbb{R}}^{n}$,
$\operatorname{Diag}({\bm{v}})$ denotes a diagonal matrix with the entries of
${\bm{v}}$ along the diagonal.

Intuitively777More formally, terms of the form
$\log\det\left(I+\frac{d}{n\epsilon^{2}}{\bm{W}}{\bm{W}}^{\top}\right)$
estimate the average number of bits needed to code $n$ i.i.d. samples
${\bm{W}}\in{\mathbb{R}}^{d\times n}$ from a zero-mean $d$-dimensional
Gaussian up to a distortion $\epsilon$ [28], hence the name coding rate. , the
$R({\bm{Z}}_{\bm{\theta}};\,\epsilon)$ term of (1) measures the volume of
${\bm{Z}}_{\bm{\theta}}$, and maximizing it would diversify features from all
samples, which we will refer to as the _expansion term_ Likewise, the
$R_{c}({\bm{Z}}_{\bm{\theta}},\bm{\Pi};\,\epsilon)$ term measures the sum of
volumes of each cluster of ${\bm{Z}}_{\bm{\theta}}$ and is referred to as the
_compression term_ , since minimizing it would push features within each
cluster to stay close. It has been shown that given $\bm{\Pi}$, the features
obtained by maximizing the rate reduction
$R({\bm{Z}}_{\bm{\theta}};\,\epsilon)-R_{c}({\bm{Z}}_{\bm{\theta}},\bm{\Pi};\,\epsilon)$
has the property that the features of each cluster spread uniformly within a
subspace (within-cluster diverse), and the subspaces from different clusters
are orthogonal (between-cluster discriminative), under relatively mild
assumptions [53].

###  2.2 Unsupervised Manifold Linearizing and Clustering via MCR2

While the MCR2 formulation is designed to learn ideal representations (§1)
when the membership ${\bm{y}}$ (or equivalently $\bm{\Pi}$) is given, here we
are interested in the unsupervised setting where one does not have access to
membership annotations. Thus, we propose to simultaneously perform both parts
1) and 2) of Problem 1 by also optimizing over the membership $\bm{\Pi}$ of
data. This naturally leads to

$\displaystyle\max_{\bm{\theta},\bm{\Pi}\in\Omega_{\circ}}R({\bm{Z}}_{\bm{\theta}};\,\epsilon)-R_{c}({\bm{Z}}_{\bm{\theta}},\bm{\Pi};\,\epsilon)\quad\text{s.t.}\quad{\bm{Z}}_{\bm{\theta}}\in{\mathcal{S}},$
(2)

where $\Omega_{\circ}:=\\{\bm{\Pi}\in{\mathbb{R}}^{n\times k}:\forall
i\in[n],\,\,\exists\hat{y}(i)\quad\mathrm{s.t.}\,\,\Pi_{i\hat{y}(i)}=1\,\,\mathrm{and}\,\,\Pi_{ij}=0\,\,\mathrm{for}\,\,j\neq\hat{y}(i)\\}$
is the set of all ‘hard’ assignments, i.e., each row of $\bm{\Pi}$ is a one-
hot vector. However, this optimization is in general combinatorial: its
complexity grows exponentially in $n$ and $k$, and it does not allow smooth
and gradual changes of $\bm{\Pi}$. Further, a second challenge is the chicken-
and-egg nature of this problem: If one already has an ideal representation
${\bm{Z}}$, then existing subspace clustering methods can be applied on
${\bm{Z}}$ to estimate the membership. Likewise, if one is given the
membership $\bm{\Pi}$ of clusters, then solving (1) would lead to an ideal
representation. However, the ${\bm{Z}}_{\bm{\theta}}$ and $\bm{\Pi}$ at the
beginning of optimization is typically far from ideal.

Doubly Stochastic Subspace Clustering. To address the combinatorial of
estimating the cluster memberships, we draw inspiration from the closely
related problem of subspace clustering, where the goal is to cluster $n$
samples assumed to lie close to a union of $k$ low-dimensional subspaces (§1).
In this case, one typically does not directly learn an $n\times k$ matrix
denoting memberships of $n$ points into $k$ subspaces. Instead, one first
learns an affinity matrix $\bm{\Pi}\in{\mathbb{R}}^{n\times n}$ signaling the
similarity between pairs of points, and then applies spectral clustering on
the learned $\bm{\Pi}$ to obtain a final clustering [11, 13, 27, 25, 16, 50].
In particular, requiring doubly-stochastic constraints on the affinity
$\bm{\Pi}$ is shown theoretically to suppress false inter-cluster connections
for clustering problems [9] along with state-of-the-art empirical performance
for subspace clustering problems [24].

Inspired by the above, we propose a constraint set $\Omega$ for matrix
$\bm{\Pi}$ to be the set of $n\times n$ doubly stochastic matrices,

$\Omega=\\{\bm{\Pi}\in{\mathbb{R}}^{n\times n}:\bm{\Pi}\geq
0,\quad\bm{\Pi}{\bm{1}}=\bm{\Pi}^{\top}{\bm{1}}={\bm{1}}\\}.$ (3)

However, this constraint alone is insufficient for strong clustering
performance: Consider the optimization of (2) with respect to
$\bm{\Pi}\in\Omega$ only, and note that the objective is strongly convex with
respect to $\bm{\Pi}$. Since we maximize a convex function with respect to
convex constraints $\Omega$, an optimal $\bm{\Pi}$ would lie at an extreme
point of $\Omega$, which for doubly stochastic matrices is a permutation
matrix. This is not ideal for clustering, as it implies that every point is
assigned to its own distinct cluster, and there is no incentive to merge
points into larger clusters. To resolve this issue, we follow the approach in
[24] and add $\ell_{2}$ regularization888Other choices of regularization are
also possible: Essentially any function which achieves its minimum over
$\Omega$ at the uniform matrix could potentially be used, e.g., the negative
entropy function $\sum_{ij}\Pi_{ij}\log(\Pi_{ij})$.
$\frac{\gamma}{2}\|\bm{\Pi}\|_{F}^{2}$ to $\bm{\Pi}$ which biases $\bm{\Pi}$
toward the uniform matrix $\tfrac{1}{n}\bm{1}\bm{1}^{\top}$, so by tuning
$\gamma$ we can also tune the sparsity level of $\bm{\Pi}$. This results in
our final proposed formulation, dubbed Manifold Linearizing and Clustering
(MLC):

$\displaystyle\max_{\bm{\theta}}\quad
R({\bm{Z}}_{\bm{\theta}};\,\epsilon)-R_{c}({\bm{Z}}_{\bm{\theta}},\bm{\Pi}_{\bm{\theta}};\,\epsilon)-\frac{\gamma}{2}\|\bm{\Pi}_{\bm{\theta}}\|_{F}^{2}$
(4)
$\displaystyle\quad\text{s.t.}\quad{\bm{Z}}_{\bm{\theta}}\in{\mathcal{S}},\,\,\bm{\Pi}_{\bm{\theta}}\in\Omega,$
$\displaystyle\text{where
}R({\bm{Z}}_{\bm{\theta}};\,\epsilon)=\log\det\left({\bm{I}}+\frac{d}{n\epsilon^{2}}{\bm{Z}}_{\bm{\theta}}{\bm{Z}}_{\bm{\theta}}^{\top}\right),\text{
and}$
$\displaystyle\medmath{R_{c}({\bm{Z}}_{\bm{\theta}},\bm{\Pi}_{\bm{\theta}};\,\epsilon)=\frac{1}{n}\sum_{j=1}^{n}\log\det\left({\bm{I}}+\frac{d}{\epsilon^{2}}{\bm{Z}}_{\bm{\theta}}\operatorname{Diag}((\bm{\Pi}_{\bm{\theta}})_{j}){\bm{Z}}_{\bm{\theta}}^{\top}\right)}.$

Note that here $\bm{\Pi}_{\bm{\theta}}=\bm{\Pi}_{\bm{\theta}}({\bm{X}})$ is
now also parameterized by a neural network. While this is constrained
optimization which may appear difficult to handle, we explain in §2.3 how we
parameterize ${\bm{Z}}_{\bm{\theta}}$ and $\bm{\Pi}_{\bm{\theta}}$ via neural
networks so that the constraints are satisfied by construction. Below, we note
a few advantages of the proposed formulation.

Figure 2: Overall architecture for optimizing the proposed manifold
linearizing and clustering (MLC) objective (4). Given $n$ input samples
${\bm{X}}$ each lying in ${\mathbb{R}}^{D}$, their $d$-dimensional
representation is given by ${\bm{Z}}_{\bm{\theta}}({\bm{X}})$, where
$\bm{\theta}$ denotes network parameters. Further, their doubly stochastic
membership matrix $\bm{\Pi}_{\bm{\theta}}({\bm{X}})$ is given by taking an
inner product kernel of the output of the cluster head
${\bm{C}}_{\bm{\theta}}({\bm{X}})$ followed by a doubly stochastic projection.

Parameterizing $\bm{\Pi}$ via a Neural Network versus Free Variables. An
alternative way to parameterize the membership would be to directly take
$\bm{\Pi}$ as decision variables in $\Omega$, as opposed to outputs of a
neural network. However, this leads to maintaining $O(n^{2})$ variables which
is prohibitive for large datasets (e.g., $n=10^{6}$ for ImageNet). In
contrast, this is not the case if one parameterizes $\bm{\Pi}$ as a neural
network, since one can do stochastic gradient descent such that for each batch
both the memory and computational complexity is at most square of the batch
size (Figure 2).

Comparison with NMCE. As mentioned in §1.2, NMCE [23] approaches Problem 1
also by optimizing MCR2 over both the representation and membership. However,
in NMCE the membership is parameterized by an $n\times k$ matrix
$\bm{\Pi}_{n\times k}$ that models the point-cluster membership, which is
different from our doubly stochastic point-point membership matrix
$\bm{\Pi}_{\bm{\theta}}$ inspired from the state-of-the-art subspace
clustering. Note further that for NMCE the initialization of
$\bm{\Pi}_{n\times k}$ is arbitrary and has nothing to do with the structures
in the initialized representation $\bm{\Pi}_{\bm{\theta}}$, and a bad
initialization of $\bm{\Pi}_{n\times k}$ could lead to the features from
different true clusters being compressed. On the other hand, the proposed
doubly stochastic membership $\bm{\Pi}_{\bm{\theta}}$ can be initialized
deterministically using structures from self-supervised initialized features
${\bm{Z}}_{\bm{\theta}}$ (§2.3). Interestingly, optimizing (4) allows an
interpretation of linearizing each point with its neighbors. Empirically as
seen in (Table 2), the proposed MLC yields a higher clustering accuracy than
NMCE [23].

###  2.3 Algorithms

As mentioned, in the MLC objective (4), we parameterize both the
representation ${\bm{Z}}_{\bm{\theta}}$ and doubly stochastic membership
$\bm{\Pi}_{\bm{\theta}}$ via a neural network. Below we elaborate on how this
is done. We summarize the network architecture in Figure 2, and the meta
algorithm in Algorithm 1.

Parameterizing ${\bm{Z}}_{\bm{\theta}}$. We follow [53] and take some existing
network architecture as the backbone. We append a few affine layers with non-
linearity as the representation head to further transform the output in
${\mathbb{R}}^{d}$, followed by a projection layer to respect the unit sphere
${\mathbb{S}}^{d-1}$ constraint.

Parameterizing $\bm{\Pi}_{\bm{\theta}}$. In subspace clustering, the
membership $\bm{\Pi}$ given data ${\bm{X}}$ often takes the form of
$g({\bm{X}})^{\top}g({\bm{X}})$ for some (linear) transformation $g$, such as
in the inner product kernel [16, 9] where $g={\bm{I}}$ or the least square
regression [27] where $g=({\bm{I}}+\lambda{\bm{X}}^{\top}{\bm{X}})^{-1/2}$.
This motivates us to parameterize $g_{\bm{\theta}}$ by a neural network, and
take
${\bm{C}}_{\bm{\theta}}^{\top}{\bm{C}}_{\bm{\theta}}\in{\mathbb{R}}^{n\times
n}$ as the membership where ${\bm{C}}_{\bm{\theta}}$ is shorthand for
$g_{\bm{\theta}}({\bm{X}})$. Nevertheless, such an $n\times n$ matrix is in
general not doubly stochastic, i.e.,
${\bm{C}}_{\bm{\theta}}^{\top}{\bm{C}}_{\bm{\theta}}\notin\Omega$. To obtain a
doubly stochastic membership, we further apply a Sinkhorn projection layer
$P_{\Omega,\eta}(\cdot)$ [39, 10], which gives our final parameterization of
the membership as
$\bm{\Pi}_{\bm{\theta}}=P_{\Omega,\eta}({\bm{C}}_{\bm{\theta}}^{\top}{\bm{C}}_{\bm{\theta}})\in\Omega$.

Initializing ${\bm{Z}}_{\bm{\theta}}$: Self-supervised Representation Learning
via MCR2. Since the proposed MCR2 clustering objective (4) is non-convex, it
is important to properly initialize both ${\bm{Z}}$ and $\bm{\Pi}$ to converge
to good (local) minimum. On the other hand, randomly initialized features are
typically far from being ideal, since they may not satisfy the idealized
properties (§1), and further may not respect the invariance to augmentation,
i.e., the augmented samples should have their representation close to each
other. Thus, we adopt the self-supervised training strategy [23]

$\displaystyle\max_{\bm{\theta}}\quad
R\Big{(}\frac{{\bm{Z}}_{\bm{\theta}}+{\bm{Z}}_{\bm{\theta}}^{\prime}}{2};\,\epsilon\Big{)}+\lambda\sum_{i=1}^{n}|{\bm{z}}_{i}^{\top}{\bm{z}}^{\prime}_{i}|,$
$\displaystyle\quad\text{s.t.}\quad{\bm{z}}_{i}^{\prime},{\bm{z}}_{i}\in{\mathbb{S}}^{d-1},\quad\forall
i\in[n],$ (5)

where for every $i$, ${\bm{z}}_{i}$ and ${\bm{z}}^{\prime}_{i}$ are features
of different augmentations of the $i$-th sample. This essentially requires
that features from different augmentations of the same sample should be as
close as possible, whereas features from different samples should be as
uncorrelated as possible.

Initializing $\bm{\Pi}_{\bm{\theta}}$. An ideal initialization of
$\bm{\Pi}_{\bm{\theta}}$ would be such that if $(\bm{\Pi}_{\bm{\theta}})_{ij}$
has a high value then points $i,j$ are likely to be from the same true cluster
and vice versa. On the other hand, after the self-supervised feature
initialization mentioned above, ${\bm{Z}}_{\bm{\theta}}$ already have some
structures which we can utilize. Thus, we propose to initialize
$\bm{\Pi}_{\bm{\theta}}$ with
$P_{\Omega,\eta}({\bm{Z}}_{\bm{\theta}}^{\top}{\bm{Z}}_{\bm{\theta}})$, which
is easily implemented by copying the parameters from ${\bm{Z}}_{\bm{\theta}}$
to ${\bm{C}}_{\bm{\theta}}$ once after the self-supervised initialization of
the former, i.e., from the feature head to the cluster head in Figure 2.

Data Augmentation. Beyond initializing ${\bm{Z}}_{\bm{\theta}}$, it is often
desirable to incorporate augmentation in optimizing (4). Specifically, from
$\\{{\bm{X}}^{(a)}\in{\mathbb{R}}^{D\times n}\\}_{a=1}^{A}$ the dataset
${\bm{X}}$ under $A$ different augmentations, one computes
$({\bm{Z}}_{\bm{\theta}}^{(a)}\in{\mathbb{R}}^{d\times
n},\bm{\Pi}_{\bm{\theta}}^{(a)}\in{\mathbb{R}}^{d\times n})$ for each
augmentation $a$, and use in (4)

${\bm{Z}}_{\bm{\theta}}=P_{{\mathbb{S}}^{d-1}}\left(\frac{1}{A}\sum_{a=1}^{A}{\bm{Z}}_{\bm{\theta}}^{(a)}\right),\quad\bm{\Pi}_{\bm{\theta}}=\frac{1}{A}\sum_{a=1}^{A}\bm{\Pi}_{\bm{\theta}}^{(a)}\in\Omega.$
(6)

Note that one can benefit from parallelization by putting
${\bm{X}}^{(a)},{\bm{Z}}_{\bm{\theta}}^{(a)},\bm{\Pi}_{\bm{\theta}}^{(a)}$ for
each augmentation $a$ on one computing device, since
$\bm{\Pi}_{\bm{\theta}}^{(a)}$ only depends on ${\bm{X}}^{(a)}$ but not from
other augmentations.

Algorithm 1 MLC: Unsupervised Manifold Linearizing and Clustering

1:${\bm{X}}\in{\mathbb{R}}^{D\times n}$, $\epsilon,\gamma,\eta,\lambda>0$,
$d,k,n_{b},T,A\in{\mathbb{Z}}_{\geq 0}$

2:initialize ${\bm{Z}}_{\theta}$ by self-supervised representation learning
via MCR2 $\triangleright$ (5)

3:initialize $\bm{\Pi}_{\theta}$

4:for $t=1,\dots,T$ do

5: $\bar{{\bm{X}}}\in{\mathbb{R}}^{D\times n_{b}}\leftarrow$ sample a batch
from ${\bm{X}}$

6: $\bar{{\bm{X}}}^{(1)},\dots,\bar{{\bm{X}}}^{(A)}\leftarrow$ apply $A$
augmentations to $\bar{{\bm{X}}}$

7: $\bar{{\bm{Z}}}_{\theta},\bar{\bm{\Pi}}_{\theta}$ $\leftarrow$ forward pass
with $\\{\bar{{\bm{X}}}^{(a)}\\}_{a=1}^{A}$ and network parameters
$\bm{\theta}$ $\triangleright$ (6)

8: $\nabla_{\bm{\theta}}$(4) $\leftarrow$ backward pass with respect objective
(4)

9: $\bm{\theta}\leftarrow$update $\bm{\theta}$ using some optimizer on
$\nabla_{\bm{\theta}}$(4)

10:end for

11:run spectral clustering on $\bm{\Pi}_{\theta}$ to estimate labels
$\hat{{\bm{y}}}$ of samples

12:${\bm{Z}}_{\theta},\hat{{\bm{y}}}$

##  3 Experiments on Real Datasets

Metrics. To evaluate the clustering quality, we run spectral clustering on
learned membership matrix $\bm{\Pi}$, and report the normalized mutual
information (NMI, [42]) and clustering accuracy (ACC, [21]), as are commonly
used in clustering tasks. To evaluate the learned representation, we define
the following metric: for a collection of points
${\bm{W}}=[{\bm{w}}_{1},\ldots,{\bm{w}}_{l}]\in{\mathbb{R}}^{d\times
l}\,\,(l>d)$ with associated singular values $\\{\sigma_{i}\\}_{i=1}^{d}$,
define the numerical rank of ${\bm{W}}$ as
$\operatorname*{arg\,min}_{r}\left\\{r:\sum_{i=1}^{r}\sigma_{i}^{2}/\sum_{i=1}^{d}\sigma_{i}^{2}>0.95\right\\}$.
Now, one can measure the numerical rank of the learned representaion
${\bm{Z}}$, as well as that of each ground-truth cluster999They are defined by
the true labels ${\bm{y}}$ (§2), so that the numerical rank metric is
decoupled from the quality of learned membership $\bm{\Pi}$. of ${\bm{Z}}$. A
low numerical rank of ${\bm{W}}$ implies that points in ${\bm{W}}$ lie close
to a low-dimensional subspace. We further report the cosine similarity of
learned representation, which is simply $|{\bm{z}}_{i}^{\top}{\bm{z}}_{j}|$
for points $i$ and $j$, since $\|{\bm{z}}_{i}\|=1$ by construction in (4).
Finally, to compare the efficiency of methods we report the training time in
§3.2, where the experiments are run on $2$ Nvidia RTX3090 GPUs.

###  3.1 Comparison with Subspace Clustering

To demonstrate the ability of MLC to cluster the samples and linearize the
manifolds, we conduct experiments on CIFAR10, which consists of RGB images
from $10$ classes such as planes, birds, and deers. As mentioned in §1
subspace clustering methods rely crucially on the assumption that data lie
close to a union of linear subspaces, which many real-world dataset may not
satisfy. To show that this is the case, we additionally compare the proposed
method with subspace clustering methods. As we shall see, applying subspace
clustering directly on self-supervised features of CIFAR10 will yield low
clustering accuracy. In contrast, MLC is able to achieve high clustering
accuracy, and moreover, produce a union-of-orthogonal-subspace representation
on which subspace clustering methods can also achieve high accuracy.

Data. We use the training split of CIFAR10 containing $50000$ RGB images, each
of size $3\times 32\times 32$. We use the augmentation specified in the
Appendix to perform self-supervised representation learning (5) and get
${\bm{Z}}_{\mathrm{self}}$. For a fair comparison, the so-learned
${\bm{Z}}_{\mathrm{self}}$ are used both as initialization for MLC (line 2 of
Algorithm 1), and as the input for subspace clustering methods101010The self-
supervised features ${\bm{Z}}_{\mathrm{self}}$ empirically exhibit some union-
of-subspace structure, and are typically used for subspace clustering, as also
seen in [53, §3.2] and [55, §4.2]. . In MLC, for each image in each batch we
randomly sample $A=2$ augmentations to apply on the image. As an additional
comparison, we also run subspace clustering methods on the features
${\bm{Z}}_{\mathrm{MLC}}$ learned by MLC.

Methods. We compare with the elastic-net subspace clustering with active-set
solver (EnSC, [50]) and sparse subspace clustering with orthogonal matching
pursuit solver (SSC-OMP, [51]), using off-the-shelf implementation provided by
the authors111111https://github.com/ChongYou/subspace-clustering. We search
the parameters of EnSC over
$(\gamma,\tau)\in\\{1,5,10,50,100\\}\times\\{0.9,0.95,1\\}$ and those of SSC
over
$(k_{\mathrm{max}},\epsilon)\in\\{3,5,10,20\\}\times\\{10^{-4},10^{-5},10^{-6},10^{-7}\\}$,
and report the run with the highest clustering accuracy for each method. We
summarize detailed parameters for MLC in the Appendix.

(a)

(b)

Figure 3: Coding rates (as loss terms in (4)) and numerical ranks (§3.1) of
the features learned by MLC on CIFAR10 as epoch varies.

Figure 4: Cosine similarity
$|{\bm{Z}}_{\mathrm{MLC}}^{\top}{\bm{Z}}_{\mathrm{MLC}}|$ of the features
${\bm{Z}}_{\mathrm{MLC}}$ learned by MLC.

Table 1: Clustering accuracy and normalized mutual information for subspace clustering (EnSC, SSC-OMP) on self-supervised features ${\bm{Z}}_{\mathrm{self}}$, features ${\bm{Z}}_{\mathrm{MLC}}$ learned by MLC, and manifold clustering (MLC) on ${\bm{X}}$, where ${\bm{X}}$ is $6\cdot 10^{4}$ images from $10$ classes of CIFAR10. Method |  |  Input Data |  ACC |  NMI
---|---|---|---|---
EnSC |  |  ${\bm{Z}}_{\mathrm{self}}$ |  72.2 |  67.9
|  ${\bm{Z}}_{\mathrm{MLC}}$ |  81.5 |  79.2
SSC-OMP |  |  ${\bm{Z}}_{\mathrm{self}}$ |  67.8 |  64.5
|  ${\bm{Z}}_{\mathrm{MLC}}$ |  78.4 |  76.3
MLC |  |  ${\bm{X}}$ |  86.3 |  78.3

Results. Figure 4 reports the coding rates (as loss terms in (4) and numerical
ranks of features learned by MLC as epoch varies. As a first note, the coding
rate $R$ of all features (the blue curve in LABEL:fig:cifar10-a) decreases
only slightly as epoch goes, indicating that the overall representation is
diverse in the feature space. Indeed, the numerical rank of all features (the
dark curve in LABEL:fig:cifar10-b) stays $118$ which is close to the dimension
$128$ of the feature space. This is in sharp contrast to the deep subspace
clustering methods where all the features collapse to a one-dimensional
subspace [15]. Moreover, as the coding rate $R_{c}$ of clustered features (the
orange curve in LABEL:fig:cifar10-a) goes down, the numerical ranks of
features from each ground-truth cluster decrease. For instance, the
representation from true cluster $3$ has a numerical rank of $37$ in the first
step and $24$ in the last step. This implies that most representation gets
linearized better and clustered more accurately, even though the MLC objective
(4) is unsupervised, i.e., it does not use ground-truth labels ${\bm{y}}$.
Last but not the least, note that the features within each ground-truth
cluster spread well in a low-dimensional subspace, e.g., the numerical ranks
for the true clusters at the last step are within $[13,23]$. This achieves the
desired within-cluster diverse property (§1), as opposed to the neural
collapse phenomenon that appears with the cross-entropy loss.

To compare MLC with subspace clustering methods, we report clustering accuracy
and normalized mutual information for EnSC, SSC-OMP on self-supervised
features ${\bm{Z}}_{\mathrm{self}}$, features ${\bm{Z}}_{\mathrm{MLC}}$
learned by MLC, and MLC on ${\bm{X}}$, where ${\bm{X}}$ is $6\cdot 10^{4}$
images from $10$ classes of CIFAR10. In addition we plot the cosine similarity
of the features learned by MLC in Figure 4. Remarkably, the highest clustering
accuracy is $86.3\%$ achieved by MLC on ${\bm{X}}$, which surpasses EnSC
($72.2\%$) and SSC-OMP ($67.8\%$) on ${\bm{Z}}_{\mathrm{self}}$ by a large
margin, even though ${\bm{Z}}_{\mathrm{self}}$ is used both as initialization
for MLC and input for EnSC and SSC-OMP. Interestingly, using instead the
features ${\bm{Z}}_{\mathrm{MLC}}$ learned by MLC, the clustering performance
of EnSC and SSC-OMP increases and even becomes comparable to MLC, e.g., EnSC
achieves $79.2\%$ normalized mutual information compared to $78.3\%$ of MLC.
This suggests that ${\bm{Z}}_{\mathrm{MLC}}$ has a union-of-subspace structure
that can be utilized by subspace clustering. Indeed, as seen in Figure 4,
features from different clusters tend to have a small similarity, i.e., being
orthogonal to each other. This demonstrates the between-cluster discrimination
(§1) as desired.

Table 2: Clustering accuracy and normalized mutual information on large scale datasets. For a fair comparison, all methods use ResNet-18 as backbone.  Method / Dataset |  | CIFAR10-10 | CIFAR100-20 | CIFAR100-100 | Tiny ImageNet-200
---|---|---|---|---|---
Metrics |  | ACC | NMI | ACC | NMI | ACC | NMI | ACC | NMI
DeepCluster ${}_{\mathrm{(ECCV^{\prime}18)}}$ |  | 37.4 | - | 18.9 | - | - | - | - | -
IIC ${}_{\mathrm{(ICCV^{\prime}19)}}$ |  | 61.7 | 51.1 | 25.7 | 22.5 | - | - | - | -
SCAN ${}_{\mathrm{(ECCV^{\prime}20)}}$ |  | 87.6 | 78.7 | 46.8 | 45.9 | 34.3 | 55.7 | - | -
RUC+SCAN ${}_{\mathrm{(CVPR^{\prime}21)}}$ |  | 90.3 | - | 53.3 | - | - | - | - | -
IMC-SWAV ${}_{\mathrm{(Arxiv^{\prime}21)}}$ |  | 89.1 | 81.1 | 49.0 | 50.3 | 43.9 | 58.3 | 28.2 | 52.6
NMCE ${}_{\mathrm{(Arxiv^{\prime}22)}}$ |  | 83.0 | 76.1 | 43.7 | 48.8 | - | - | - | -
MLC |  | 86.3 | 78.3 | 52.2 | 54.6 | 49.4 | 68.3 | 33.5 | 67.5
Table 3: Running time in minutes and clustering accuracy on CIFAR100-100. For a fair comparison, all methods use ResNet-18 as backbone. Method / Metric |  | Running Time | ACC
---|---|---|---
Stage |  | I | II | III | Total |
SCAN ${}_{\mathrm{(ECCV^{\prime}20)}}$ |  | 308.3 | 33.3 | 54.7 | 396.3 | 34.3
IMC-SWAV ${}_{\mathrm{(Arxiv^{\prime}21)}}$ |  | 529.4 | - | - | 529.4 | 43.9
MLC |  | 266.7 | 17.7 | - | 284.4 | 48.3
Table 4: Clustering accuracy on imbalanced datasets: (a) Imb-CIFAR10, (b) Imb-CIFAR100-100. For a fair comparison, all methods use ResNet-18 as backbone.  Method / Dataset |  |  (a) |  (b)
---|---|---|---
IMC-SWAV ${}_{\mathrm{(Arxiv^{\prime}21)}}$ |  |  65.7 |  38.2
SCAN ${}_{\mathrm{(ECCV^{\prime}20)}}$ |  |  62.9 |  31.1
MLC |  |  80.0 |  46.1

###  3.2 Comparison with Deep Clustering Methods

We further compare the proposed MLC with state-of-the-art deep clustering
methods. Note that most methods reported (all except NMCE which is discussed
in §2.2) do not aim to learn a union-of-orthogonal-subspace representation, in
contrast to MLC. As we will see, MLC achieves clustering accuracy comparable
to state-of-the-art methods on large scale datasets with faster computational
time, and further surpasses them on extreme yet realistic cases like datasets
of imbalanced clusters.

Compared Methods. We conduct experiments with MLC, SCAN [45], and IMC-SWAV
[31].121212The authors are aware of a preprint [30] which approaches image
clustering via a combination of self/semi-supervised learning and pseudo-
labeling. However, to the best of our effort we are unable to reproduce the
numbers reported in this paper using the implementation provided by the
authors. We discuss the details in the Appendix and thus do not report their
numbers here. Training details can be found in the Appendix. In addition we
include the numbers reported from DeepCluster [7], IIC [19], RUC [33] and NMCE
[23]. For a fair comparison, all methods reported use ResNet-18 as the
backbone, which is also commonly adopted by other methods.

Datasets. Beyond CIFAR10 (§3.1), we further use CIFAR100-20, CIFAR100-100 and
Tiny Imagenet-200 to evaluate the performance of our method. Both CIFAR100-100
and CIFAR100-20 contain the same $50000$ train images and $10000$ test images
with size $32\times 32\times 3$, while the former are split into $100$
clusters and the latter $20$ super clusters. Finally, Tiny ImageNet contains
$100000$ train images and $10000$ test images with size $64\times 64\times 3$
split into $200$ clusters.

Results on Large-scale Datasets. We report clustering accuracy and normalized
mutual information on CIFAR10, CIFAR100-20, CIFAR100-100, and TinyImageNet in
Table 2, and we further report running time in minutes for CIFAR100-100 in
Table 3. As seen, the highest clustering performance on CIFAR10 is achieved by
RUC+SCAN ($90.3\%$ ACC) and IMC-SWAV ($81.1\%$ NMI), where MLC yields a
slightly lower ACC of $86.3\%$ and NMI of $78.3\%$. We note some interesting
semantic interpretation for the clustering obtained by MLC in the Appendix. On
the other hand, MLC performs comparably with other methods on CIFAR100-20 by
achieving an ACC of $52.2\%$ and NMI of $54.6\%$. Notably, MLC outperforms
SCAN and IMC-SWAV on CIFAR100-100 and TinyImageNet-200 by a large margin,
while using lower running time: E.g., on CIFAR100-100, MLC yields an accuracy
of $49.4\%$ in $291$ minutes, whereas IMC-SWAV has $43.9\%$ using $529$
minutes, and SCAN has $34.3\%$ in $396$ minutes.

Imbalanced Clusters. Note that for CIFAR10 or CIFAR100 each cluster contains
approximately the same number of samples. On the other hand, natural images
are typically imbalanced, i.e., the clusters have unequal number of samples.
To mimic this setting, we take a naive approach to construct the following
imbalanced datasets. For the $10$ clusters of CIFAR10, we remove half of the
samples from odd-numbered clusters (i.e., clusters $1,3,\dots,9$) from both
the training and test split. We refer to the reduced dataset Imb-CIFAR10.
Likewise we construct Imb-CIFAR100-100. We run two state-of-the-art methods
IMC-SWAV and SCAN as well as the proposed MLC on Imb-CIFAR10 and Imb-
CIFAR100-100.

Table 4 shows clustering accuracy on the imbalanced datasets Imb-CIFAR10 and
Imb-CIFAR100-100. As a first observation, the clustering accuracy of all
methods is lower on the imbalanced datasets than on the balanced counterparts,
as expected. Notably, MLC suffers from the least performance drop, e.g., when
moving from CIFAR10 to Imb-CIFAR10 the accuracy of MLC drops from $86\%$ to
$80\%$, whereas that of SCAN and IMC-SWAV decreases from above $87\%$ to below
$66\%$.

##  4 Conclusion

This paper studies the problem of simultaneously clustering and learning an
union-of-orthogonal-subspace representation for data, when data lies close to
a union of low-dimensional manifolds. To address the problem we propose an
objective based on maximal coding rate reduction and doubly stochastic
membership inspired by the state-of-the-art subspace clustering results. We
provide an efficient and effective parameterization of the membership
variables as well as a meta-algorithm to optimize the representation and
membership jointly. We further conduct experiments on datasets with larger
number of clusters and imbalanced clusters and show that the proposed method
achieves state-of-the-art performance. We believe that our work provides a
general and unified framework for unsupervised learning of structured
representations for multi-modal data.

## References

  * [1] Mahdi Abavisani and Vishal M Patel.  Deep multimodal subspace clustering networks.  IEEE Journal of Selected Topics in Signal Processing, 12(6):1601–1614, Apr. 2018.
  * [2] David Arthur and Sergei Vassilvitskii.  k-means++: The advantages of careful seeding.  In the Eighteenth Annual ACM-SIAM Symposium on Discrete Algorithms. Society for Industrial and Applied Mathematics, June 2006.
  * [3] Bahman Bahmani, Benjamin Moseley, Andrea Vattani, Ravi Kumar, and Sergei Vassilvitskii.  Scalable K-Means++.  Proceedings VLDB Endowment, 5(7), Mar. 2012.
  * [4] Adrien Bardes, Jean Ponce, and Yann LeCun.  VICReg: Variance-Invariance-Covariance regularization for Self-Supervised learning.  In International Conference on Learning Representations, 2022.
  * [5] Paul Bradley, Olvi Mangasarian, and W Street.  Clustering via concave minimization.  In Advances in neural information processing systems, 1996.
  * [6] Joan Bruna and Stéphane Mallat.  Invariant scattering convolution networks.  IEEE transactions on pattern analysis and machine intelligence, 35(8):1872–1886, Aug. 2013.
  * [7] Mathilde Caron, Piotr Bojanowski, Armand Joulin, and Matthijs Douze.  Deep clustering for unsupervised learning of visual features.  In European conference on computer vision, pages 132–149, July 2018\.
  * [8] Ting Chen, Simon Kornblith, Kevin Swersky, Mohammad Norouzi, and Geoffrey E Hinton.  Big self-supervised models are strong semi-supervised learners.  Advances in neural information processing systems, 33:22243–22255, 2020.
  * [9] Tianjiao Ding, Derek Lim, Rene Vidal, and Benjamin D Haeffele.  Understanding doubly stochastic clustering.  In the 39th International Conference on Machine Learning, volume 162 of Proceedings of Machine Learning Research, pages 5153–5165. PMLR, 2022.
  * [10] Marvin Eisenberger, Aysim Toker, Laura Leal-Taixé, Florian Bernard, and Daniel Cremers.  A unified framework for implicit sinkhorn differentiation.  In IEEE/CVF Conference on Computer Vision and Pattern Recognition, pages 509–518, 2022.
  * [11] Ehsan Elhamifar and Rene Vidal.  Sparse subspace clustering.  In IEEE Conference on Computer Vision and Pattern Recognition, pages 2790–2797, June 2009.
  * [12] Ehsan Elhamifar and René Vidal.  Sparse manifold clustering and embedding.  Advances in neural information processing systems, 24, 2011.
  * [13] Ehsan Elhamifar and Rene Vidal.  Sparse subspace clustering: Algorithm, theory, and applications.  IEEE transactions on pattern analysis and machine intelligence, 35(11):2765–2781, 2013.
  * [14] Edward Forgey.  Cluster analysis of multivariate data: Efficiency vs. interpretability of classification.  Biometrics, 1965.
  * [15] Benjamin D Haeffele, Chong You, and René Vidal.  A critique of Self-Expressive deep subspace clustering.  In International Conference on Learning Representations, 2020.
  * [16] Reinhard Heckel and Helmut Bölcskei.  Robust subspace clustering via thresholding.  IEEE transactions on information theory, 61(11):6320–6342, 2015\.
  * [17] R C Jancey.  Multidimensional group analysis.  Australian Journal of Botany, 14:127–130, 1966.
  * [18] Pan Ji, Tong Zhang, Hongdong Li, Mathieu Salzmann, and Ian Reid.  Deep subspace clustering networks.  In Advances in Neural Information Processing Systems, volume 2017-Decem, pages 24–33, 2017.
  * [19] Xu Ji, João F Henriques, and Andrea Vedaldi.  Invariant information clustering for unsupervised image classification and segmentation.  In IEEE/CVF International Conference on Computer Vision, July 2018\.
  * [20] Mohsen Kheirandishfard, Fariba Zohrizadeh, and Farhad Kamangar.  Deep low-rank subspace clustering.  In 2020 IEEE/CVF Conference on Computer Vision and Pattern Recognition Workshops. IEEE, June 2020.
  * [21] Minsik Lee, Jieun Lee, Hyeogjin Lee, and Nojun Kwak.  Membership representation for detecting block-diagonal structure in low-rank or sparse subspace clustering.  Proceedings of the IEEE Computer Society Conference on Computer Vision and Pattern Recognition, 07-12-June:1648–1656, 2015.
  * [22] Chun-Guang Li, Chong You, and René Vidal.  On geometric analysis of affine sparse subspace clustering.  IEEE J. Sel. Top. Signal Process., 12(6):1520–1533, Dec. 2018.
  * [23] Zengyi Li, Yubei Chen, Yann LeCun, and Friedrich T Sommer.  Neural manifold clustering and embedding.  arXiv [cs.LG], Jan. 2022.
  * [24] Derek Lim, René Vidal, and Benjamin D Haeffele.  Doubly stochastic subspace clustering.  arXiv [cs.LG], Nov. 2020.
  * [25] Guangcan Liu, Zhouchen Lin, Shuicheng Yan, Ju Sun, Yong Yu, and Yi Ma.  Robust recovery of subspace structures by low-rank representation.  IEEE transactions on pattern analysis and machine intelligence, 35(1):171–184, Jan. 2013.
  * [26] Stuart Lloyd.  Least squares quantization in PCM.  Technical report, Bell Laboratories, 1957.
  * [27] Can-Yi Lu, Hai Min, Zhong-Qiu Zhao, Lin Zhu, De-Shuang Huang, and Shuicheng Yan.  Robust and efficient subspace segmentation via least squares regression.  In European conference on computer vision, pages 347–360. Springer, 2012.
  * [28] Yi Ma, Harm Derksen, Wei Hong, and John Wright.  Segmentation of multivariate mixed data via lossy data coding and compression.  IEEE transactions on pattern analysis and machine intelligence, 29(9):1546–1562, 2007.
  * [29] James B McQueen.  Some methods for classification and analysis of multivariate observations.  In Fifth Berkeley Symposium on Mathematical Statistics and Probability, pages 281–297, 1967.
  * [30] Chuang Niu, Hongming Shan, and Ge Wang.  SPICE: Semantic pseudo-labeling for image ClustEring.  arXiv [cs.CV], Mar. 2021.
  * [31] Foivos Ntelemis, Yaochu Jin, and Spencer A Thomas.  Information maximization clustering via Multi-View Self-Labelling.  arXiv [cs.CV], Mar. 2021.
  * [32] Vardan Papyan, X Y Han, and David L Donoho.  Prevalence of neural collapse during the terminal phase of deep learning training.  Proceedings of the National Academy of Sciences of the United States of America, 117(40):24652–24663, Oct. 2020.
  * [33] Sungwon Park, Sungwon Han, Sundong Kim, Danu Kim, Sungkyu Park, Seunghoon Hong, and Meeyoung Cha.  Improving unsupervised image clustering with robust learning.  In IEEE/CVF Conference on Computer Vision and Pattern Recognition, pages 12278–12287, 2021.
  * [34] Vishal M Patel and Rene Vidal.  Kernel sparse subspace clustering.  In 2014 IEEE International Conference on Image Processing (ICIP). IEEE, Oct. 2014.
  * [35] Xi Peng, Jiashi Feng, Shijie Xiao, Jiwen Lu, Zhang Yi, and Shuicheng Yan.  Deep sparse subspace clustering.  arXiv [cs.CV], Sept. 2017.
  * [36] Herbert Robbins and Sutton Monro.  A stochastic approximation method.  Ann. Math. Stat., 22(3):400–407, 1951.
  * [37] Daniel P Robinson, Rene Vidal, and Chong You.  Basis pursuit and orthogonal matching pursuit for subspace-preserving recovery: Theoretical analysis.  arXiv [cs.LG], Dec. 2019.
  * [38] Aurko Roy, Mohammad Saffar, Ashish Vaswani, and David Grangier.  Efficient content-based sparse attention with routing transformers.  Transactions of the Association for Computational Linguistics, 9:53–68, Feb. 2021.
  * [39] Michael E Sander, Pierre Ablin, Mathieu Blondel, and Gabriel Peyré.  Sinkformers: Transformers with doubly stochastic attention.  In International Conference on Artificial Intelligence and Statistics, Oct. 2021.
  * [40] Mahdi Soltanolkotabi and Emmanuel J Candés.  A geometric analysis of subspace clustering with outliers.  Annals of statistics, 40(4):2195–2238, 2012.
  * [41] Mahdi Soltanolkotabi, Ehsan Elhamifar, and Emmanuel J Candes.  Robust subspace clustering.  Annals of statistics, 42(2):669–699, 2014.
  * [42] Strehl and Ghosh.  Cluster ensembles – a knowledge reuse framework for combining multiple partitions.  Journal of machine learning research, 2002.
  * [43] Tom Tirer and Joan Bruna.  Extended unconstrained features model for exploring deep neural collapse.  In International Conference on Machine Learning, Feb. 2022.
  * [44] Manolis Tsakiris and Rene Vidal.  Theoretical analysis of sparse subspace clustering with missing entries.  In Jennifer Dy and Andreas Krause, editors, International Conference on Machine Learning, pages 4975–4984, 2018.
  * [45] Wouter Van Gansbeke, Simon Vandenhende, Stamatios Georgoulis, Marc Proesmans, and Luc Van Gool.  SCAN: Learning to classify images without labels.  In European conference on computer vision. Springer, 2020.
  * [46] Ulrike von Luxburg.  A tutorial on spectral clustering.  Statistics and computing, 17(4):395–416, 2007.
  * [47] Yining Wang, Yu-Xiang Wang, and Aarti Singh.  A deterministic analysis of noisy sparse subspace clustering for dimensionality-reduced data.  In International Conference on Machine Learning, 2015.
  * [48] Yu-Xiang Wang and Huan Xu.  Noisy sparse subspace clustering.  Journal of Machine Learning Research, 2016.
  * [49] Chong You, Chun-Guang Li, Daniel Robinson, and Rene Vidal.  Is an affine constraint needed for affine subspace clustering?  In International Conference on Computer Vision (ICCV). IEEE, Oct. 2019.
  * [50] Chong You, Chun Guang Li, Daniel P Robinson, and Rene Vidal.  Oracle based active set algorithm for scalable elastic net subspace clustering.  In the IEEE conference on computer vision and pattern recognition, pages 3928–3937, 2016.
  * [51] Chong You, Daniel P Robinson, and René Vidal.  Scalable sparse subspace clustering by orthogonal matching pursuit.  In IEEE Conference on Computer Vision and Pattern Recognition, June 2016.
  * [52] Yang You, Igor Gitman, and Boris Ginsburg.  Large batch training of convolutional networks.  arXiv [cs.CV], Aug. 2017.
  * [53] Yaodong Yu, Kwan Ho Ryan Chan, Chong You, Chaobing Song, and Yi Ma.  Learning diverse and discriminative representations via the principle of maximal coding rate reduction.  In Neural Information Processing Systems, June 2020.
  * [54] Junjian Zhang, Chun-Guang Li, Chong You, Xianbiao Qi, Honggang Zhang, Jun Guo, and Zhouchen Lin.  Self-Supervised convolutional subspace clustering network.  In IEEE/CVF conference on computer vision and pattern recognition, 2019.
  * [55] Shangzhi Zhang, Chong You, René Vidal, and Chun-Guang Li.  Learning a self-expressive network for subspace clustering.  In IEEE Conference on Computer Vision and Pattern Recognition, Oct. 2021.
  * [56] Jinxin Zhou, Xiao Li, Tianyu Ding, Chong You, Qing Qu, and Zhihui Zhu.  On the optimization landscape of neural collapse under MSE loss: Global optimality with unconstrained features.  In International Conference in Machine Learning, Mar. 2022.

(a) Learned cluster 1

(b) Learned cluster 2

(c) Learned cluster 3

(d) Learned cluster 4

(e) Learned cluster 5

(f) Learned cluster 6

(g) Learned cluster 7

(h) Learned cluster 8

(i) Learned cluster 9

Figure 5: Principal images (defined in §A) of clusters learned by MLC on
CIFAR10.

##  Appendix A Semantic Interpretability of the Learned Representation and
Clusters on CIFAR10

Recall that MLC is designed to perform clustering while learning a union-of-
orthogonal-subspace representation (§1), where each cluster defines a low-
dimensional subspace. Therefore, we further visualize the different directions
within each learned cluster or subspace. Specifically, after a final
clustering is obtained (line 11 of Algorithm 1), we take the features from
each learned cluster and apply Principal Component Analysis (PCA) to them to
obtain the first $8$ principal components. These correspond to the $8$ rows
for each cluster in Figure 5. Recall that the principal components are
mutually orthogonal, indicating uncorrelated directions within one cluster. To
visualize those directions or principal components in images, we take the
features that are closest to the principal components and visualize the
corresponding original images.

Interestingly, the rows of images corresponding to principal components appear
to exhibit some semantic ‘concepts’. For example in Figure 5(b), row $1$ and
$8$ are respectively white and red trucks, while row $4$ are the trucks that
ship sand or mud; row $1$ of Figure 5(d) are deers with trees as background.
This further suggests that the learned embedding seems to preserve distance
within each cluster (as desired in §1), i.e., images that are close/far in
semantic meaning will be close/far in the feature space. Note however, that
some learned clusters do not align fully with the ground-truth labels. For
instance, rows $1$ and $3$ of Figure 5(h) are cats while all other rows in
this cluster are dogs. On the other hand, one may argue that Figure 5(h) are a
cluster of cats and dogs of light colors, whereas Figure 5(c) is a cluster of
those of brown colors, which could be a semantically meaningful clustering
even though it does not align with the ground-truth labels. We believe it
would be an interesting future work to use MLC to discover new semantics that
are not present in the given labels.

Table 5: Ablation study on the roles of different parts of Algorithm 1 and on
using augmentation.

Ablation Study on CIFAR-10 |  | Clustering
---
Accuracy
Full Algorithm 1 | 86.3%
| Replacing self-supervised initialization (line 2)
---
with random initialization
20.0%
| Replacing updating MLC objective (4) (lines 4-10)
---
with subspace clustering (EnSC)
73.4%
Not using augmentation in lines 4-10 | 80.0%

##  Appendix B Role of Augmentation

Recall that data augmentation was used both in the self-supervised
initialization (line 2 of Algorithm 1, see ‘Initializing
${\bm{Z}}_{\bm{\theta}}$’ in §2.3) and in updating the MLC objective (lines
4-10, ‘Data Augmentation’ in §2.3). Below we give additional clarification on
the role of augmentation therein.

###  B.1 Augmentation for Initializing the Features

Since the proposed MLC objective (4) is highly non-convex, the (local)
solution that a first-order optimizer converges to in general depends on the
initialization. However, before line 2 of Algorithm 1 is executed, the
features ${\bm{Z}}$ at initialization could be very far from union-of-
orthogonal-subspace (as desired by Problem 1), since the neural network has an
arbitrary architecture and initialization. To at least promote some ideal
structures in the features, we conduct line 2 of Algorithm 1 so that the
features from an original sample and its augmented copy are close, while
features from different samples spread out in the feature space. This is a
common idea used in contrastive learning, and more related in [53, §3.2] [23,
§3.6] that are both based on MCR2 as in this paper (even though the
formulations are different, as argued in §1.2). Empirically, initializing the
features using augmentation (line 2 of Algorithm 1) is important for the
following-up steps: as seen in Table 5, on CIFAR-10, if one uses random
initialization to replace this step, then the final clustering accuracy is
$20\%$, in sharp contrast to $86.3\%$.

###  B.2 Augmentation for Updating MLC Objective (4)

In optimizing MLC (lines 4-10 of Algorithm 1), augmentation empirically
improves clustering performance. As one can see in Table 5, on CIFAR-10 using
the sample self-supervised initialization of the features, MLC achieves only
$80\%$ clustering accuracy without augmentation, in contrast to $86.3\%$ with
augmentation. We attribute this difference to the fact that augmentation
enriches the diversity of samples the algorithm sees.

##  Appendix C Role of Different Parts of Algorithm 1

###  C.1 Initialization of the Features (Line 2)

Please kindly refer to §B.1.

###  C.2 Updating the MLC Objective (4) (Lines 4-10)

The main novelty of this paper lies in updating the MLC objective that learns
both the representation ${\bm{Z}}_{\bm{\theta}}$ and a doubly stochastic
membership $\bm{\Pi}_{\bm{\theta}}$. Note that in this step, clustering is
pursued by modeling the membership $\bm{\Pi}_{\bm{\theta}}$, as opposed to the
self-supervised feature initialization step where no membership is explicitly
pursued. This step is indeed important for clustering: as seen in Table 5, on
CIFAR-10, the clustering accuracy on the self-supervised initialized features
${\bm{Z}}_{\bm{\theta}}$ is only $73.4\%$, in contrast to $86.3\%$ obtained
after updating the MLC objective (4).

###  C.3 Spectral Clustering (Line 11)

Since the proposed MLC learns a doubly stochastic membership that signals
pair-wise similarity between points, it is standard to run spectral clustering
[46] to compute a final set of clusters from the learned membership. This is
done once at the very end of Algorithm 1, and is rather efficient compared to
the other parts of Algorithm 1: for instance, using an unaccelerated
implementation from SciPy, it takes less than $30$ seconds to perform spectral
clustering on a $10000\times 10000$ matrix.

##  Appendix D Details on Experiment Settings

###  D.1 Synthetic Union-of-Manifold Data

We perform simulations to visualize the properties of the proposed manifold
learning and clustering method. As seen in Figure 1(a), we generate data
${\bm{X}}$ from two manifolds on the sphere ${\mathbb{S}}^{2}$, each
consisting of $200$ samples. The points from the first manifold (green) take
the form

$\displaystyle{\bm{x}}_{i}=\begin{bmatrix}\cos\big{(}A\sin(\omega\phi_{i})\big{)}\cos\phi_{i}\\\
\cos\big{(}A\sin(\omega\phi_{i})\big{)}\sin\phi_{i}\\\
\sin\big{(}A\sin(\omega\phi_{i})\big{)}\end{bmatrix}+\bm{\epsilon}_{i},$ (7)

where $A=0.2$ and $\omega=5$ sets the curvature of the manifold,
$\bm{\epsilon}_{i}\sim{\mathcal{N}}({\bm{0}},0.05{\bm{I}}_{3})$ is the
additive noise, and we take $\phi_{i}=\frac{2\pi i}{100}$ for $i=1,\dots,100$
to generate $100$ points. On the other hand, the points from the second
manifold (blue) are simply $100$ samples from
${\mathcal{N}}([0,0,1]^{\top},0.05{\bm{I}}_{3})$. We take the feature
dimension $d=3$ to be equal to he input dimension $D=3$. We paramterize both
the feature head $f_{\bm{\theta}}$ and the cluster head $g_{\bm{\theta}}$ to
be a simple fully-connected network with $100$ hidden neurons, followed by a
Rectified Linear Unit as non-linearity and a projection operator onto the
sphere ${\mathbb{S}}^{2}$. Figures 1(b), 1(c) and 1(d) report the features
${\bm{Z}}_{\bm{\theta}}$ with random initialization (i.e., before line 2 of
Algorithm 1), with self-supervised initialization, and at convergence of MLC.
Notably, despite ${\bm{Z}}_{\bm{\theta}}$ being noisy and only approximately
piece-wise linear, as epoch goes ${\bm{Z}}_{\bm{\theta}}$ gradually transform
to two linear subspaces: the green points converge to a $2$-dimensional
subspace (intersected with ${\mathbb{S}}^{2}$) and the blue points converge to
a $1$-dimension subspace.

###  D.2 Training Details on Real Datasets

MLC. As said, we use ResNet-18 as the backbone for experiments on CIFAR10,
CIFAR100-20, CIFAR100-100 and Tiny-ImageNet-200, and the imbalanced
counterparts Imb-CIFAR10, Imb-CIFAR100-100. We also fix the batch size to be
$1024$ in all experiments. In self-supervised initialization of
${\bm{Z}}_{\bm{\theta}}$ (line 2 of Algorithm 1), we use the precision (§2.1)
parameter $\epsilon^{2}=0.2$, a LARS optimizer [52] (as is also done in [8,
23]) with a learning rate of $0.3$ and trained MLC for $1000$ epochs. On the
other hand, in the training of MLC objective, we use $\epsilon^{2}=0.1$,
$\gamma=0.05$, and $\eta=0.175$ for the entropy regularization in the Sinkhorn
projection [10] layer $P_{\Omega,\eta}(\cdot)$. We fix the backbone and for
each batch, we perform one update for parameters in the feature head
${\bm{Z}}_{\bm{\theta}}$ and one update for parameters in the cluster head
${\bm{C}}_{\bm{\theta}}$. For each head we use one SGD optimizer [36] with a
learning rate of $10^{-2}$, momentum of $0.9$, and weight decay of $5\cdot
10^{-4}$. Finally, for all experiments, we use the augmentation from [4]
detailed below in PyTorch code.

Augmentation 1 Augmentations for real datasets

import torchvision.transforms as t

t.Compose([

t.RandomResizedCrop(32,scale=(0.04, 1.0)),

t.RandomHorizontalFlip(p=0.5),

t.RandomGrayscale(p=0.2),

t.RandomApply([t.ColorJitter(0.4, 0.4, 0.4, 0.1)], p=0.8),

GaussianBlur(p=0.1)

])

SCAN and IMC-SWAV. Recall that we conduct experiments on CIFAR100-100, Imb-
CIFAR10, and Imb-CIFAR100-100 with SCAN [45], IMC-SWAV [31] and MLC, and
report clustering and running time in Tables 3 and 4. We use off-the-shelf
implementation131313https://github.com/wvangansbeke/Unsupervised-
Classification, https://github.com/foiv0s/imc-swav-pub provided by the
authors. For a fair comparison, SCAN, IMC-SWAV and MLC all use ResNet-18 as
the backbone. Finally, the hyper-parameters of SCAN and IMC-SWAV are set to be
the ones optimally chosen for CIFAR10 and CIFAR100 respectively provided by
the authors.

SPICE. As mentioned, the preprint [30] proposed a method SPICE that appears to
achieve state-of-the-art performance in image clustering. We tried to
reproduce their results on CIFAR-100-20 using the official
implementation141414https://github.com/niuchuangnn/SPICE, commit 5eba538..
However, the provided implementation ran into a few errors, which are also
observed151515https://github.com/niuchuangnn/SPICE/issues/27,https://github.com/niuchuangnn/SPICE/issues/31.
Despite our best effort to fix those issues, the experiments yield only $14\%$
clustering accuracy on CIFAR100-20 as opposed to the $53\%$ reported in the
paper [30]. Therefore, we note this observation and do not include SPICE in
the main text.

Table 6: Clustering accuracy and normalized mutual information of MLC and NMCE
on CIFAR10 over $10$ random seeds, using the same self-initialized features.
For the purpose of comparison, both methods use the same optimization strategy
and hyper-parameters optimally tuned for NMCE. Consequently, the clustering
performance of MLC reported here is lower than that in Table 2.

Method |  | Metric | Seed |  | Mean | Std.
---|---|---|---|---|---|---
| 0 | 1 | 2 | 3 | 4 | 5 | 6 | 7 | 8 | 9 |
MLC |  | ACC | 84.5 | 84.8 | 84.8 | 84.6 | 84.4 | 84.4 | 84.0 | 84.3 | 84.4 | 84.6 |  | 84.5 | 0.24
| NMI | 76.6 | 77.1 | 76.8 | 76.8 | 76.5 | 76.4 | 76.1 | 76.4 | 76.4 | 76.5 |  | 76.6 | 0.28
NMCE |  | ACC | 83.7 | 82.1 | 81.6 | 73.7 | 80.4 | 77.9 | 81.7 | 81.4 | 72.7 | 80.9 |  | 79.6 | 3.69
| NMI | 74.4 | 71.2 | 70.4 | 65.2 | 70.0 | 68.1 | 72.7 | 70.8 | 69.2 | 69.8 |  | 70.2 | 2.49

(a) CIFAR100-20

(b) CIFAR100-100

(c) Tiny ImageNet-200

Figure 6: Cosine similarity
$|{\bm{Z}}_{\mathrm{MLC}}^{\top}{\bm{Z}}_{\mathrm{MLC}}|$ of the features
${\bm{Z}}_{\mathrm{MLC}}$ learned by MLC on more complicated datasets:
CIFAR100-20, CIFAR100-100, Tiny ImageNet-200.

##  Appendix E Additional Comparison of MLC and NMCE on Stability with Respect
to Random Seeds

As detailed in §2.2, one of the advantages of the proposed MLC over NMCE [23]
is that MLC has a more stable performance with respect to random seeds, since
MLC is able to initialize the membership deterministically using structures
from the self-supervised initialized features. Below we conduct extra
experiments to provide empirical evidence. We first fix a self-supervised
initialization of features that is in turn used for both NMCE and MLC. Then,
using this very same initialization of features, we update NMCE and MLC
objective respectively with $5$ different seeds: recall that NMCE initializes
the membership randomly whereas MLC initializes the membership
deterministically using the initialized features. To make a valid comparison,
for both methods we further use the same optimization strategy and hyper-
parameters that are optimally161616For NMCE, we use the implementation as well
as the parameters provided in https://github.com/zengyi-li/NMCE-release. tuned
for NMCE (which are not optimal for MLC): precision $\epsilon^{2}=0.2$, #
epochs $100$, LARS optimizer for ${\bm{Z}}_{\bm{\theta}}$ with an initial
learning rate $0.3$ decayed to $0$ in a cosine annealing manner. Table 6
reports clustering accuracy and normalized mutual information of MLC and NMCE
over $10$ random seeds. As expected, MLC has a more stable clustering
performance by having a standard deviation of clustering accuracy and
normalized mutual information less than $0.28$, in contrast to more than
$2.49$ achieved by NMCE. Further, MLC achieves higher mean clustering
performance than NMCE, as also observed in Table 2. Last but not least, we
note that the numbers in Table 6 are not comparable to those in Table 2, since
for MLC the hyper-parameters and optimizers are different, and for NMCE an
additional step that fine tunes the backbone is used in Table 2.

##  Appendix F Addtional Visualization on Learned Representation and Clusters

Figure 6 presents the cosine similarity (as defined in the preamble of §3) of
the representation learned by MLC on CIFAR100-20, CIFAR100-100 and
TinyImageNet-200 (for the counterpart on CIFAR10 see §3.1). As seen, the
cosine similarity maps form approximately block diagonal structures, showing
that the features from different clusters are roughly orthogonal to each
other. This is desired by the between-cluster discrimination (§1).

Finally, we provide additional visualization of principal images on
CIFAR100-20 (see §A for definition) in Figure 8.

(a) Learned cluster 1

(b) Learned cluster 2

(c) Learned cluster 3

(d) Learned cluster 4

(e) Learned cluster 5

(f) Learned cluster 6

(g) Learned cluster 7

(h) Learned cluster 8

(i) Learned cluster 9

(j) Learned cluster 10

(k) Learned cluster 11

(l) Learned cluster 12

(m) Learned cluster 13

(n) Learned cluster 14

(o) Learned cluster 15

(p) Learned cluster 16

(q) Learned cluster 17

(r) Learned cluster 18

(s) Learned cluster 19

(t) Learned cluster 20

Figure 8: Principal images (defined in §A) of clusters learned by MLC on
CIFAR100-20.