Fusion over the Grassmann Manifold for Incomplete-Data Clustering
Authors: Jeremy Scott Johnson
Huanran Li
Chun Gan
Zheng Lu
Matthew Malloy
Daniel L. Pimentel-Alarcón
Venue: Submitted to ICLR 2023 (2023)
Keywords: high-rank matrix completion
subspace clustering
manifold learning

Abstract: This paper presents a new paradigm to cluster incomplete vectors using subspaces as proxies to exploit the geometry of the Grassmannian. We leverage this new perspective to develop an algorithm to cluster and complete data in a union of subspaces via a fusion penalty formulation. Our approach does not require prior knowledge of the number of subspaces, is naturally suited to handle noise, and only requires an upper bound on the subspaces’ dimensions. In developing our model, we present local convergence guarantees. We describe clustering, completion, model selection, and sketching techniques that can be used in practice, and complement our analysis with synthetic and real-data experiments.

Under review as a conference paper at ICLR 2023
FUSION OVER THE GRASSMANN MANIFOLD FOR INCOMPLETE-DATA CLUSTERING
Anonymous authors Paper under double-blind review
ABSTRACT
This paper presents a new paradigm to cluster incomplete vectors using sub- spaces as proxies to exploit the geometry of the Grassmannian. We leverage this new perspective to develop an algorithm to cluster and complete data in union of subspaces via fusion penalty formulation. Our approach does not require prior knowledge of the number subspaces, is naturally suited to handle noise, and only requires upper bound the subspaces* dimensions. In developing our model, we present local convergence guar- antees. We describe clustering, completion, model selection, and sketching techniques that can be used in practice, and complement our analysis with synthetic and rcal-data experiments.
INTRODUCTION
Supposc we observe subset of entries in a data matrix whose columns lie near union of subspaces, for example:
4
6
9
16
8
*
5
4
*
16
5
*
*
*
7
6
12
18
-1
2
&
5
*
*
1
0
-1
5
6
19
9
5
2
-3
8
1
* 7
14
*
13
8
where the unobserved entries are marked with *. Our goals are (i) to complete the unobserved entries, (ii) to cluster the columns according to the subspaces, and (zit) to learn the underlying subspaces. In the example above, we should (i) obtain the following (ground truth) completion:
4
16
-1
-7
1
5
4
14
16
5
-7
-18
10
8
7
-12
2
18
-1
28
- -18
X
2
1
11
4
5
1
9
19
9
-3
1
-1
3
14
13
we should also (ii) cluster the columns of X into two groups, X2, x7) and (xx, X4, Xs, Xs, Xg), and (iei) obtain bases for two 2-dimensional subspaces (given by any subset of linearly independent columns from each group).
This problem is often known as high-runk matrir completion (HRMC) [24, 21| or as subspace clustering with missing dala, and it has a wide range of applications, including tracking moving objects in computer vision [13, 14, 30, 31, 33, 35, 42]. predicting target interactions for drug discovery [28, 44, 45, 49), and identifying groups in recommender systems [37, 66, 77]. While there exists theory detailing conditions under which the HRMC goals above are feasible (e.g., sufficient sampling and subspaces genericity) [51], existing algorithms present a variety of shortcomings (more details in Section below).
The fundamental difficulty that all HRMC approaches face lies in assessing distances (e.g., euclidean, in the form of inner products) between partially observed vectors, for

---

Under review as a conference paper at ICLR 2023
the simple reason that this requires overlapping observations which become increasingly unlikely in low-sampling regimes [24]. To circumvent this problem, we introduce a new paradigm to cluster incomplete vectors, using subspaces as proxies, thus avoiding the need to calculate distances or inner products or other notions of similarity between incomplete vectors, as other methods require. To this end we assign each (incopplteeddata) point its own (full-data) subspace, and simultaneously minimize over the Grassmann manifold: (a) the chordal distance between each point and its assigned subspace, to guarantee that the subspace stays near the observed entries, and (b) the geodesies between subspaces of all data, to encourage the subspaces from points that belong together to fuse (i.e, represent the same space). At the end of this minimization, clustering the proxy subspaces using standard procedures like k-means or spectral clustering [9, 23, 56, 61, 67, 71, 73] can be done as proxy for clustering the incomplete-data (goal ii). The ability to cluster the subspaces rather than the incomplete-data is the key strength we gain by moving to the Grassmannian After clustering, the missing entries can be filled (goal i) using low-rank matrix completion. Once the data is clustered and completed, the underlying subspaces can be trivially inferred (goal dii) with singular value decomposition. Local convergence guarantees follow casily from known manifold optimization results. We complement our theoretical results with experiments on both synthetic and real data that show the potential of the foundational fusion-over -eee--vvvreeheeeee wwweerrte-Grrssmaae Grassmann formulation.
2 RELATED WORK
Number of Subepaczes
K
Fall-Data
SC
Due to its broad applicability, HRMC has attracted consid- erable attention in recent years. Existing approaches can be divided in three main groups: generalizations from low-rank matrir completion (LRMC), generalizations from subspace clustering (SC), and methods specifically tailored for HRMC (see [38] for a recent survey).
Missing Data
LRMC
HRMC (This paper)
Figure 1: HRMC a generaliza- tion of principal component anal- ysis (PCA), LRMC, and SC.
HRMC LRMC, LRMC seeks to exactly recover the missing entries of data matrix X whose columns lie in a single low-dimensional subspace [11]. One can view HRMC as a generalization of LRMC, where the columns of X are known to lie in a union of subspaces (UoS), each of low dimension, but it is not known to which subspace each column belongs (see Figure 1). Research in LRMC over the last decades has resulted in theory and algorithms that guarantec perfect recovery under reasonable assumptions (c.g., random sampling and bounded-coherence of the data) [10, 11, 15, 16, 29, 53). Hence, given a HRMC problem, if the number of underlying subspaces, say K, and the maximum of their dimension, say r, are low, one could be tempted to cast HRMC as a LRMC problem. In such case, the single subspace containing all the columns of X would have dimension no larger than r' ;=r K. This would, however, completely ignore the union structure present in the data, and therefore require more observed entries in order to complete X. We can see this by noting that each column must have more observed entries than the subspace containing it [51]. This means that even in the fortunate case where r' is low enough, using LRMC would require K times more observations than HRMC. This is especially prohibitive in applications such as Metagenomics or Drug Discovery, where data is extremely sparse and costly to acquire. In general, may be too large to even allow the use of LRMC.
HRMC vs SC. SC aims to cluster the columns of full-data matrix X according to a UoS that is not known a priori [22]. One can thus view HRMC as the generalization of SC to the case where data is missing (see Figure 1). There exists a vast repertoire of theory and algorithms that guarantee perfect clustering under reasonable assumptions (e.g., sufficient sampling and subspace separation) [68, 43, 65, 2, 59, 19). Hence, a natural approach to HRMC is thus to fill missing entries naively (with zeros, means, or LRMC) prior to clustering with a full-data method, like sparse subspace clustering [22, 40, 75]. Unfortunately, this approach may work if data is missing at rate inversely proportional to the dimension of the subspaces [64], but fails with moderate volumes of missing data, as data filled naively no longer lies in union of subspaces [21].
2

---

Under review as a conference paper at ICLR 2023
Tailored HRMC algorithms. Algorithms specifically designed to solve the HRMC prob- lem can be further divided in the following subgroups: (1) neighborhood methods that cluster points according to their overlapping coordinates [24], (2) alternating methods, like EM [50], k-subspaces [7], group-lasso [52, 64], SR [39], or MCOS [41] (3) liftings, which exploit the second-order algebraic structure of unions of subspaces [68, 70, 48, 26, 27/. and (4) integer programmiing [57]. Neighborhood methods require either abundant observations or à super-polynomial number of samples (to produce enough overlaps). Liftings require squaring the dimension an already high aigh-dimensional problem, which severely limits their applicability. Integer programming approaches are similarly restricted to small data.
To summarize, while much research has been devoted HRMC, current algorithms have shortcomings, and little is known regarding their theoretical guarantees.
Our work in context. Among the methods discussed above, the approach of this paper is perhaps closer in principle to [47]. which uses a similar Grassmannian optimization model to study the single-subspace problem of LRMC. This paper generalizes these ideas to the much harder multiple-subspaae problem of HRMC, while maintaining the local convergence guarantees of Proposition 5.1 [47]. The main difference between [47] and our formulation is that the former only considers predefined subset of geodesic distances (see equations (17)-(19) in [47]), which determine the Grassmannian points that must be matched. In [47], these subsets of geodesics
be chosen somewhat arbitrarily, because in LRMC all points belong to the same subspace. so-called gossip protocol is therefore suitable in the easier problem of LRMC. In contrast, HRMC requires that only certain subsets of the Grassmannian be matched (the points corresponding to the unknown clustering to be learnt). Without knowing priori the correct clusters, one cannol. utilize the gossip method and must therefore use all pairwise geodesics so as to not introduce bias.
Note: Appendices A-E contain review on the mathematical background involved in our formulation.
3 MODEL AND MAIN RESULTS
Let X XI....XX, E Im lie near a union of subspaces with dimension upper bounded by r. Let xfe RiBil denote the observed entries of x;, indexed by Ω (l...., m). We propose assigning to each observed vector xf proxy subspace U; span(U,). Our goal is to estimate the true subspace 24" to which x; belongs by (a) enforcing that the proxy space Ui contains a possible completion of x and (b) minimizing the distance between individual proxy spaces Uú and U₁ to build consensus. This is done via the following optimization problem, where the first term achieves goal (a) and the second term achieves goal (b):
n
A n
min U.....U, U.E.U.
Σ (x
U,)
Σ &(U,U U,). 9
(1)
where
do(xf. U,) 1
and d,(U,,
arccos2 0r (r(U.U
Here S(m,r) r) denotes the Stiefel manifold of m x orthonormal matrices, À 2 0 is regular- ization parameter, oz(-)) denotes the gth largest singular value, and X, is the orthonormal matrix spanning all the possible completions of a non-zero xf. The space of all possible completions of xf2 is therefore X;" := span(X.)), which clearly contains the true data x₁. The matrix x can be casily constructed as follows. If x = o, then x I, the identity matrix. Otherwise, X is the m x (m - [2,1 + 1) matrix formed with xf normalized and filled with zeros in the unobserved rows, concatenated with the (m - ((:)) canonical vectors indicating
3

---

Under review as a conference paper at ICLR 2023
du(U,, Up
24
d. Up
24
A.0
44
The semi-spheres represent the Grassmannian G(m, where each point 24 represents subspace (in the particular case of G(3, 1), the line going from the origin to (4). Left: Intuitively, the chordal distance (x(x., U1) is informal measure of distance between the subspace 24 and an incomplete point x² The left image should only be taken as intuition since A;" may not live on the same Grassmannian and the chordal distance should not be thought of as a geodesic distance. Center: The geodesic distance (U, Uj) measures the distance over the Grassmannian between Ui and Uj. Right: The Euclidean gradient vector V. falls out of the Grassmann manifold; account for the Grassmannian curvature, each geodesic step needs to be adjusted according to (4).
Figure 2:
the unobserved rows of xf. For example, if xff + O is observed in the first [22] rows, then
} |
X
0
m [0 +1
When xf fully observed, x simplifies to x normalized. Recall that Grassmannians G(m, r) are quotient spaces of Stiefel manifolds 8(m, r) by action of the orthogonal group of r x orthonormal matrices. Since both terms (xf",, U:) and dg (U; Uj) are invariant under this quotient, the objective function in (1) does not depend on the choice of basis, and descends to function on the Grassmannian.
Why should this work? The chordal distance (xft, U,), as defined in [72], is not a formal distance on the Grassmannian, but rather measures how far U4, is from containing a possible completion of xf. More precisely, dx(x?" Uj) is the cosine of the angle between the nearest completion of xf² and the r-plane 24, If the top singular value o (X) U;) is then X, and U4 intersect at least a line, meaning that the proxy space Ui contains a possible completion of x While merely forcing Uj to contain a possible completion offers no way to distinguish one possible completion to another, consensus among data is built different proxies Ui and 24 are forced towards one another by the geodesic term. In other words, U4, and U are allowed to be near each other, and hence form clusters, only if they both contain possible completions of both x and xfe The term chordal distance used in this way is adopted from [72) and should not be confused with the more common chordal distance between points on the Grassmannian [18]. See Figure? 2 to build some intuition.
Solving (1) The gradients of (1) with respect to U; over the Grassmannian are given by:
vdl(ff U,)
-X U: U1,
(2)
v dU U,) 2 Carccoss:(
(I (I- u.
(3)
V 1 (
where Vj and w, the leading left. and right singular vectors of X U, and v/. w/ are
the gth left and right singular vectors of U, U,U,TU, U, The key behind these expressions is that tangent vectors on the Grassmannian can be computed as projections of gradient vectors in Euclidean space [1, 20). In fact, that is exactly what the gradient expressions in (2) and (3) are: -20, (X!" !! in (2) is the gradient of de (xff U;) with respect to the entries of U; in Fuclidean space. The multiplication by I I- UU takes the horizontal direction of the
4

---

Under review as a conference paper at ICLR 2023
tangent vector with respect the quotient, thus mapping the gradient from Euclidean space to the Grassmannian. The same is true for the term I - U in (3). To summarize, (2) and (3) are gradient directions in the manifold of subspaces, rather than matrices: (2) is the steepest direction along which the subspace Ui be descended to potential subspace that contains xf" and (3) is the steepest direction along which the subspacc 24 can be descended to the subspace Uz The derivations of these gradients in Appendix E. Putting together the chordal and geodesic gradients, our overall descent direction for U; is given by
n
V = V U1)
V d (U, Uj).
Observe that U: V, falls out of the Grassmannian for every step size j (see Figure 2). To adjust for the curvature of the manifold, the update after taking a geodesic step of size over the Grassmannian in the direction of -V; is given by equation (2.65) in [20], which in our context reduces to:
diag cos (1 (n)r) i) diag sin (1/ Mi)
U + U.E.
r
E
where r E is the compact singular value decomposition V, In our implementation, we use Armijo step sizes, given by 7 = 8 where пo > 0, and 8, 7 E (0, 1) are the Armijo tuning parameters related to the initial step size and step scarch granularity [1], and v is the smallest non-negative integer such that
Σ ((U.) ((Ru, no 2 -7 Σ Vi.
where fi(U,) is the component of the objective function (1) holding i fixed and ranging j. and Ru, (A) performs the geodesic step described by (4) in the direction of A.
Convergence guarantees. One advantage of approach is that techniques in [47] to obtain local convergence guarantees like the following:
use standard
Un)} be a sequence of iterates generated by the geodesic steps given by equation (4) with Armijo steps sizes 1) as defined above. Then the sequence will converge to critical point of (1).
Proposition 1.
Let
Proof. It suffices to show the (rather technical) fact that the gradient steps in (4) instance of Accelerated Line Search (ALS) given in [t] and outlined in Algorithm 1, where the product manifold n will serve the Riemmannian manifold M, the tangent space of which is the cartesian product of tangent spaces of each constituent G. To see this, let T T' denote the tangent bundle of (set of all tangent vectors to) M, and let T&M denote the tangent space of M at N E .M. In our case, U is the tuple (UI...., Un). and equation (4) serves as the retraction! Ru, on each component, so that Ry = (Ru, is indeed retraction by recognizing (4) as the exponential map Exp TuG - G and noting that, on Riemannian manifold, the exponential map is a retraction, and that the product of exponential maps is again exponential map [1]. For our sequence of gradient--elaten tangent vectors, we use the negative gradient, which is clearly always gradient-related. The
Ru,).. One can verify that this
' mapping R from TM to M such that its restriction to T,M, denoted Ru, satisfies local rigidity condition which preserves gradients U; see the rightmost illustration in Figure build some intuition, or Chapters 3 and
[1] for a more careful treatment of these definitions. 2Given cost function f on Riemannian manifold M, sequence of tangent vectors {At) A: E M, is reeteed if, for any sequence (4ee that converges to non-critical point of f. the corresponding subsequence (At)ek is bounded and satisfies lim sup, 45,, lek (V ((4)). A.) <
5

---

Under review as a conference paper at ICLR 2023
Algorithm 1: Accelerated Line Search (ALS)
Require: Riemannian manifold M; continuously differentiable scalar field f on M; retraction R from T' to M; scalars do 0, c, 6,8,75 (0, 1).
Input: Initial iterate U E MM
Output: Sequence of iterates {U.).
for t 0. 1, do
Pick A, E M such that the sequence of tangent vectors {A.} } is gradient related.
Select Ui-1 such that
S(4) - f(4+1) 2 c(((4) - [(Ru, (mA.)))
(5)
where th is the Armijo step size for the given Yo, 3, 7, A..
end
A
Small A
Ideal A
Large
Figure 3: À 2 0 in (1) regulates how clusters fuse together. Tf À = O, each point is assigned to subspace that exactly contains it (overfitting.. The larger X, the more penalize subspaces being apart, which results in subspaces getting closer to form fewer clusters. The extreme case À the special case of PCA and LRMC, where only one subspace is allowed to explain all data.
gradient on the product manifold is the cartesian product of the gradients on each constituent manifold, i.e. V(f) (V(f)..... V (fa)). Moreover, the inner product on the tangent space is the sum of the inner products on the constituent tangent spaces. Therefore, if jA,
A E Tu, M, gradicnt-related for each Mi, then {(A Ana)] is gradient-related on the product manifold. Furthermore, setting Uu+1 = Rı,, (m A:) satisfies the bound in (5) with e = Thus Proposition 1 follows as consequence Theorem 4.3.1 and Corollary 4.3.2 in [1].
FUSION IN PRACTICE
Clustering, completion, and subspace inference. Recall that the solution to (1) pro- vides an estimator Ut of U1, the true subspace from which is drawn. After solving (1), one may form the matrix D whose (t entry is given by d, (U, Uj) and use it as input to any distance-batt based clustering method, such as k-means [9, 23, 56]. spectral clustering [61, 67, 71, 73], or DBSCAN [25, 32, 55]. While prior knowledge of the number of subspaces K may be required for some clustering methods (e.g., k-means, or spectral clustering), at all to solve (1). Hence, by choosing clustering method that doesn't require knowing K (e.g., DBSCAN), our approach can be applied to situations where K is unknown. After clustering, one agglomerate all the data points corresponding to the tth cluster in the same matrix X and run any low-rank malrix completion (LRMC) algorithu (e-g., [6, 8, 10, 11, 12, 36, 53, 69]) estimate its completion X₁. Finally, one can run principal component analysis (PCA) |46, 74] on Xg to recover an estimate basis Üg of the t underlying subspace U:
not required
Penalty parameter and model selection. Intuitively, the chordal term (1) forces each subspace to be close to its assigned data point, and the geodesic term forces subspaces from different data points to be close to one another. The tradeoff between these two quantities is determined by the penalty parameter X > 0. If À 0, then the geodesic term is ignored and there a trivial solution where each subspace exactly contains its assigned data point (thus attaining the minimum, zero, for the chordal distance). If X > 0, the geodesic term forces subspaces from different data points to get closer, even they no longer contain exactly their assigned data points. As À grows, subspaces get closer and closer (see Figure 3). The

---

Under review as a conference paper at ICLR 2023
extreme case (A oo) forces all subspaces to fuse into one (to attain zero in the second term), allowing only one subspace to explain all data, which is the equivalent of PCA the complete-data case, and LRMC if data is missing. In other words, PCA and LRMC are the special cases of our formulation with X = xc.
Ultimately, the effect of À will be reflected in the distance matrix D which in turn determines the number of clusters. The smaller X, the more clusters, up to the extreme where each data point is in its own cluster. Conversely, the larger X, the fewer clusters, up to the extreme point where all points are clustered together. The more subspaces, the more accuracy, but the more degrees of freedom (overfitting). To determine the best X, one can compute a goodness of fit test, like the minimum effective dimension [33, 68), that quantifies the tradeoff between accuracy and degrees of freedom. Similarly, we can iteratively increase r in (1) to find all the data points that lie in 1-dimensional subspaces, then all the data points that lie in 2-dimensional subspaces, and so on (pruning the data at each iteration). This will result in an estimate of the number of subspaces K, and their dimensions.
Initialization. In our implementation we initialize (1) with a solution to the problem when A= 0, Le., when each subspace perfectly contains the observed entries of its assigned data point, To this end, for each we first construct an m x r matrix whose first column is equal to x² in its observed entries, and whose remaining entries are filled with standard normal entries, known to produce incoherent and uniformly distributed subspaces [24]. This matrix is then orthonormalized to produce the initial estimate Uj, which, by construction, contains xf, thus producing do(x;" U;) =0 O.
Computational complexity. We point out that the main caveat of our approach is its quadralic complexity in the number of samples. Fortunately, subspace clustering allows a simple approach to sketching both samples and features [62]. That is, one may solve (1) with a subset of n' < n columns, and a subset of m' < m rows (e.g., those with most observations), resulting in an improved complexity, quadratic in n' as opposed to n. With the solution of (1), one can use a clustering method, a LRMC algorithm, and PCA, as described above, to produce subspace estimates OI...., ÖK, with K' < K. Each of the remaining n - n' incomplete data points x,² that were not used to solve (1) and that have more than r observations (a fundamental requirement of subspace clustering [51]) can be trivially assigned (
to the subspace estimate producing the largest projection coefficient o where 0 E RR denotes the restriction of Üg to the observed rows of x8 (notice that GgT 0 is invertible for almost every rank-r UR whenever |Sa]] > [51]). If xf is assigned to UR its completion can be trivially estimated as X, Ü2 All the data points xf that are too far from all of the subspace estimates (equivalently, the data points whose coefficients are smaller than pre-determined parameter) can be used to solve (1) again for refined clustering.
EXPERIMENTS
In this section we present series of experiments on real and synthetic data, in particular the Hopkins 155 dataset [63], and the Smartphone dataset for Human Activity Recognition in Ambient Assisted Living (AAL) [4]. Rather than establishing new state-of-the-art, these experiments have the intention to serve as proof of concept, showing the potential of our approach, which in this first introduction and basic formulation performs comparable to prominent methods [75]. In our experiments we initialize (1) as described in Section 4, with r fixed, as known a priori in both, the simulations, and the real datasets. We do not specify K, and we make special adjustments to handle noise, as is nol required by our approach. The attained solution to (1) is used as input to spectral clustering [61, 67, 71, 73) (though, as described in Section 4, other clustering algorithms, such as k-means [9, 23, 56] or DBSCAN [25, 32, 55) could be used). We measure accuracy in terms of clustering error, given by mins n Σ (()(()))y; where I denotes the indicator function, and M is a function that maps the estimated cluster labels & E (I,...,K) assigned to x xf....,xff to the true labels y E (1,...,K)'
7

---

Under review as a conference paper at ICLR 2023
m
100, n 100, r 5
n
100, - 5, K
m 100, n
100, K
0.6
papert
Clustering Error
0.5
2F-55C 2=
2F+S5G K
0.5
25153C 150
SSC EWZF,
0.4
1th papert 2
jihis papert K
0.4
SSC BWZF: 150
#1+555: m 203
0.4
0.3
0.3
SSC EWZF: 200
0.3
GF (thi papar): m 150
0.2
0.2
OF (thi paport: m 200
0.2
0.1
0.1
0.1
0.0
0.0
0.0
0.2
0.4
0.6
0.8
1.0
0.2
0.4
0.6
0.8
1.0
0.2
0.4
0.6
0.8
1.0
Sampling Rate (p)
Sampling Rate (p)
Sampling Rate (p)
Figure Clustering error (average over 10 trials) as a function of sampling rate for different synthetic settings.
Baseline comparisons and full-data sub-optimality. A recent survey [38] shows that most state-of-the-ar algorithms for HRMC (including MC+ [75], EM [50], GSSC [52], k- subspaces [7], and more) have similar performance, with varying winners on specific scenarios ambient dimension gap, fraction of missing data, and number of subspaces. Based on this recent survey [38], and others [52], we chose ZF-SSC and SSC-EWZF as baselines, which has been seen in [38] to have nearly identical performance as MC+SSC, EM, GSSC, k-subspaces and k- k-GROUSE [7] in the scenarios discussed in our paper.
depending on subspace
Synthetic data, In all our simulations we first generate K matrices U E m with ii.d. ((00 entries, to use as bases of the true subspaces. For each k we generate a matrix Ø E Rexi, also with i.i.d. N(0, 1) entries, to use as coefficients of the columns in the kth subspace. We then form X as the concatenation [U10], U205. missing data we sample each entry independently with probability p. Figure shows the clustering results as a function of the sampling rate for a variety of settings, tuning the parameter manually. Notice that, even with this first formulation, we perform comparable to existing methods, and even better in some cases, especially in low-sampling regimes.
0k To induce
Object tracking in Hopkins 155. This dataset contains 155 videos of K or K moving objects, such as checkerboards, vehicles, and pedestrians. In each video, collection of n mark points are tracked through all frames. The locations over time of the ith point are stacked to produce x; E Rm, so that the points corresponding to the same object lie near a low-dimensional subspace [60, 35] (r varies from video to video, from to 3). In all cases we fixed the penalty parameter À to 1. To induce missing data (e.g., produced by occlusions) we sample each entry independently with probability p. Figure shows the clustering results.
Human activity recognition in Smartphone AAL Dataset. This dataset contains n = 5744 instances, each with m = 561 features related to pre- pre-processed accelerome- ter and gyroscope time series and summary statistics (3], related to K = 2 activities: walking, and other movements, each approximated by a subspace of dimension r = 4. Recall that the complexity of (1) is quadratic n, so if solved directly, this dataset that would produce an unmanageable computational complexity. However, using the sketching techniques described in Section it can be solved quite efficiently. In particular, we only used m' 158 features (related to the accelerometer's and gyroscope's minimum, maximum, standard deviation, and mean parameters over time), and n' 100 samples selected at random, evenly distributed among classes. In all cases we fixed the penalty parameter X to 10-5 The results are summarized in Figure 5. Notice that. our approach outperforms existing methods in the low-sampling regime.
Lastly, we note the disparity in performance between our model and existing algorithms when the missing data rate is low. The main motivation for our approach is incomplete data. While our formulation can certainly be used with full data, we acknowledge that it would be over-kill, and consequently suboptimal in that scenario, which has been extensively studied. Hence, it is not surprising that methods tailored for full-data outperform ours in such setting. However, no full-data method outperforms ours when data is missing.
8

---

Under review as a conference paper at ICLR 2023
Hopkins155
SmartphoneALL
0.30
Clustering Error
0.4
(thi paperi
0.20
2F-55C
0.10
0.1
sapert
papers K
0.5
0.6
0.8
0.9
1.0
0.2
0.4
0.6
0.8
1.0
Sampling Rate (p)
Sampling Rate (p)
Figure 5: Clustering error function of sampling rate for real datasets. Left: average over 120 videos with 2 objects, and 35 videos with K objects. Right: average over 20 trials.
FUTURE DIRECTIONS AND CHALLENGES
The main formulation presented in this paper is a non-convex optimization that relics on the simultaneous interactions of many terms. A proper analysis of the model is therefore challenging. One difficulty confronted is the complex geometry of the zero-set of the chordal distance lerm, more precisely, the intersections of many of these zero-sets, one for each column data. While intuition suggests that the model is encouraged by the geodesic term to find regions of "dense" intersection, and therefore build consensus, more precise formulation this intuition has evaded us. This is further highlighted in Figures and by the fact. that the performance of our model decreases as K grows, indicating that it is not currently understood how the combination of the chordal and geodesic terms encourage consensus amongst so many cross-cluster terms. It is our belief that well-designed weighted version of (1), such as
min U1 U.E...
D Σ Lxxx
A 2
Σ (U, 1J=1
(U, (J)
where the weights Wij 2 O quantify how much attention is given to cach penalty, is key to unlocking better performance and understanding of the model. Our immediate future work will focus on investigating options for these weights, such as inverse distance functions, or k-nearest neighbors, known to dramatically improve the performance, computational complexity, and tolerance to K of fusion formulations in Euclidean space [5, 17, 58, 76).
7 CONCLUSIONS
This paper presents new paradigm for clustering incomplete datapoints using subspaces as praxies leverage the geometry of the Grassmannian. This new perspective enables clustering and completion of data in a union of subspaces. This work should be understood the idea of fusion penalties in the Grassmann manifold, for
as the first introduction the problem of high-rank matrix completion. Rather than establishing our approach as the state-of-the-art our experiments have the intention to serve as proof concept, showing that there is potential in our approach, in the hopes to ignite future work on several directions, such as the study of weighted versions described in (6), the choice of penalty parameters, and variants robust to outliers.
9

---

Under review as a conference paper at ICLR 2023
REFERENCES
[1] P-A Absil, Robert Mahony, and Rodolphe Sepulchre. Optimization algorithms on matria manifolds. Princeton University Press, 2009.
[2] Pankaj K Agarwal and Nabil H Mustafa. "K-means projective clustering". In: Proceed- ings of the twenty- third A ACM of databose systems. 2004, pp. 155-165.
symposium on Principles
[3] Davide Anguita et al. "A public domain dataset for human activity recognition using smartphones." In: Esann. Vol. 3. 2013, p. 3.
[4] Davide Anguita al. UCI Machine Learning Repository. 2013. URL: http://archive. ics. uci. edu/ml.
[5] Anestis Antoniadis and Jianqing Fan. 'Regularization of wavelet approximations", In: Journal of the American Statistical Association 96.455 (2001), pp. 939-967.
[6] Laura Balzano, Robert Nowak, and Benjamin Recht. "Online identification and tracking of subspaces from highly incomplete information". In: 2010 48th Annual allerton conference on communication, control, and computing (Allerton). IEEE. 2010, pp. 704- 711.
[7] Laura Balzano et al. "K-subspaces with missing data". In: 2012 IEEE Statistical Signal Processing Workshop (SSP). IEEE. 2012, pp. 612-615.
[8] Jonathan Bauch, Boaz Nadler, and Pini Zilber. "Rank iterative least squares: efficient recovery of ill-conditioned low rank matrices from few entries". In: SIAM Journal on Mathematics of Data Science 3.1 (2021), pp. 439-465,
[9] Leon Bottou and Yoshua Bengio. "Convergence propertics of the k-means algorithms". In: Advances in neural information processing systems. 1995, pp. 585-592.
[10] Emmanuel J. Candes and Terence Tao. "The Power of Convex Relaxation: Near- Optimal Matrix Completion". In: JEEE Trunsactions on Information Theory 56.5 (2010), pp. 2053-2080. DOI: 10. 1109/TIT. 2010. 2044061.
[11] Emmanuel Candès and Benjamin Recht. "Exact matrix completion via convex optimization". In: Foundations of Computational mathemotics 9.6 (2009), pp. 117-772.
[12] Emmanuel Candès et al. "Robust. principal component analysis?" In: Journal of the ACM (JACM) 58.3 (2011), pp. 1-37.
[13] Jiuwen Cao et al. "Extreme learning machine and adaptive sparse representation for image classification". In: Neural networks 81 (2016), pp. 91 102.
[14] Guangliang Chen and Gilad Lerman. "Spectral curvature clustering (SCC)". In: Inter- national Journal of Computer Vision 81.3 (2009), pp. 317-330
[15] Y. Chen. "Incoherence-optimal matrix completion". In: IEEE Transactions on Infor- mation Theory 61.5 (2013), pp. 2909 2923.
[16] Y. Chen al. "Coherent matrix completion". In: International Conference on Machine Learning 61.5 (2014), pp. 674-682.
[17] Eric C Chi and Kenneth Lange. "Splitting methods for convex clustering" In: Journal of Compulational and Graphical Statistics 24.4 (2015), pp. 994-1013.
[18] John H. Conway, Ronald H. Hardin, and Ncil J. A. Sloane. "Packing lines, plancs, etc., packing in Grassmannian spaces". In: Exper. Math. 5.2 (1996), pp. 139-159.
[19] Harm Derksen et al. "Segmentation of multivariate mixed data via lossy coding and compression" Tn: Visual Communications and Image Processing 2007. Vol. 6508. SPIE. 2007, pp. 170-181.
[20] Alan Edelman, Tomás A Arias, and Steven T Smith, "The geometry of algorithms with hogoognny constraints". In: SIAM journal on Malrir Analysis and Applications 20.2 (1998), pp. 303-353.
[21] Ehsan Elhamifar. "High-rank matrix completion and clustering under self-cxpressive models". In: Advances in Neural Information Processing Systems 29 (2016), pp. 73-81. [22] Ehsan Elhamifar and René Vidal. "Sparse subspace clustering: Algorithm, theory, and applications" In: IEFF transactions on pattern analysis and mochine intelligence 35.11 (2013), pp. 2781-
10

---

Under review as a conference paper at ICLR 2023
[23] Charles Elkan. "Using the triangle inequality to accelerate k-means". In: Proceedings of the 20th international conference on Machine Learning (ICML-03). 2003, pp. 147-153. Brian Eriksson, Laura Balzano, and Robert Nowak. "High-rank matrix completion". In: Artificial Intelligence and Statistics. PMLR. 2012, pp. 373 381.
[24]
[25]
Martin Ester et al. "A density-based algorithm for discovering clusters in large spatial databases with noise." " In: Kdd. Vol. 96. 34. 1996, pp. 226-231.
[26] Jicong Fan and Tommy WS Chow. "Non-lincar matrix completion". In: Pattern Recog- nition 77 (2018), pp. 378 394.
[27] Jicong Fan and Madeleine Udell. "Online high rank matrix completion". In: Proceedings of the IEEE/CC /CVF Conference on Computer Vision and Pattern Recognition. 2019, pp. 8690-8698,
[28] Sonny Gan. "The Application of Spectral Clustering in Drug Discovery". PhD thesis. University of Sheffield, 2013.
[29] D. Gross. "Recovering low-rank matrices from few coefficients in any basis", In: IEEE Transactions on Information Theory 57.3 (2011), pp. 1548-1566.
[30] Jeffrey Ho et al. "Clustering appearances of objects under varying illumination condi- lions". In: 2003 IEEE Compuler Sociely Conference on Computer Vision and Pallern Recognition, 2003. Proceedings. Vol. 1. IEEE. 2003, pp. I-I.
[31] Wei Hong et al. "Multiscale hybrid linear models for lossy image representation". In: IEEE Transactions on Image Processing 15.12 (2006), pp. 3655-3671.
[32] Jian Hou, Huijun Gao, and Xuelong Li. DDSANs-DBSCAN: A Parameter- Free Clustering Algorithm". In: IEEE Transactions on Image Processing 25.7 (2016), pp. 3182 3193. DoI: 10. 10.11091 1109/TIP. 2016. 2559903.
[33] Kun Huang, Yi Ma, and René Vidal. "Minimum effective dimension for mixtures of subspaces: A robust. GPCA algorithum and its applications" In: Proceedinggs of the 2004 IEEE Compuler Suciely Conference on Compuler Vision and Pollern Recognilion, 2004. CVPR 2004. Vol. 2. IEEE. 2004, pp. II-II.
[34] Lee M. John. Introduction to Riemannian Manifolds, 2nd. Springer International Publishing, 2018
[35] Ken-ichi Kanatani. "Motion segmentation by subspace separation and model selection". In: Proceedings Eighth IEEE International Conference on computer Vision. 1CCV 2001. Vol. 2. IEEE. 2001, pp. 586 591.
[36] Raghunandan H. Keshavan, Andrea Montanari, and Sewoong Oh. "Matrix Completion From Few Entries". In: IEEE Transuclions on Information Theory 56.6 (2010), pp. 2980 2998. DOI: 10. 1111/TTTTTTT.. 2010. 2046205.
[37] Hamidreza Koohi and Kourosh Kiani. "A new method to find neighbor users that im- proves the performance of collaborative filtering". In: Expert Systems with Applications 83 (2017), pp. 30-39.
[38] Connor Lane et al. "Classifying and comparing approaches to subspace clustering with missing data". In: 2019 JEEE/CT /CVF International Conference on Computer Vision Workshop (ICCVW). IEEE. 2019.
[39] Chun-Guang
Guang Li and René Vidal. "A structured sparse plus structured low-rank frame- work for subspace clustering and completion". In: IEEE Transactions on Signal Pro- cessing 64.24 (2016), pp. 6557-6570.
[40] Chun-Guang Li, Chong You, and René Vidal. "Structured sparse subspace clustering: A joint affinity lcarning and subspace clustering framework". In: IEEE Transactions Image Processing 26.6 (2017), pp. 2988-3001
[41] Ziheng Li et al. "Matrix completion with column outliers and sparse noise". In: Infor- mation Sciences 573 (2021), pp. 125-140,
[42] Le Lu and René Vidal. "Combined central and subspace clustering for computer vision applications". In: (2006), pp. 593-600.
[43] Yi Ma et al. "Estimation of subspace arrangements with applications in modeling and segmenting mixed data". In: SLAM review 50.3 (2008), pp. 413 458.
[44] Mohamed G Malhat, Hamdy M Mousa, and Ashraf B El-Sisi. "Clustering of chemical data sets for drug discovery". In: 2014 9th international conference on informatics and systems. IEFE. 2014, DEKM-11.
11

---

Under review as a conference paper at ICLR 2023
[45] Korosh Mashayekh et al. "Clustering and Sampling of the c-Met Conformational Space: A Computational Drug Discovery Study". In: Combinatorial chemistry high throughput screening 22.9 (2019), pp. 635-648.
[46] Thomas Minka. " "Automatic choice of dimensionality for PCA" In: Advances in neural information processing systems 13 (2000), pp. 598-604.
[47] B. Mishra et al. "A Riemannian gossip approach to subspace learning on Grassmann manifold". In: Machine Learning 108.10 (2019), pp. 1783-1803)
[48] Greg Ongie et al. "Algebraic variety models for high-rank matrix completion". In: International Conference on Machine Learning. PMLR. 2017, pp. 2691-2700.
[49] Nolen Joy Perualila- Tan et al. "Weighted similarity- basied clustering of chemical struc- tures and bioactivity data in early drug discovery". In: Journal of bioinformatics and computational biology 14.04 (2016), p. 1650018.
[50] Daniel Pimentel, R Nowak, and Laura Balzano. "On the sample complexity of subspace clustering with missing data". In: 2014 IEEE Workshop on Stotistical Signal Processing (SSP). IEEE 2014, pp.
[51] Daniel Pimentel-Alarcon and Robert Nowak. "The informattion-theorllie requirements of subspace clustering with missing data". In: International Conference on Machine Learning. PMLR. 2016, pp. 802 -10.
[52] Daniel Pimentel-Alarcee et al. "Group-sparse subspace clustering with missing data". In: 2016 JEEE Statistical Signal Processing Workshop (SSP). IEEE. 2016, pp. 1-5.
[53] Benjamin Recht. "A Simpler Approach to Matrix Completion." In: Journal of Machine Learning Research 12.12 (2011).
[54]
Budhaditya Saha et al. "Sparse subspace clustering via group sparse coding". In: Proceedings of the 2013 SIAM International Conference on Data Mining. SIAM, 2013, pp. 130-138.
[55] Erich Schubert et al. "DBSCAN revisited, revisited: why and how you should (still) use DBSCAN". In: ACM Transactions on Database Systems (TODS) 42.3 (2017), pp. 21. [56] David Sculley. "Web-scle k-means clustering". Tn: Proceedings of the 19th international conference on World wide web. 2010, pp. 1177 1178.
[57]
Akhilesh Soni el al. "Integer Programming Approaches To Subspace Clustering With
Missing Data". In: ()-
[58] Kean Ming Tan and Danicla Witten. "Statistical propertics of convex clustering". In: Electronic journal of statistics 9.2 (2015), 2324.
[59] Michael E. Tipping and Christopher M Bishop. "Probabilistic principal component
analysis". Tn: Journal of the Royal Statistionl Society: Series B (Statistical Methodology) 61.3 (1999), pp. 611-622
[60] Carlo Tomasi and Takeo Kanade. "Shape and motion from image streams under ort orthography: a faclorization method"" In: International journal of compuler vision 9.2 (1992), pp. 137-154.
[61]
Tao Tong et al. "Ono-step spectral clustering based on self-paced lcarning". In: Pattern Recognition Letters 135 (2020), pp. 8 14.
[62]
Panagiotis A Traganitis and Georgios B Giannakis. "Sketched subspace clustering". In: JEFE Transactions on Signal Processing 66.7 (2017), pp. 1663-1675.
[63] Roberto Tron and René Vidal. "A benchmark for the comparison of 3-d motion segmentation algorithmss". In: 2007 IEEE conference on compuler vision and pullern recognilãon. IEEE. 2007, pp. 1-8.
[64] Manolis Tsakiris and René Vidal. "Thcoretical analysis of sparse subspace clustering with missing entrics". In: International Conference on Machine Learning. PMLR. 2018, pp. 4975 4984.
[65]
Paul Tseng. "Nearest q-flat to m points". In: Journal of Optimization Theory and Applications 105.1 (2000), pp. 249-252.
[66] Farman Ullah, Ghulam Sarwar, and Sungchang Lee. "N-screen aware multicriteria. hybrid recommender system using weight based subspace clustering" In: The Scientific World Journal 2014 (2014).
12

---

Under review as a conference paper at ICLR 2023
[67] Patrick Veenstra, Colin Cooper, and Steve Phelps. "Spectral clustering using the kNN- MST similarity graph". In: 2016 8th Computer Science and Electronic Engineering (CEEC). IEEE. 2016, pp. 222-227.
[68] Rene Vidal, Yi Ma, and Shankar Sastry. "Generalized principal component analysis (GPCA)". In: IEEE transactions on pattern analysis and machine intelligence 27.12 (2005), pp. 1945-1959.
[69] René Vidal and Paolo Favaro. "Low rank subspace clustering (LRSCY", In: Pattern Recognition Letters 43 (2014), pp. 47-61.
[70] René Vidal, Roberto Tron, and Richard Hartley. "Multiframe motion segmentation with missing data using PowerFactorization and GPCA". In: International Journal of Computer Vision 79.1 (2008), pp. 85 105.
[71] Ulrike Von Luxburg. "A tutorial on spectral clustering". Tn: Stalistics and computing
17.4 (2007), pp.
[72] Ely Kerman Wei Dai and Olgica Milenkovie. "A Geometrie Approach to Low-Rank Matrix Completion". In: IEEE Transaction on Information Theory 58.1 (2012), pp. 237- 247.
[73] Guoqiu Wen. "Robust self-tuning spectral clustering". In: Neurocomputing 391 (2020),
pp. 243 248.
[74] Svante Wold, Kim Esbensen, and Paul Geladi. "Principal component analysis". In: Chemometrics and intelligent laboratory systems 2.1-3 (1987), pp. 37-52.
[75] Congyuan Yang, Daniel Robinson, and Rene Vidal. "Sparse subspace clustering with missing entries", In: International Conference on Machine Learning. PMLR. 2015, pp. 2463-2472.
[76] Ming Yuan and Yi Lin. "Model selection and estimation in regression with grouped variables", In: Journal of the Royal Statistical Society: Series B (Statistical Methodology)
68.1 (2006), pp. 49-67.
[77] Wen Zhang et al. "RP-LGMC: rating prediction based on local and global information with matrix clustering". In: Computers & Operations Research 129 (2021), p. 105228.
A STIEFEL AND GRASSMANN MANIFOLDS
The primary mathematical object involved in this work is the Grassmannian G(m,r) This smooth compact manifold of dimension r(m- - -r). A full expository on the Sticfel and Grassmann manifolds is given [20]. Here, we record the most basic necessary ideas needed in order to have a working understanding of the tools used in the above. To describe this, it necessary to define precursor objects; the orthogonal group O(m) and the Stiefel manifold S(m,r). The objects of interest are thus:
1. The orthogonal group O(m) consisting of m x orthogonal matrices;
2. The Stiefel manifold s(m,r) consisting of m r orthonormal matrices;
3. The Grassmann manifold G(m, r) obtained by identifying those matrices in S(m, whose columns span the same subspace (a quotient manifold).
In this setting, the Sticfel manifold defined as a quotient space of the orthogonal group. Here, two orthogonal matrices U and V are identified if their first columns are identical or, equivalently, if U = V where Q is an orthogonal (m - r) (m -r) block. Therefore S(m,r) = (mm) - r). Going further, the Grassmannian is defined a quotient space of the Stiefel manifold where two Stiefel elements are identified if their columns span the same r-dimensional subspace. Therefore G(m, = 8(m, (r))
Given the above, it is clear that we may describe elements of the Stiefel and Grassmann manifolds using concrete representatives that can be stored on a computer. A point on the x orthonormal matrix. A point the Grassmann
Stiefel manifold may be stored manifold, however, being a linear subspace, does not have a unique representative and can be stored arbitrary orthonormal matrix so long it spans the correct subspace.
13

---

Under review as a conference paper at ICLR 2023
B PRINCIPAL ANGLES AND SINGULAR VALUES
Recall the notion of principal angles between subspaces: let U E S(m,pp) p) and V E S(m, ?) be orthonormal bases for two arbitrary subspaces of IR" Assume, without loss of generality, that S < S m. The principal angles between Span(U) and span(V) are defined via the following construction. Let un E sp) and V1 E Vpan(V) be unit vectors such that lu] v1 is maximal. Inductively, let that u u; and v vj = for all S j and jug V&l maximal. The principal angles are defined as og arccos up va for all k = ..P.....P.
span (U) and Va span(V) be unit vectors such
This constructive definition is too cumbersome to use in practice. We opt for the following alternative computation via the singular value decomposition. Let uI....: up and Vi..... Va be the columns of U and V respectively. Compute the singular value decomposition U"V Udiagl..... Jj .j.... .."..." Set U' UÚ and V VV and denote their columns by uf and vj respectively. Observe that span = (U) and span(V) = span(V) and furthermore that U V = diagl..... .i.......] that is
Sai i=j i
vj
The vectors uf and vi correspond to those in the constructive definition. We therefore have that the i-th principal angle a; relates to the i-th singular value via d; coS az.
C CHORDAL DISTANCE
The chordal distance between points on the Grassmannian G(m,r). introduced and studied in [18], is defined via
sin
where the oj are the principal angles between the points described above. The authors of [72] introduce a notion distance between a partially observed vector x9 1Rm and a subspace
U Gm,, via formulation closcly related to the chordal distance, which they give the same name
Let X denoted the orthonormal matrix spanning all possible completions of x²: If = then X = the identity matrix. Otherwise, X is the m x (m - | + 1) matrix formed with x² normalized and filled with zeros in the unobserved rows, concatenated with the (m - [R]) canonical vectors indicating the unobserved rows of x², Let (x"T denote the
largest singular value of x U. Then
(x U) = sin V xx
This metric is studied [72]. Of particular importance is the following fact stated as Theorem 2 in [72]: the preimage of 0 under dx(x".)) is the closure of the preimage of 0 under fr(x", ,)) where
fr 2 U) = min J| - PaP
weR*
and Pp denotes projection onto the entries indexed by Q. That is, fF is the Frobenious norm, which is often used to search for subspaces U consistent with data. The Frobcnius norm may not be continuous, whereas the chordal distance continuous and differentiable.
D GEODESIC DISTANCE ON THE GRASSMANNIAN
The geodesic distance d, (U, Uj) is derived from the intrinsic geometry the Grassmann manifold and depends on the metric which defines the manifold structure. Let 7 (a, b) M
14

---

Under review as a conference paper at ICLR 2023
be eurve general Riemannian manifold (M.g) with metric g. Then the length of y is defined as [34]
L(r) - L." V ((t((),
The canonical metric for the Grassmann manifold coincides with the Euclidean metric inherited from (mm): Se, U,)
(U,, U₁) tr(Uū Uj [20]. To compute the geodesic distance, we therefore require knowledge the geodesic segment connecting U; and Uj with respect to the metric This is described in Lemma 1 of [72]: Let V,EV; be the singular value decomposition of the matrix U,UU; and denote the f-th singular value by de = cos ag. Set U, = U,V, and U V and note that U, Σ. Then the geodesic with respect to ge from U, to Uj given by U(t). 0 S & 1, where the path U(t) given by diag (L... cos act, act.....)
U0 Gl
sin diag (I... sin ast, ant......))
V.
where the columns of G [..., -.....] . E S(r, m) are defined as 02. U..e
if Xe #1 if A =
1
9:
Here, the subscript l denotes the f-th column of the corresponding matrix. We therefore have
diag (l... . - ot sin agt.....)
ú(t) [U, Gl
V.T.
diag (I..., a cos art,
Denote = -dag ([..., -a sin agt....) and C = diag ([
cos ast, .). Then
T
S
Gl
v.
G
c
Vi[S.]]
T.ūū
U.G
s
V.
G
G
c
=v v
V.diag((............
where we use the fact that U, G E 8(mm r), and that U G 0 [72]. Recall that tr(AB) = tr( (BBA),, hence tr(ú Ú)= = tr(diag([.... ()) az We therefore have L(U(!))
So a oet V. az. Recalling that ag cos ae, we finally have
(GU,,U))
Σ arccos2 (,)
E GRADIENTS ON THE GRASSMANNIAN
In this section, we derive the expressions in (2) and (3) that govern the fusion steps our formulation. For a function F(U) defined on the Grassmannian, the graduate of F at U is given by equation (2.70) in [20], which we record here:
VF Fu - UUT Fu.
where Fu the matrix whose entries are given by [Fu]y S0,
Chordal gradient. To obtain the gradient of the chordal distance dx(xf, U;)
presented in
(2), consider the partial derivative with respect to the (a, b)th element U,
OdL (x., U.) U.lbbb 2 U,) = -2 -x O U;)
do, (x U
0U:
0 U.lbbb
15

---

Under review as a conference paper at ICLR 2023
To obtain the partial derivative of the leading singular value di, observe that x U and X²UU; U; share singular values: if U; VEW!, then X (X(V)EW' with X s(m,r), the result is compact singular value docomposition. Recall that vi and
wi denote the leading left. and right singular vectors of x U,. Since xx U, 01 we have that
0o
:v x U.w; al.
OU: wi + x U
Ow:
x
ali
a(Uijah
Uiiijb
The first and third terms are zero, because w, is the leading right singular vector of X U,, so (X 0ıvi, which implics
vv! a
Ov
U;)
Va.
allbb
and because dv. a. Vi = 0, as seen by differentiating both sides viv, = (and similarly for the third term). To compute the second term, note that V; E span(X.) since it is column of V,V, and the space spanned by X invariant under multiplication by V. Now, (v. x XxxxTT;; V: = Vi, since Xx acts on vectors the projection onto span
Hence v, x = v,T. The second term then becomes vl a0. JU..f wi = [vila[w]]s It follows
that
001 = Mvlalwwlb. 8
(7)
From this, we have dxf projects onto the tangent space of the Grassmannian at U, as described before [1, 20).
U,) 2 - U,U Joi,,w where multiplication by I U.
Geodesic gradient. For the gradient of the geodesic distance d'(U;, Uj) in (3) let use de as shorthand for (²U Uj), and recall that ví and wf denote the eth left and right singular vectors UJU, U. Then the partial derivative with respect to the (a, b)th element of U,
(U, Uj)
-2 arccos OF do V1- aii
-2 arccos of
Wij
aU:
E= V1-
where the first equality follows because U; = second equality follows by parallel arguments as the derivation of (7) for the leading singular value. The last equation is the Euclidean gradient. Projecting onto the tangent space al U;, as described before [1, 20), we obtain the following gradient on the Grassmannian
= U))
(U,U))
U:), and the
2 arccos of (I U.U.)
Vd(UU
V1 of
16


Reviews:
[{"review_id": "lJq0TfTshs", "reviewer": "Reviewer_MfUS", "summary": "", "questions": "", "limitations": "", "rating": 8, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "The paper tackles the high rank matrix completion problem in which only partial observations are present and the underlying data comes from the union of subspaces.  \n\nThe fundamental difficulty with all existing approaches is that they rely on assessing distances (e.g., euclidean, or in the form of inner products) between partially observed vectors To circumvent this problem, this paper introduce a new paradigm to cluster incomplete vectors, using subspaces as proxies, thus avoiding the need to calculate distances or inner products or other notions of similarity between incomplete vectors, as other methods require. In this approach, each (incomplete-data) point is assigned its own (full data) subspace, and the objective is to simultaneously minimize over the Grassmann manifold two terms: (1) the chordal distance between each point and its assigned subspace, to guarantee that the subspace stays near the observed entries, and (2) the geodesics between subspaces of all data, to encourage the subspaces from points that belong together to fuse (i.e, represent the same space). These two terms (1) and (2) constitute the 'fusion' between two different types of distances in the non-convex optimizaiton problem.\n\nIn a nutshell, this work introduces the idea of fusion penalties in the Grassmann manifold, for the problem of high-rank matrix completion.", "summary_of_the_review": "In summary, I find the paper interesting and novel, hence  worthy of a conference publication.", "strength_and_weaknesses": "The paper is quite interesting and well written. The idea of fusion (specifically, in this example combining the chordal and geodesic terms in the optimization function) is interesting.\n\n\n", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"correctness": "4: All of the claims and statements are well-supported and correct.", "technical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work.", "empirical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work."}, "extra_text": {}}, {"review_id": "P-53WCzFC9", "reviewer": "Reviewer_vG2m", "summary": "", "questions": "", "limitations": "", "rating": 5, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "This paper introduces a new paradigm to cluster incomplete vectors using subspaces as proxies to exploit the geometry of the Grassmannian. The authors leverage this new perspective to develop an algorithm to cluster and complete data in a union of subspaces via a fusion penalty formulation. The analysis with synthetic and real-data experiments is conducted to serve as proof of concept.", "summary_of_the_review": "This paper introduces a novel paradigm that optimizes on the Grassmannian to complete and cluster incomplete data in a union of subspaces. The core idea is novel and neat, and the paper is well-written. However, there is some unnecessarily technical writeup, some parts appear redundant. And the experimental results section is not sufficient.", "strength_and_weaknesses": "Strength:\n\n1. The core idea that to cluster incomplete data using subspaces as proxies is novel. \n\n2. This paper is well-written and clearly organized.\n\n3. Detailed theoretical derivations and analyses are provided.\n\nWeaknesses\n\n1. The scenario to which the article applies, and the problem it solves, are unclear. What is incomplete data, and what is the relationship between incomplete data and HRMC?\n\n2. The article provides detailed mathematical proofs, but I think they are too much at the expense of readability. I carefully read the article twice but still had difficulties following it.\n\n3. There are not enough experiments to prove the efficiency of the proposed method, and the paper did not provide codes or other supplementary materials for replication.\n", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"correctness": "4: All of the claims and statements are well-supported and correct.", "technical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work.", "empirical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work."}, "extra_text": {}}, {"review_id": "Fcp0gevPJo", "reviewer": "Reviewer_6KMp", "summary": "", "questions": "", "limitations": "", "rating": 1, "confidence": 5, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "The paper proposes to solve the high-rank matrix completion problem by learning a fusion (or a collection) of subspaces. Since the subspaces form the Grassmann manifold, the paper poses the problem as a learning problem on the Grassmann manifold. It proposes a particular cost function that minimizes the distances between all the n subspaces, which correspond to n data points. No technical results are presented. Convergence of the standard stochastic gradient method on manifolds is discussed (which is a known result). Empirical results show some benefits but they are not clear to me. ", "summary_of_the_review": "I consider the paper requires much improvements in order to be considered a viable submission. ", "strength_and_weaknesses": "Pros:\nThis is an interesting problem.\n\nCons:\nI found a number of arguments in the paper very problematic. This is apart from the lack of theoretical basis of the arguments presented in the paper (which by the way appear very hand wavy at best).\n\nIn Page 1, it is mentioned “the geodesics between subspaces of all data, to encourage the subspaces from points that belong together to fuse (i.e, represent the same space)…” as a contribution. However, no results are shown to see whether the modeling actually allows all n points to converge to the same. There seems to be a mis-match between the text statement and the actual work.\n\nNo motivation for the loss function, especially the fusion cost term. Why is it a good choice? Ablation studies are missing.\n\nAs the paper alludes to “high-rank” completion problems, literature survey should be done on what are precise difficulties and what works already exist. \n\nThe baselines are not properly explained. Furthermore, the legend of those are mismatching to those in the text. Also, it would make sense to include [39] as a baseline by considering all the pairs as the paper considers [39] it to be “tightly related”. \n\n\n", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"correctness": "2: Several of the paper’s claims are incorrect or not well-supported.", "technical_novelty_and_significance": "1: The contributions are neither significant nor novel.", "empirical_novelty_and_significance": "1: The contributions are neither significant nor novel."}, "extra_text": {}}, {"review_id": "93m2O-5oUMY", "reviewer": "Reviewer_CCam", "summary": "", "questions": "", "limitations": "", "rating": 6, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "This paper proposes an approach for \"high-rank matrix completion\", a combination of matrix completion with clustering on a union of subspaces. It designs an optimization problem that pushes the (incomplete) observation's subspace to contain a possible completion of the observation, while pushing for different subspaces to get close to each other. Optimization is done by gradient descent on the product of grassmanian manifold, and there is a local convergence theorem. Several numerical experiments are conducted with promising results.", "summary_of_the_review": "While I still think this paper might be a better fit for a journal, the problem formulation is novel and it has interesting theoretical and experimental content that tackles a difficult problem with many applications.", "strength_and_weaknesses": "Strength: This paper builds on well-grounded results in Riemannian manifold optimization. Subspace clustering and matrix completion is solved together, and it doesn't require the number of clusters as an input. It doesn't have low-rank assumption on matrix completion. Experimental results in the low sampling regime are promising.\n\nWeakness: The paper seems to have crammed a lot of material to fit the page limit. It may have been a better fit for a journal where they can define, derive, and present the results much more rigorously and carefully.", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"correctness": "4: All of the claims and statements are well-supported and correct.", "technical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work.", "empirical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work."}, "extra_text": {}}]

Author rebuttal:

Decision: Reject