Connected Subspace Clustering: Hardness, a Scalable Heuristic, and an Application to Sea Level Geodesy

Connected Subspace Clustering: Hardness, a Scalable Heuristic, and an Application to Sea Level Geodesy Johanna Hillebrand∗

Jan Höckendorff†

arXiv:2608.14215v1 [cs.LG] 14 Aug 2026

Melanie Schmidt∗

Jürgen Kusche‡

Christian Sohler†

Kelin Luo§

Heiko Röglin‡

Bernd Uebbing‡

Abstract. Constrained optimization extends classical optimization by integrating side information, making it widely applicable across scientific and engineering domains. Consider a setting where we are measuring a set of variables at different physical locations. When grouping these measurements, it is natural to try to satisfy two different objectives. On the one hand, we want our groups (clusters) to contain similar measurements. On the other hand, we want to have physically coherent clusters. Thus, we have a constrained clustering problem where the constraint models the coherence. Motivated by an application in geodesy, where contiguous regions of the sea surface must be identified for subsequent principal component analysis of sea level anomaly data, we introduce the Connected Subspace Clustering problem. Given highdimensional points and a connectivity graph on them, the objective is to partition the points into k connected clusters while minimizing their total squared distance to the clusters’ best-fit m′ -dimensional affine subspaces. We prove that, even for m′ = 0 and when the connectivity graph is a grid graph with holes, this problem is NP-hard to approximate within a factor of Ω(n1/2−ε ) for every ε > 0, where n is the number of measurements in the input. We then introduce an efficient Lloyd-style heuristic that alternates subspace fitting with an iterative component-merging procedure to enforce connectivity. Our method returns exactly k connected regions by construction, whereas unconstrained subspace clustering methods leave up to 1,966 disconnected fragments and substantially higher reconstruction cost. In an experimental study of 160 configurations on global sea level time series, our merging-based connectivity repair is the strongest of four strategies in 73.75% of configurations, and consistently outperforms competitors such as (connected) Ward’s method across all tested cluster counts. The resulting regions isolate signals that align with established climate indices such as those describing the El Niño–Southern Oscillation and Indian Ocean Dipole. Although developed for sea level geodesy, the approach applies in principle to other spatially embedded multivariate time series, such as climate fields, remote sensing imagery, neuro-imaging, and sensor networks.

∗ Heinrich Heine University Düsseldorf, Universitätsstraße 1, 40225 Düsseldorf, Germany (johanna.hillebrand@hhu.de). † University of Cologne, Albertus-Magnus-Platz, 50923 Köln, Germany (hoeckendorff@cs.uni-koeln.de). ‡ University of Bonn, Regina-Pacis-Weg 3, 53115 Bonn, Germany. § University at Buffalo, 12 Capen Hall, Buffalo, New York 14260-1660, USA (kelinluo@buffalo.edu).

1 Introduction Clustering is a widely used machine learning tool, applied in fields as diverse as medical imaging and astronomy. In practice, additional domainspecific structure often imposes constraints on how clusters can be formed. Consequently, various forms of constrained clustering, such as capacitated or fairnessaware clustering, have been discussed extensively in the recent literature. Connectivity-constrained clustering naturally arises when handling spatial datasets, which are common in many fields including remote sensing, neuro-imaging and sensor networks, and it can be described as the requirement of identifying contiguous regions with coherent dynamics. More formally, it is defined as a family of constrained clustering problems where a connectivity graph G = (V, E) is given in addition to the input (the standard input contains the points, the target number of clusters, and the metric) [19]. The connectivity graph is unweighted and undirected. The constraint is that for every cluster, the induced subgraph in G has to be a connected graph. The graph has no guaranteed correspondence with the metric, so points u, v that share an edge {u, v} can still be very far away from each other with respect to the metric. This sort of constraint has a long history. In a survey from 1985, Murtagh [35] reviews applications from social sciences, earth sciences, pattern recognition and image processing and describes contiguity-constrained agglomerative clustering methods. The idea was rediscovered in 2008 by Ge et al. [19] under the name connected clustering, in the context of community detection which spurred subsequent work on connectivity-constrained clustering [11, 14, 21]. Another important line of work is subspace clustering which addresses the fact that high-dimensional data often contains meaningful patterns that exist only within specific subsets of dimensions. In particular, with modern sensing technologies, the number of measured dimensions grows quickly. Traditional clustering methods may fail to detect these patterns due to the curse of dimensionality [26, 38]. Restricting attention to fewer, more relevant directions improves accuracy and interpretability in complex datasets. Thus, subspace clustering has become especially valuable in fields like bioinformatics, image analysis, motion segmentation and market research. Formally, it is defined as follows: given points P ⊂ Rm , a number k and a target dimension m′ , find k affine subspaces of dimension m′ that minimize the sum of the squared distances of points to their closest subspace. Representative approaches include CLIQUE, using grid-based dense units [2], PROCLUS, using projected k-medoids with per-cluster relevant dimensions [1], SUBCLU a DBSCAN-style algo-

rithm, which uses density-based clusters in subspaces defined via density-reachability [23] and numerous other methods [12]. When datasets are both structured and highdimensional, the ideas of subspace clustering and connectivity naturally meet and thus, we introduce Connected Subspace Clustering, a problem that combines the graph-topological constraint of connectivity with the statistical objective of minimizing the subspace reconstruction error in one optimization problem. 1.1 Our Contributions (i) We introduce Connected m′ -Subspace k-Clustering, which combines affinesubspace reconstruction error with connectivity in a given graph. (ii) We prove that, even for m′ = 0 and grid graphs with holes, the problem is NP-hard to approximate within a factor of Ω(n1/2−ε ) for every ε > 0. (iii) We design a generic (heuristic) merging subroutine that converts fragmented cluster assignments into connected regions and can be incorporated into Lloyd-style k-subspaces algorithms. On bounded-degree connectivity graphs, the routine runs in O(n(k + log n)) time. (iv) We provide extensive experimental evaluation for an important application in sea level geodesy, comparing our algorithm to methods from the areas of classical clustering, subspace clustering and hyperspectral imaging. 1.2 Application The motivation for our clustering problem stems from sea level geodesy, in particular from the analysis of climate-driven sea level variability. Contemporary sea level rise is predominantly driven by climate change via glacier/ice sheet melt and volumetric expansion (ocean warming) [18]. Since 1993, satellite altimetry has provided (near) global records of sea level anomaly (SLA) with centimeter-level accuracy on an almost gap-free ocean grid, excluding only land areas and parts of the high latitudes. This data enables the study of sea level variability on global as well as regional scales. A standard tool for analyzing such data is principal component analysis (PCA), which decomposes the 3-dimensional SLA field into spatial patterns known as empirical orthogonal functions (EOFs) and associated temporal components, both ordered by explained variance. On the global domain this ordering by variance favors large-scale, dominant signals such as the longterm sea level trend and the annual cycle (e.g., seasonal land-water exchange). Regionally important but lowervariance phenomena, such as the El Niño–Southern Oscillation (ENSO), therefore tend to appear only in higher-order modes (Fig. B.1). Thus, analyses performed on the global field can obscure spatially localized

(a) Greedy k-means++

(b) Ward’s method

(c) SSC-OMP

(d) EGCSC

(e) Connected Ward’s method

(f) Connected Subspace Clustering

Figure 1.1: Clustering of the sea level anomaly dataset for k = 25. Unconstrained clustering and subspace clustering methods produce spatially fragmented clusters (a-d), whereas connectivity-constrained methods (ef) produce spatially coherent clusters. Connected Subspace Clustering (f) produces spatially coherent and interpretable clusters, while also achieving the lowest subspace-clustering cost among the displayed methods. signals. For this reason, regional analyses are standard in sea level research [30]. Such analyses require partitioning the ocean surface into subregions that exhibit internally similar SLA variability and then can be studied separately, for example by applying PCA on each region. Partitions are often defined manually or based on prior geographic knowledge. In this work, we replace such a partitioning with data-driven clustering of altimetry time series. The resulting clusters group ocean areas with similar SLA behavior and provide a basis for downstream tasks such

as reconstructing pre-satellite-era sea level from tide gauges [8] and separating steric and mass contributions [37]. Because the downstream analysis applies PCA separately within each region, a natural partitioning objective is the total rank-m′ reconstruction error over all regions. This is exactly the sum of squared distances to the regions’ best-fit m′ -dimensional affine subspaces. Requiring each cluster to be geographically contiguous therefore leads directly to our Connected Subspace Clustering problem. Existing methods do not simultaneously optimize

subspace reconstruction error and enforce connectivity in a prescribed spatial graph. Figure 1.1 shows the output of six methods, including our own. Greedy kmeans++ and Ward’s method (scikit-learn [36]) yield spatially fragmented regions, but spatial contiguity is crucial for oceanographic analysis (Figure 1.1(a,b)). Typical subspace clustering algorithms also do not enforce connectivity and thus can produce highly fragmented clusters (Figure 1.1(c)). Related work in hyperspectral imaging (HSI) also handles high-dimensional data under spatial coherence constraints e.g. [6, 42, 20]. While it shares these challenges with our setting, there are key differences in goals and problem formulation. Most importantly, it encourages spatial coherence but does not strictly enforce it while our problem formulation requires strictly connected clusters. Applying HSI on ocean data does not produce connected clusters (Figure 1.1(d)). A well-established connected clustering method is connected Ward’s method, which is a rather old agglomerative clustering method ([35, 29]). It provides coherent clusters, but yields uneven patches that do not align well with oceanographic phenomena (Figure 1.1(e)). The regions resulting from our method, in comparison, are both contiguous and interpretable (Figure 1.1(f)). 1.3 Related Work We are not aware of any algorithm specifically designed for our problem formulation of connected subspace clustering. For subspaces of dimension zero, i.e., for k-means, and generally for point-based clustering problems, some results on connected clustering are known. Connected k-center: The radius-based k-center problem is somewhat easier to tackle than sum-based objectives like k-means. Surprisingly, it is not O(1)-approximable: The connected kcenter problem is hard to approximate in o(log∗ k) [15]. On trees, the connected k-center problem is solvable optimally [11, 19]; for general graphs, there exists a O(log2 k)-approximation [11]. Connected k-median/kmeans: For general graphs, connected k-median is NPhard to approximate within Ω(n1−ϵ ), but is solvable on trees [14]. The hardness extends to Euclidean k-means and thus to subspace clustering. On star-like graphs, connected k-means is NP-hard to approximate within Ω(log n) and admits an O(log n)-approximation [21]. To our knowledge, the case where the connectivity graph is a grid has not been studied before this work. ILP formulations: The point-based connected clustering problem has also been modeled as an integer linear program in the area of cartography in the context of political redistricting (e.g. [41]). While ILP solvers nowadays are quite fast, this approach does not scale to our data. As mentioned above, we are not aware of any algorithms proposed for subspace clustering where the con-

nectivity is explicitly enforced but now review related algorithms and heuristics. The approaches need to be able to handle larger data sizes (the ocean dataset includes 532,783 time series with 365 monthly time steps). Iterative improvement heuristics: Ge et al. [19] propose a heuristic called NetScan for connected k-center and adapt it to connected k-means. The broad idea of this algorithm is: initialize centers, assign while enforcing connectivity via simultaneous breadth-first search (BFS) with a look-ahead for “bridge”-nodes, then recompute better centers for each subgraph and iterate this process until a convergence criterion is reached. This approach does not scale well, and there is no implementation available. We implement a related Lloyd-style scheme with a more efficient assignment step. Connected Ward’s method: Agglomerative clustering is a general approach that starts with singleton clusters and then subsequently merges a pair of clusters into one cluster until the desired number of clusters is reached. This can naturally be adapted to respect connectivity by only allowing merges of clusters that are neighboring. We use the implementation of this method from scikit-learn [36] that is based on Murtagh [35] and Li [29]. Subspace clustering methods: Subspace clustering is a large and active research area. We choose two well-established methods to investigate the connectivity of the produced solutions: Sparse subspace clustering with orthogonal matching pursuit (SSC-OMP) [45], which addresses subspace clustering without any connectivity constraints, and elastic net subspace clustering (EnSC) [44], which encourages well-connected affinity graphs within recovered subspaces, but it does not enforce connectivity in an externally specified spatial graph. Hyperspectral imaging (HSI): Hyperspectral images typically are 3D tensors (two spatial dimensions and one spectral dimension). Many spatial-spectral clustering methods developed for HSI data encourage spatial coherence via a regularization term, which encourages neighboring pixels to have the same label. Such regularization promotes local coherence but generally does not require every cluster to induce a single connected component. This is the method of choice for, e.g., identifying landscape types. We nevertheless compare our method with two popular HSI approaches [6], efficient graph convolutional subspace clustering (EGCSC) and efficient kernel graph convolutional subspace clustering (EKGCSC). 2 Problem Definition, Hardness and Algorithmic Overview We define the connected subspace clustering problem as follows: Let P ⊂ Rm be a set of m-dimensional points with |P | = n (corresponding to n time series of length m). We also think of this point set as an n × m-matrix where the i-th point is in the i-th row. Let d : Rm × Rm → R denote the Euclidean

used, with no further structural assumptions. In our experiments, N is a grid graph with holes as dictated by our application. The routine first decomposes every cluster into its connected components. It then repeatDefinition 2.1. (Feasible clustering) For a num- edly selects a smallest component and merges it into an ber of target clusters k, a partition of P into clusters adjacent cluster, until k connected components remain. C1 , . . . , Ck is feasible if the subgraph of G induced by Ci The target cluster is chosen according to which neighboring cluster subspace best explains the points of the is connected for all i ∈ {1, . . . , k}. selected component. The routine can be plugged into The purpose of the clusters from the point of various clustering algorithms to establish connectivity, the application is that computing a PCA with m′ see Section 3. components on a cluster shall induce a smaller error 2.3 Algorithmic Framework: Our heuristic than computing it on the global dataset, in order to enable better analysis of sea level behavior for each local combines an initial clustering with a Lloyd-style iterative improvement procedure. In each iteration, we comarea. ′ To this end, we want to minimize the sum of squared pute a best-fit m -dimensional affine subspace for every distances of points to the best-fit subspace of their cluster and reassign each point according to its reconcluster. The best-fit subspace of a set of points Ci ⊂ Rm struction error. Since these reassignments may violate is the affine subspace of a given dimension m′ that connectivity, we consider several strategies for mainminimizes the sum of squared Euclidean distances from taining or restoring connected clusters, including iterative component merging, postprocessing, spatial filterthe points in Ci to the subspace. P It is obtained by computing the mean µi = |C1i | x∈Ci x, the centered ing, and locally restricted reassignment. The framework is intentionally modular: different initialization methdata matrix X = [x − µi : x ∈ Ci ] ∈ Rm×|Ci | , and then ods and connectivity strategies can be combined within the singular value decomposition of X. Then the bestthe same iterative procedure. Section 3.1 presents these fit subspace is µi + span(u1 , . . . , um′ ), where u1 , . . . , um′ methods in detail. Section 4 compares the resulting conare the top m′ left singular vectors of X, which are the figurations on sea level anomaly data. m′ top eigenvectors of XX T . 3 An Algorithm for Connected Subspace Definition 2.2. (Connected m′ -subspace kClustering We begin with the general subspace clus′ clustering) Given P , G, k, and m , find a featering algorithm (Subsection 3.1), followed by two key sible clustering C1 , . . . , Ck of P with correspondcomponents: the initialization (Subsection 3.2) and a ′ ing best-fit Pk Psubspaces Si2 of dimension m such method to (re)establish connectivity (Subsection 3.3). that d(x, S ) is minimized, where distance which we use as a distance measure for time series. Let G = (P, E) denote an unweighted and undirected connectivity graph. We want to partition G into connected subgraphs.

i=1

x∈Ci

i

d(x, Si ) = d(x − µi , span(u1 , . . . , um′ )). ′

2.1 Hardness We prove that the connected m subspace k-clustering problem is hard to approximate to a constant even for m′ = 0. This is even true if the connectivity graph is a grid graph with holes which is the situation present in the application that we study in Section 4. We show this via a non-trivial reduction from the vertex-disjoint paths problem on grid graphs. The proof of the following theorem is given in the supplementary material, see Appendix A.2. Theorem 2.3. For any k ∈ N, connected m′ subspace k-clustering with m′ = 0 is NP-hard to approximate within factor Ω(|V ′ |1/2−ε ) for any ε > 0, on grid graphs G′ = (V ′ , E ′ ) with holes in Euclidean tdimensional space, where t ≥ k. 2.2 Establishing Connectivity: We design a generic merging subroutine for clusterings with fragmented clusters. It operates on a neighborhood graph N over the input points; N may be any finite, simple, undirected graph, as only adjacency information is

3.1 Basic Algorithm The general structure of our algorithm follows a Lloyd-style iterative improvement strategy: Alternately, assign points to clusters and update the clusters’ optimal representation. Lloyd-type methods strictly decrease the objective each step but can take exponentially many iterations in the worst case (e.g., even for k-means [4]). Stronger convergence guarantees require additional assumptions in the input [43]. We make no such assumptions and cap the number of iterations to ensure polynomial running time. This works ′ well for subspace clustering P because the m -dimensional subspace S minimizing x∈Ci d(x, S) is just the best-fit subspace of dimension m′ and it can be computed via principal component analysis. This procedure is given in Alg. 3.1. Several steps need further specification. First and foremost, we need a good initial clustering to start with. It is known that for k-means, the initialization is of utmost importance. We discuss several options in Section 3.2. Second, we need to specify a convergence criterion. We use the following simple criterion: Either terminate after a fixed number of iterations or when as-

Algorithm 3.1 Connected Subspace Clustering

tering, but with Ward’s objective function. This algo1: Require: P = {x1 , . . . , xn }, number of clusters k, rithm was publicly available and we use the implementation of scikit-learn [36]. subspace dimension m′ 2: Ensure: Clusters C1 , . . . , Ck with best-fit sub- 4) Greedy k-means++ with connectivity (ConnKMeans++) As a final initialization method, we use spaces S1 , . . . , Sk of dimension m′ scikit-learn’s implementation of greedy k-means++ [36], 3: Compute an initial clustering C1 , . . . , Ck a variant of k-means++ by [5], reported to deliver very 4: while stopping criterion not met do good results, e.g. [7]. For completeness, pseudocode is in 5: for i = 1 to k do Subsection C.1. Since k-means does not enforce spatial 6: Compute the best-fit subspace connectivity, we post-process the output with Alg. 3.2 7: Si ← µi + span(v1 , . . . , vm′ ) ′ to obtain a connected clustering. 8: C ←∅ i

for x ∈ P do Let j ∈ arg min1≤i≤k d(x, Si ) 11: Cj′ ← Cj′ ∪ {x} 12: for i = 1 to k do 13: Set Ci ← Ci′ 14: Optional: Establish connectivity for Ci 15: Optional: Establish connectivity for all Ci

3.3 Methods for (re)Establishing Connectivity Assigning points to best-fit subspaces in Step 11 of Algorithm 3.1 can violate the connectivity constraint. We evaluate the following options: (i) Reestablish connectivity after each iteration via Algorithm 3.2 (IterMerge), (ii) establish connectivity once at the end via the same algorithm (PostMerge), (iii) apply a filtering step that improves, but does not strictly enforce connectivity (SmoothMerge), and finally (iv) restrict Step 11 signments no longer change. Third, it is unclear how to maintain connectivity throughout (IntegratedConn). to maintain or (re)establish connectivity for the clus- For an overview of the different options, see Table 3.1. ters. This can be done after each iteration or once after We discuss each option in detail below. the procedure terminated. There is a big caveat to this Method During iteration Final repair question, since enforcing connectivity can increase the IterMerge Component merging No objective function (Section 3.3). 9: 10:

3.2 Initial Clustering We employ five initialization methods based on four clustering approaches. Two methods use the same connectivity-constrained agglomerative procedure with different dissimilarity measures. 1) Agglomerative clustering without a connectivity constraint (Agglo-ST) Thompson and Merrifield [40] apply average-linkage clustering with a spatiotemporal distance function D(Xi , Xj ) = 1 − exp(−d(xi , xj )/c0 )r(ti , tj ), where Xi = (xi , ti ) combines location and time series. The exponential term penalizes spatial separation, while Pearson’s r rewards temporal similarity. Connectivity is only implicit and thus, connected clusters appear only when the data are smooth in space and time. This is explained in greater detail in the appendix (see Subsection B.2). 2) Agglomerative connected clustering The main difference to standard agglomerative clustering is the restriction of merges to pairs of neighboring clusters. We defer the pseudocode to Appendix C.1. We implemented Algorithm C.1 P and use P it with (i) average 1 linkage D(A, V ) = |A||V x∈A y∈V D(x, y) (Conn| Agglo-Euc) and (ii) with the spatiotemporal distance function described in 1) (Conn-Agglo-ST ). 3) Ward’s method with connectivity (Conn-Ward) This follows the same idea as agglomerative connected clus-

PostMerge None SmoothMerge Label smoothing IntegratedConn Local reassignment

Yes Yes Yes

Table 3.1: Connectivity strategies. “Final repair” indicates whether Algorithm 3.2 is applied after the iterative procedure terminates. 3.3.1 Merging Connected Components To restore connectivity (Alg. 3.2), we form the induced subgraph on the neighborhood graph N by retaining only edges whose endpoints share a cluster label. The connected components of this subgraph are contiguous subregions of clusters (see Figure 3.1(a)). We then repeatedly select the smallest component and merge it with an adjacent cluster whose fitted subspace is preferred by most of its points, until k components remain (see Figure 3.1). If a reassignment reduces the number of components below k, the connected components can be split along edges of a spanning tree until exactly k components remain. All ties in arg min and arg max expressions are broken arbitrarily but consistently. The procedure is motivated by the expectation that an unconstrained assignment typically consists of several substantial contiguous regions together with smaller disconnected fragments. More details are in Subsection C.2.1, and the proof of Theorem 3.1 is

(a) Before merging

(b) After merging

Figure 3.1: Application of Alg. 3.2 to an assignment of grid points for 15 clusters. The merging step enforces connectivity by reassigning disconnected subregions, prioritizing spatial coherence over original cluster labels. in Subsection C.2.2. An optimized implementation comAlgorithm 3.2 Merging Connected Components 1: Require: Possibly

disconnected clusters C1 , . . . , Ck , subspaces S1 , . . . , Sk , the connectivity graph N 2: Ensure: Connected clustering C̃1 , . . . , C̃k 3: Construct the graph GC = (V, E ′ ) where E ′ = {{u, v} | {u, v} ∈ E and ∃i : u, v ∈ Ci } 4: Compute set Z of connected components of GC 5: Compute the component graph GCC = (Z, E ′′ ) where E ′′ = {{A, B} | ∃u ∈ A, ∃v ∈ B, A, B ∈ Z, {u, v} ∈ N } 6: while |Z| > k do 7: Compute Amin ∈ arg minA∈Z |A| 8: for x ∈ Amin do 9: Compute the winner w(x) = arg min{d(x, Si ) | i ∈ [k] s. t. ∃B ⊂ Ci , B ∈ Z, {Amin , B} ∈ E ′′ } 10: Compute iw = arg maxi |{x ∈ Amin | w(x) = i}| 11: Reassign all points in Amin to cluster Ciw 12: Update graphs GC and GCC 13: For each A ∈ Z create a cluster C̃ and output these

3.3.2 Improving Connectivity via Filtering Enforcing connectivity at every step has the side effect of a fluctuating cost (the reassignment decreases cost, while enforcing connectivity raises it). For SmoothMerge, we encourage spatial coherence after each reassignment step using Gaussian weighted label voting. For each point x, nearby grid points contribute to the support of their current cluster label, with weights decreasing with their Haversine distance from x. We assign x the label with maximum weighted support. This operation generally reduces fragmentation but does not guarantee connectivity, so Algorithm 3.2 is applied after the iterative procedure terminates. Full details are given in Subsection C.2.3. 3.3.3 Integrating Connectivity Constraints We also evaluate incorporating local connectivity information directly into the reassignment step, rather than relying exclusively on connectivity repair afterward. For this we replace Step 10 and 11 of Alg. 3.1 with Alg. 3.3.

Algorithm 3.3 Integrated connectivity for each point x ∈ P do for each i = 1 to k such that there exists y ∈ Ci with y ∈ δN (x) ∪ {x} do Compute d(x, Si ) = d(x − µi , span(v1 , . . . , vm′ )) Assign x to the cluster Ci′ that minimizes d(x, Si )

putes the initial connected components by BFS/DFS, precomputes point-to-subspace rankings and component score vectors, and maintains the smallest component in a min-heap. Component merges use union–find with union-by-rank and path compression, while winner The sole change is a locality restriction: x may lookups and score-vector updates cost O(k) per iteraswitch from Ci to Cj only if a neighbor of x in N lies tion; these ingredients yield the bound below. in Cj and d(x, Sj ) < d(x, Si ). Among these candidate Theorem 3.1. Algorithm 3.2 can be implemented labels, it chooses the cluster whose fitted subspace is in O(|E| + n(k + log n)) time. Consequently, its run- closest. This locality restriction encourages spatially ning time is O(n(k + log n)) on bounded-degree graphs, coherent assignments but does not guarantee that every including the grid graphs used in our application, and cluster remains connected. We therefore apply the final O(nk) whenever k = Ω(log n). merging step in Step 15.

4 Experimental Evaluation We evaluate our Connected Subspace Clustering framework on global sea level anomaly data, with the goal of identifying contiguous regions exhibiting internally similar temporal variability. 4.1 Dataset, Preprocessing and the Connectivity Graph We use the Copernicus Marine Service (CMEMS) global gridded sea level anomaly product [13]. We choose the CMEMS product because its enhanced spatial sampling and time-step accuracy are well suited to identifying regional spatiotemporal patterns. This dataset contains monthly sea level anomaly (SLA) values on a 0.25◦ × 0.25◦ regular, rectilinear grid (720 latitudes and 1440 longitudes). We restrict the analysis to grid cells with complete SLA observations throughout the study period. Owing to the dataset’s validity mask, this excludes land, many coastal cells, and parts of the high latitudes, leaving 532,783 valid ocean grid points. For tractability, all reported clustering experiments use this grid coarsened to a 2◦ × 2◦ resolution, yielding 8,160 grid points. Following the preprocessing strategy of [40, 16], we apply a temporal boxcar filter (half-width seven months; approximately 90% gain at a 24-month period) and a spatial Gaussian filter based on Haversine distance, thereby accounting for the Earth’s curvature. The spatial filter is applied independently to each monthly field using normalized Gaussian weights with half-width h = 500, truncated at a radius of 3σ, where σ = h/1.178. Masked grid cells are excluded from the weighted average. All methods are evaluated on both unfiltered and filtered data. More on the dataset can be found in Subsection B.1 and further preprocessing details are described in Subsection D.1.

(a) k-means cost function

(b) Subspace clustering cost for m′ = 15 Figure 4.1: Comparison of our five initial clustering methods evaluated with two different cost functions. README provides complete environment setup, configuration, and execution instructions for reproducing the reported experiments. The experiments were conducted on a single machine with an AMD Ryzen 7 PRO 7840U processor (8 cores, 16 threads) and 30 GiB of RAM.

4.1.3 Method Names We use the abbreviations Agglo-ST, Conn-Agglo-Euc, Conn-Agglo-ST, Conn-KMeans++, and Conn-Ward for the five initialization methods from Subsection 3.2. The four con4.1.1 Connectivity Graph We map the com- nectivity strategies are IterMerge, PostMerge, Smoothbinatorial structure of the data to a connectivity graph Merge, and IntegratedConn; see Table 3.1. N = (V, E) with one vertex per grid point and edges for 4.2 Initializer Quality under k-Means and each of the four neighbors in the grid. Because masked Subspace Objectives The k-means problem is a specells are omitted, N is a grid graph with holes. cial case of subspace clustering with m′ = 0, and we To account for the periodicity of longitude, retained use k-means as an initial clustering. When evaluating cells in the westernmost and easternmost grid columns the quality of the initial clusterings, we made an interare treated as adjacent when they occur at the same esting observation: For Figure 4.1, we ran all five inilatitude. Thus, N can be viewed as a subgraph of a twotialization methods with an increasing number of clusdimensional grid on a cylinder. Optionally, diagonal ters. For each initialization method, we evaluate the reedges can be added; all algorithms described in this sulting partition before applying the iterative subspacepaper support either choice. In the evaluation the fourclustering procedure. We score each fixed partition once neighbor graph is used. The set of neighbors of a vertex using the k-means objective and once using the rank-15 x ∈ V is δN (x) = {y | {x, y} ∈ E}. subspace-reconstruction objective. Under the k-means 4.1.2 Hardware and Source Code We imple- objective, Conn-Ward and Conn-KMeans++ outpermented all algorithms in Python. The anonymized form the other approaches, see Figure 4.1(a). However, source code is available via https://anonymous.4open. under the subspace clustering objective, the ranking science/r/subspace_clustering-C00F. The repository flips (Figure 4.1(b)). Notably, agglomerative connected

clustering outperforms all other methods, and the tailored method by Thompson and Merrifield (Agglo-ST ) is now competitive with Conn-KMeans++. Conn-Ward has comparatively high subspace-reconstruction cost, particularly for smaller values of k. This trend is consistently observed across multiple configurations (see Subsection D.4). These results show that initializer quality is strongly objective-dependent: methods that perform well under the k-means objective do not necessarily provide good starting partitions for the rank-m′ reconstruction objective used in our downstream regional PCA analysis. We examine the geographic and oceanographic interpretation of the resulting regions in Subsection 4.6.

Table 4.1: Initial and final subspace-reconstruction objectives for IterMerge on filtered data with k = 8 and m′ = 15. For this representative setting, IterMerge improves the objective value for all five initialization methods. Initial objective Final objective Agglo-ST Conn-Agglo-Euc Conn-Agglo-ST Conn-KMeans++ Conn-Ward

6058.8 5892.1 5925.4 6019.4 6250.7

5681.5 5656.5 5716.8 5742.0 5728.4

Even for larger k, Agglo-ST consistently underperforms (Tables D.6, D.8, and D.10); e.g., at k = 20, it is at least 8.98% higher than the next worst clustering method, Conn-Ward. The other four initialization methods produce usable starting partitions on both filtered and unfiltered data for all tested configurations.

(a) Clustering on filtered input data.

(b) Clustering on unfiltered input data. Figure 4.2: Agglo-ST initialization for k = 25 on (a) filtered and (b) unfiltered SLA data. Without filtering, Agglo-ST produces one dominant cluster and several small regions, some of which are too small to fit the requested subspace dimension. 4.3 Robustness / Filtering We next examine how spatial and temporal filtering affect the performance of the initialization methods. Agglo-ST is highly sensitive to preprocessing. On unfiltered data, it frequently degenerates into one dominant cluster and several very small clusters (Figure 4.2). Some of these clusters contain too few points for our PCA implementation to fit an m′ -dimensional subspace, so the resulting initialization cannot be used by Conn-Subspace for those parameter settings.

4.4 Comparison of Connectivity Strategies We investigate how different strategies for reestablishing connectivity (see Subsection 3.3) during subspace clustering impact the objective. Therefore, we compare the four connectivity strategies across five initializations, four subspace dimensions m′ ∈ {5, 10, 15, 30}, four cluster counts k ∈ {8, 15, 20, 25}, and both filtered and unfiltered data. This yields 160 base configurations. Each base configuration is evaluated with all four connectivity strategies, resulting in 640 runs. For each base configuration, we compare the final subspace-reconstruction objective attained by the four strategies. Overall, IterMerge attains the lowest final objective value in 118 of the 160 base configurations (73.75%). SmoothMerge and IntegratedConn follow, PostMerge is worst. In two of the 160 base configurations (1.25%), Agglo-ST on unfiltered data produces clusters that are too small to fit the requested subspace dimension; no refined solution is therefore available for these configurations (Subsection 4.3). In addition to outperforming the other connectivity strategies most frequently, IterMerge lowers the objective value relative to the corresponding initialization in every valid configuration. Table 4.1 shows one representative setting; complete results are given in Subsection D.4. As illustrated in Figure 4.3, because connectivity is enforced only at the end, SmoothMerge and especially PostMerge show a final cost spike. The iterative smoothing step in SmoothMerge mitigates this effect, but does not remove it. IterMerge and IntegratedConn behave more like Lloyd’s method, where costs improve through the iterations. Even though the strict connectivity enforced by Alg. 3.2 in IterMerge can lead to a temporary increase in the cost function after

(a) IterMerge

(b) IntegratedConn

(c) SmoothMerge

(d) PostMerge

Figure 4.3: Cost function trajectories for the four connectivity heuristics over algorithm iterations. IterMerge and IntegratedConn show relatively continuous improvements, while SmoothMerge and PostMerge show late-stage cost increases due to delayed enforcement of connectivity. the merging, it is the approach that most successfully from suboptimal starts; the largest reduction is 20.3% balances the connectivity constraint with the subspace (Agglo-ST initialization, unfiltered data, k = 15, clustering objective. m′ = 30). Among the configurations compared in Figure 4.4, Conn-Agglo-Euc followed by Conn-Subspace with IterMerge attains the lowest objective value for every tested k, on both filtered and unfiltered data. Conn-Ward is the connected clustering method most readily available prior to this work. While Conn-Ward relies on a hierarchical agglomeration based on local variance, it lacks the ability to revise cluster assignments after initial formation. In contrast, Conn-Subspace can iteratively correct misassignments. 4.5.1 Comparison with SubspaceClustering and Spatial-Spectral Methods We additionally compare against two unconstrained subspace-clustering methods, SSC-OMP and EnSC, and two spatial-spectral methods from hyperspectral image segmentation, EGCSC and EKGCSC. None of these methods guarantees that each output cluster 4.5 Comparison with Connected and induces a single connected component in the prescribed Subspace-Clustering Baselines Across the valid graph. Implementation and parameter details are configurations, Conn-Subspace consistently lowers given in Subsection D.2. SSC-OMP and EnSC both the objective value relative to its initialization, even yield higher subspace-reconstruction objectives and Figure 4.4: Subspace-reconstruction objective as a function of the targeted number of clusters k for ConnWard and Conn-Subspace initialized with Conn-AggloEuc and using IterMerge. Results are shown for filtered and unfiltered data.

more fragmented partitions than our method. For example, on unfiltered data with k = 25, SSC-OMP attains an objective value of 14,473.1, compared with 2,947.9 for Conn-Agglo-Euc followed by Conn-Subspace with IterMerge at m′ = 30. Despite only 25 targeted clusters, SSC-OMP produces up to 1,966 connected components; EnSC is less fragmented (at most 222) but its objective value stays consistently above ours (see Table C.1 and subsection D.4). EGCSC and EKGCSC yield less fragmentation than SSC-OMP and EnSC and attain lower objective values, yet in every tested configuration they still produce more connected components than requested clusters and remain worse than our method (Table D.2 and subsection D.4). In contrast, our method returns exactly k connected clusters.

(a) Signal from cluster 8 (red) together with the ENSO index (blue).

(b) Signal from cluster 9 (red) together with the 4.6 Implications for Sea Level Research ENSO index (blue). Combining subspace clustering with Empirical Orthogonal Function (EOF) analysis isolates regionally coherFigure 4.6: ENSO index comparison ent sea level variability that global EOFs can blur. We illustrate this using the best-performing configuration, Conn-Agglo-Euc followed by Conn-Subspace with Iter- the second EOF. The leading component of cluster 13 Merge, on filtered data with k = 15 (see Figure 4.5). similarly shows correspondence with the Indian Ocean For example, cluster 1 directly relates to the Agulhas Dipole (IOD) index, although the cluster extends beyond the IOD region and contains additional variability. This regional separation supports targeted analysis of sea level drivers and geodetic applications in which regional altimetry EOFs can improve tide-gauge reconstruction [24, 8]. Potential inconsistencies at cluster boundaries warrant further study.

Figure 4.5: Clustering achieved by Conn-Subspace with IterMerge, for k = 15. western boundary current, where water is transported southward along the South African coast and then diverted eastward back into the Indian Ocean. Similarly, cluster 12 is related to the Gulf Stream region in the Atlantic Ocean. While the trend and annual cycle dominate the global EOF decomposition (Appendix, Figure B.1(a),(b)), regional decomposition reveals locations in which other modes dominate. Cluster 8 spans much of the Atlantic basin, whereas cluster 9 covers the eastern tropical Pacific. The leading principal component of cluster 9 shows clear visual correspondence with the El Niño–Southern Oscillation (ENSO) index, wheras that of cluster 8 shows no comparable correspondence (Figure 4.6). This illustrates how the clustering distinguishes regions according to their dominant temporal variability. In cluster 9, the annual cycle shifts to

5 Conclusion We introduced Connected Subspace Clustering, established strong inapproximability even for zero-dimensional subspaces and grid graphs with holes, and developed a scalable Lloyd-style heuristic with graph-based connectivity guarantees. Experiments on global sea level anomaly data show that the bottom-up merging (IterMerge) routine is strongest of the four connectivity strategies and that Conn-Subspace consistently improves its initialization. Unlike the evaluated baselines, it returns exactly k connected regions while attaining substantially lower reconstruction cost. The resulting regions also expose interpretable regional modes, including ENSO-related variability. Limitations of the proposed approach include sensitivity to the neighborhood graph and to the initialization; promising directions for the future include further study of the simultaneous optimization of the objective function and connectivity (bi-criteria optimization) and the extension to other datasets or settings. Overall, combining the subspace clustering objective with topological information offers a practical route to achieve structure-aware clusterings.

References [1] C. C. Aggarwal, J. L. Wolf, P. S. Yu, C. Procopiuc, and J. S. Park, Fast algorithms for projected clustering, ACM SIGMoD record, 28 (1999), pp. 61–72.

[12] E. Elhamifar and R. Vidal, Sparse subspace clustering: Algorithm, theory, and applications, TPAMI, 35 (2013), pp. 2765–2781.

[13] E.U. Copernicus Marine Service Information (CMEMS). Marine Data Store (MDS), Global Ocean Gridded L 4 Sea Surface Heights And [2] R. Agrawal, J. Gehrke, D. Gunopulos, and Derived Variables Reprocessed 1993 Ongoing, 2024. P. Raghavan, Automatic subspace clust. of high acc. on 2024-09-11. dim. data for data mining appls., in Proc. SIGMoD, [14] J. Eube, K. Luo, D. Reineccius, H. Röglin, 1998, pp. 94–105. and M. Schmidt, Connected k-median with disjoint and non-disjoint clusters, in Proc. ESA, 2025. [3] D. Aloise, A. Deshpande, P. Hansen, and P. Popat, Np-hardness of euclidean sum-of- [15] J. Eube and H. Röglin, New algorithms and squares clustering, Machine Learning, 75 (2009), hardness results for connected clustering, CoRR, pp. 245–248. abs/2511.19085 (2025). [4] D. Arthur and S. Vassilvitskii, How slow is the [16] Y. Ewerdwalbesloh, Finding coherent regions k-means method?, in Proc. SoCG, 2006, pp. 144– in the ocean in cmip6 model simulations, master’s 153. thesis, Rheinische Friedrich-Wilhelms-Universität Bonn, 2022. [5] D. Arthur and S. Vassilvitskii, k-means++: The advantages of careful seeding, in Proc. SODA, [17] M. Formann and F. Wagner, The VLSI layout problem in various embedding models, in WG, 2007, pp. 1027–1035. Springer, 1990, pp. 130–139. [6] Y. Cai, Z. Zhang, Z. Cai, X. Liu, X. Jiang, and Q. Yan, Graph convolutional subspace clus- [18] B. Fox-Kemper, H. Hewitt, C. Xiao, G. Aðalgeirsdóttir, S. Drijfhout, T. Edtering: A robust subspace clustering framework for wards, N. Golledge, M. Hemer, R. Kopp, hyperspectral image, TGRS, 59 (2020), pp. 4191– G. Krinner, A. Mix, D. Notz, S. Nowicki, 4202. I. Nurhati, L. Ruiz, J.-B. Sallée, A. Slangen, and Y. Yu, Ocean, Cryosphere and Sea [7] M. E. Celebi, H. A. Kingravi, and P. A. Level Change, Cambridge University Press, 2021, Vela, A comparative study of efficient initializap. 1211–1362. tion methods for the k-means clustering algorithm, ESWA, 40 (2013), pp. 200–210. [19] R. Ge, M. Ester, B. J. Gao, Z. Hu, B. K. Bhattacharya, and B. Ben-Moshe, Joint clus[8] J. A. Church, N. J. White, R. Coleman, ter analysis of attribute data and relationship data: K. Lambeck, and J. X. Mitrovica, Est. of the The connected k-center problem, algorithms and reg. dist. of sea level rise over the 1950–2000 period, applications, TKDD, 2 (2008), pp. 7:1–7:35. J. Clim., 17 (2004), pp. 2609–2625. [9] V. Cohen-Addad, H. Esfandiari, V. S. Mirrokni, and S. Narayanan, Improved approximations for euclidean k-means and k-median, via nested quasi-independent sets, in Proc. STOC, 2022, pp. 1621–1628. [10] Copernicus Climate Change Service, Climate Data Store, Sea level gridded data from satellite observations for the global ocean from 1993 to present., 2018. acc. on 2025-07-01.

[20] R. Guan, Z. Li, W. Tu, J. Wang, Y. Liu, X. Li, C. Tang, and R. Feng, Contrastive multiview subspace clustering of hyperspectral images based on graph convolutional networks, TGRS, 62 (2024), pp. 1–14. [21] N. Gupta, A. Pancholi, and Y. Sabharwal, Clustering with Internal Connectedness, in WALCOM: Algorithms and Computation, N. Katoh and A. Kumar, eds., Berlin, Heidelberg, 2011, Springer, pp. 158–169.

[11] L. Drexler, J. Eube, K. Luo, D. Reineccius, [22] M. Inaba, N. Katoh, and H. Imai, Applications H. Röglin, M. Schmidt, and J. Wargalla, of weighted voronoi diagrams and randomization to Connected k-center and k-diameter clustering, Alg., variance-based k-clustering, in Proc. SoCG, ACM, 86 (2024), pp. 3425–3464. 1994, pp. 332–339.

[23] K. Kailing, H.-P. Kriegel, and P. Kröger, [36] F. Pedregosa, G. Varoquaux, A. GramDensity-connected subspace clustering for highfort, V. Michel, B. Thirion, O. Grisel, dimensional data, in Proc. SDM, SIAM, 2004, M. Blondel, P. Prettenhofer, R. Weiss, pp. 246–256. V. Dubourg, J. VanderPlas, A. Passos, D. Cournapeau, M. Brucher, M. Perrot, [24] A. Klos, J. Kusche, L. Fenoglio-Marc, M. S. and E. Duchesnay, Scikit-learn: Machine learnBos, and J. Bogusz, Introducing a vertical land ing in python, JMLR, 12 (2011), pp. 2825–2830. motion model for improving estimates of sea level rates derived from tide gauge records affected by [37] R. Rietbroek, S.-E. Brunnabend, J. Kusche, J. Schröter, and C. Dahle, Revisiting the conearthquakes, GPS Solutions, 23 (2019), p. 102. temporary sea-level budget on global and regional [25] M. R. Kramer and J. van Leeuwen, The comscales, Proc. PNAS, 113 (2016), pp. 1504–1509. plexity of wire routing and finding minimum area layouts for arbitrary VLSI circuits, Adv. Comput. [38] K. Sim, V. Gopalkrishnan, A. Zimek, and G. Cong, A survey on enhanced subspace clusRes, 2 (1984), pp. 129–146. tering, Data mining and knowledge discovery, 26 [26] H.-P. Kriegel, P. Kröger, and A. Zimek, (2013), pp. 332–397. Clustering high-dimensional data: A survey on subspace clustering, pattern-based clustering, and [39] R. E. Tarjan and J. Van Leeuwen, Worst-case analysis of set union algorithms, JACM, 31 (1984), correlation clustering, TKDD, 3 (2009), pp. 1–58. pp. 245–281. [27] S. Lattanzi and C. Sohler, A better kmeans++ algorithm via local search, in Proc. [40] P. R. Thompson and M. A. Merrifield, A unique asymmetry in the pattern of recent sea level ICML, vol. 97, 2019, pp. 3662–3671. change, Geophysical Research Letters, 41 (2014), [28] P. Le Traon, F. Nadal, and N. Ducet, An impp. 7675–7683. proved mapping method of multisatellite altimeter [41] H. Validi, A. Buchanan, and E. Lykhovyd, data, JTECH, 15 (1998), pp. 522–534. Imposing contiguity constraints in political district[29] J. Li, Agglomerative connectivity constrained clusing models, OR, 70 (2022), pp. 867–892. tering for image segmentation, Stat. Anal. Data [42] J. Wang, C. Tang, X. Zheng, X. Liu, Min., 4 (2011), pp. 84–99. W. Zhang, and E. Zhu, Graph regularized [30] Y. Li, J. Feng, X. Yang, S. Zhang, G. Chao, spatial–spectral subspace clustering for hyperspecL. Zhao, and H. Fu, Analysis of sea level variabiltral band selection, Neural Netw., 153 (2022), ity and its contributions in the bohai, yellow sea, pp. 292–302. and east china sea, Frontiers in Marine Science, 11 [43] P. Wang, H. Liu, A. M.-C. So, and (2024), p. 1381187. L. Balzano, Convergence and recovery guarantees [31] S. P. Lloyd, Least squares quantization in pcm. of the k-subspaces method for subspace clustering, Bell Laboratories Technical Memorandum, later in ICML, 2022, pp. 22884–22918. published in IEEE Trans. Inf. Theory, 1982, 1957. [44] C. You, C.-G. Li, D. P. Robinson, and R. Vidal, Oracle based active set algorithm for scal[32] S. P. Lloyd, Least squares quantization in pcm, able elastic net subspace clustering, in Proc. CVPR, IEEE Trans. Inf. Theory, 28 (1982), pp. 129–137. 2016, pp. 3928–3937. [33] J. F. Lynch, The equivalence of theorem proving and the interconnection problem, ACM Sigda [45] C. You, D. P. Robinson, and R. Vidal, Scalable sparse subspace clustering by orthogonal Newsl., 5 (1975), pp. 31–36. matching pursuit, in CVPR, 2016, pp. 3918–3927. [34] M. Mahajan, P. Nimbhorkar, and K. R. Varadarajan, The planar k-means problem is nphard, TCS, 442 (2012), pp. 13–21. [35] F. Murtagh, A survey of algorithms for contiguity-constrained clustering and related problems, Comput. J., 28 (1985), pp. 82–88.

A

Theory

A.1 Further related work Partitional k-clustering and k-means. Our research belongs to the field of partitional metric clustering where we are given a point set P from a metric space, and the goal is to partition these points into clusters while optimizing a cost function and possibly obeying side constraints. Famous partitional clustering problems are k-median, k-means or k-center. For us, kmeans is specifically relevant: Given P ⊂ Rm , find centers c1 , . . . , ck ∈ Rm to minimize the squared error when the points P are replaced by their closest center, i. e., minimize x∈P mini=1,...,k ||x − ci ||2 . This problem is NPhard, both for constant k [3] and for constant m [34] (if both are simultaneously constant then there is a polytime algorithm for the problem [22]). The best known approximation algorithm for k-means achieves a factor of 5.912 [9]. Practical algorithms for k-means. The widely known heuristic for k-means by Lloyd was developed in 1957 [31, 32]. A milestone result is due to [5] who presented the k-means++ algorithm that outperforms Lloyd’s algorithm empirically and achieves an O(log k)approximation guarantee. State-of-the-art algorithms for practical approximation improve this to an O(1)worst-case guarantee [27]. Exact and approximate algorithms for connected k-clustering [19] and [11] consider connected k-center clustering, where the objective is to minimize the maximum radius, rather than a sum of squared errors. [19] show that connected k-center can be solved optimally when the neighborhood graph is a tree. [11] provide an O(log2 k) approximation for connected k-center with general connectivity graphs. They also consider the special case of Euclidean metrics in Rm where their result for connected k-center improves to an O(m2.5 ) approximation guarantee. The connected k-means problem is a special case of the connected subspace clustering problem where m′ = 0. It is sum-based rather than the above radiusbased k-center problem, and it uses squared Euclidean distances. [21] study the connected k-means problem under the restriction that the connectivity graph is a starlike graph. They show that the problem is NP-hard to approximate better than Ω(log n) even for this restricted graph class and provide an O(log n)-approximation for this special graph class. [14] show that the related metric connected kmedian problem can be solved optimally on trees. For general connectivity graphs, they show that the problem is NP-hard to approximate better than Ω(n1−ϵ ) and they give an (n·log2 k)-approximation algorithm for this

case. There are no results about the special case of connectivity graphs that are grids. A.2 Proofs for hardness of connected (k, p)clustering problem In this section we present our theoretical result on the computational hardness of the following problem. Problem A.1. (connected (k, p)-clustering). Given a metric space (Rm , d), V ⊂ Rm , a connected graph G = (V, E) and parameters k, p ∈ N. Find a partition of V into k clusters C1 , . . . , Ck with corresponding centers c1 , . . . , ck ∈ Rm s.t. for [k], Ci induces a P i ∈P connected subgraph in G and i∈[k] x∈Ci d(x, ci )p is minimal. Note that this is a m′ -subspace k-clustering problem with m′ = 0. Concretely we will prove the following theorem. Theorem A.2. Given k ∈ N and a constant p ∈ N. The connected (k, p)-clustering problem is NP-hard to approximate within factor Ω(|V ′ |1/2−ε ), for any ε > 0 on grid graphs G′ = (V ′ , E ′ ) with holes and Euclidean t-dimensional space, where t ≥ k. To do so we first consider a discrete variant of our clustering problem that restricts the centers to originate from the set of input vertices. Problem A.3. (discrete connected (k, p)clustering). Given a connected graph G = (V, E), a distance metric d : V × V → V and parameters k, p ∈ N. Find a partition of V into k clusters C1 , . . . , Ck with corresponding centers c1 , . . . , ck ∈ V s.t. for i ∈ [k], P Ci induces P a connected subgraph in G that minimizes i∈[k] x∈Ci d(x, ci )p . The cases p = 1 and p = 2 are commonly known as k-median and k-means clustering, respectively. We show the hardness of discrete connected (k, p)clustering by reducing from the vertex disjoint paths problem in grid graphs with holes, which is known to be NP-complete [25, 17, 33]. Problem A.4. (vertex-disjoint paths problem). Given an undirected and unweighted graph G = (V, E) and k pairs of vertices (s1 , t1 ), (s2 , t2 ), . . . , (sk , tk ). Decide if there exist vertex-disjoint paths P1 , . . . , Pk s.t. Pi is an (si , ti )-path. Let G = (V, E) with (s1 , t1 ), . . . , (sk , tk ) be the input instance for the vertex-disjoint paths problem, where G is a grid and si , ti ∈ V for all 1 ≤ i ≤ k. Since G is a grid, each vertex in V can be represented by a pair of integer coordinates. For v ∈ V , we write v = (x(v), y(v)), where x : V → N and y : V → N. We

also assume that G might contain holes, i.e., G is not necessarily a complete grid. Based on G = (V, E) and the k node pairs we construct an instance for discrete connected (k, p)clustering. Let m denote some even integer parameter to be chosen later. Then the connectivity graph Gm = (V ′ , E ′ ) is defined as follows: We start with (V ′ , E ′ ) = (V, E) and then we apply the following modifications. For each horizontal edge {u, v} ∈ E (i.e., y(u) = y(v) and |x(u) − x(v)| = 1) in the grid, we check whether there exists a pair of vertices (si , ti ) with x(si ) ∈ {x(u), x(v)} or x(ti ) ∈ {x(u), x(v)}, respectively. If this is the case then we replace in E ′ the edge (u, v) by a path (in the same row y(u) = y(v)) ′ with m intermediate vertices. Let V(u,v) be the set of vertices on this path, including u and v itself. If there does not exist a node si or ti in the same column as ′ u or v then we define V(u,v) = {u, v}. Analogously, we replace each vertical edge {u, v} ∈ E ′ by a path if there exists a node si or ti in one of the corresponding rows y(u) and y(v). For i ∈ [k], we add to V ′ an (2m + 1) × (2m + 1) grid centered at si . Let Bsi be the vertex set consisting of si and the vertices of this grid. The edges on the grid Bsi are chosen such that the grid consists of 4 disjoint ((m + 1) × m) complete grids Xs1i , Xs2i , Xs3i , Xs4i and there is no edge between two grids Xsai , Xsbi with a ̸= b and each such grid is connected to si . Additionally if there exists an edge (u, si ) in E then one subgrid is connected to u. We do the same construction for ti and let Bi = Bsi ∪ BtS i. S ′ Then V ′ = ( i∈[k] Bi ) ∪ ( (u,v)∈E V(u,v) ) ∪ V . See Figure A.1 for an example of the construction. Note that by construction, if G is a grid then Gm is a grid with holes. ′ ′ We now define the S distance function d : V × V → ′ R≥0 . Let A = V \ i∈[k] Bi . We set d(x, y) = 0 with x, y ∈ A and d(x, y) = 0 with x, y ∈ Bi for i ∈ [k]. Additionally d(x, y) = 1 with x ∈ A, y ∈ Bi and finally d(x, y) = λ with λ ∈ (1, 2] and x ∈ Bi , y ∈ Bj for i ̸= j. The constructed metric d can be viewed as a star metric, where all vertices in A are located in the center and all vertices in Bi are located at a leaf with distance 1. We first argue that if there exists a solution for the vertex-disjoint paths problem in G, then there exists a solution for the discrete connected (k, p)-clustering in Gm with cost at most n + n2 m.

s1 t1

(a)

(b)

Figure A.1: Figure (a) shows a grid graph G = (V, E) and vertex pair (s1 , t1 ). Figure (b) shows the graph Gm with m = 3. The subgraphs added at s1 and t1 are indicated in color blue. Each edge (u, v) ∈ E that crosses the row or column of s1 , t1 is replaced by a path as indicated by the red dots. solution for discrete connected (k, p)-clustering on Gm that has the claimed cost. We choose the vertices s1 , . . . , sk as the k centers and assign all vertices in Bi to si . Additionally, for all edges (u, v) of path Pi we assign ′ the vertices V(u,v) to si to ensure that all points assigned to si are connected. Next we deal with vertices in V ′ that are not yet assigned. While there are unassigned vertices, we pick an arbitrary vertex v ∈ V ′ that is assigned to some cluster but has unassigned neighbors. We then assign each neighbor of v to the same cluster as v. By assigning this way, we make sure that all clusters in Gm are connected. Now we bound the cost of this solution. Every vertex in A has distance 1 to its assigned center. For all i ∈ [k] the vertices in Bi have distance 0 to their assigned center since we have a center located in each Bi . S ′ Therefore we have costs |A| ≤ |V | + (u,v)∈E V(u,v) ≤
|V | 2 2 |V | + |V | m, where we used that |E| ≤ 2 ≤ |V | /2 ′ and for all (u, v) ∈ E, |V(u,v) | ≤ m + 2 ≤ 2m. Next we argue that if there does not exist a solution for the vertex-disjoint paths problem in G, then there does not exist a solution for discrete connected (k, p)clustering in Gm with cost smaller than m2 . Lemma A.6. If there does not exist a solution for the vertex-disjoint paths problem in G, then for all p ∈ N there does not exist a solution for discrete connected (k, p)-clustering in Gm with cost smaller than m2 .

Proof. Assume there is a solution for discrete conLemma A.5. If there exist vertex-disjoint (si , ti )nected (k, p)-clustering in Gm with cost smaller than paths Pi for all i ∈ [k] in G, then for all p ∈ N there 2 m . Then there needs to be a center located in each set exists a solution for discrete connected (k, p)-clustering 2 B , since otherwise each vertex in Bi incurs distance at i in Gm with cost at most |V | + |V | m. least 1 to its assigned center and therefore the total cost Proof. To show the lemma, we construct a feasible would be at least |Bi | = 2 · (2m + 1)2 > m2 . Therefore,

each set Bi contains one center ci and the set A contains none. Additionally, at most m2 vertices of the set Bsi can be assigned to a center not equal to ci , otherwise they would incur cost at least m2 . This implies that at least |Bsi | − m2 = 3m2 + 4m + 1 vertices in Bsi are assigned to ci . Note that Bsi consists of 4 disjoint parts with size m2 + m each. Since at least 3m2 + 4m + 1 vertices in Bsi are assigned to ci it holds by the pigeonhole principle that at least one vertex from each part Xs1i , . . . Xs4i is assigned to ci . That implies that the cluster corresponding to ci contains connected subsets of all sets Xs1i , . . . , Xs4i . Since si is enclosed by Xs1i , . . . Xs4i and each Xs1i , . . . , Xs4i has at most one edge that leaves Bsi , it must be that si is assigned to ci as otherwise the cluster containing si would not be connected. Analogously we argue that ti is assigned to ci . The cluster that contains si and ti also needs to contain an (si , ti )-path in Gm to ensure that the cluster is connected. Since the clusters are disjoint, this would imply that there are vertex-disjoint (si , ti )-paths in G, which is a contradiction to the assumption.

Lemma A.8. Given metric space (Rm , d) and a connected graph G = (V, E) with V ⊂ Rm , n = |V | and k, p ∈ N. Let C1∗ , . . . , Ck∗ be an α-approximation for the connected (k, p)-clustering problem on G. Then one can compute in O(kn) time a solution C1 , . . . , Ck with centers c1 , . . . , ck ∈ V that is an (2p α)-approximation for the discrete connected (k, p)-clustering problem on G. Proof. Let ∆ be the optimal cost for the connected (k, p)-clustering problem on G and ∆′ be the optimal cost for the discrete connected (k, p)-clustering problem on G. Let c∗1 , . . . , c∗k be the corresponding centers of C1∗ , . . . , Ck∗ . For i ∈ [k], let c′i = arg minx∈Ci∗ d(x, c∗i )p . By the relaxed triangle inequality we get k X X

d(x, c′i )p ≤

i=1 x∈Ci∗

k X X

2p−1 (d(x, c∗i )p + d(c∗i , c′i )p )

i=1 x∈Ci∗

≤

k X X

2p−1 (d(x, c∗i )p + d(c∗i , x)p )

i=1 x∈Ci∗

≤ 2p α∆ ≤ 2p α∆′ .

We are now ready to prove our hardness result for discrete (k, p)-connected clustering.

Note that computing c′1 , . . . , c′k can be done in O(kn) time. Since C1∗ , . . . , Ck∗ with centers c′1 , . . . , c′k is a Theorem A.7. Given k, p ∈ N. Discrete connected feasible solution for discrete connected (k, p)-clustering, (k, p)-clustering is NP-hard to approximate within facthe proof is complete. ′ 1/2−ε ′ tor Ω(|V | ), for any ε > 0, on grid graphs G = ′ ′ (V , E ) with holes. We will again show the hardness of the connected Proof. Given a graph G = (V, E) with |V | = n and k pairs (si , ti ) as input instance for the vertex-disjoint paths problem. Let m = na for some a ≥ 2 to be chosen later and consider the graph Gm . By Lemma A.5 and Lemma A.6 there cannot be an approximation algorithm for discrete connected (k, p)clustering on grids with holes that has an approximation m2 m2 m na−2 factor smaller than n+n 2 m ≥ 2n2 m = 2n2 = 2 . for all k, p ∈ N. Notice that the number n′ = |V ′ | of nodes in the instances of connected clustering satisfies n′ ∈ O(k · m2 + m · n2 ). Since k ≤ n and m ≥ n, this implies n′ ≤ cnm2 = cn2a+1 for a sufficiently large constant c. Therefore, n ≥ (n′ /c)1/(2a+1) and the a−2 bound on the approximation factor becomes n 2 ≥

(k, p)-clustering by reducing from the vertex-disjoint paths problem. Note that we cannot directly apply the reduction used for the discrete setting, since the continuous variant gets the metric space as input. To circumvent this issue, we focus on the Euclidean metric space and tailor an instance in k-dimensions that incurs the same distances as utilized in the construction for the discrete variant, which then allows us to use the discrete variant to solve the vertex-disjoint paths problem. Theorem A.2. Given k ∈ N and a constant p ∈ N. The connected (k, p)-clustering problem is NP-hard to approximate within factor Ω(|V ′ |1/2−ε ), for any ε > 0 on grid graphs G′ = (V ′ , E ′ ) with holes and Euclidean t-dimensional space, where t ≥ k.

Proof. Let G = (V, E) with (s1 , t1 ), . . . , (sk , tk ) be 1 5 = Ω((n′ ) 2 − 4a+2 ). Now it suffices to an instance for the vertex-disjoint paths problem. Let choose a large enough so that 5/(4a + 2) ≤ ε. n = |V | and consider the graph Gm = (V ′ , E ′ ) with Next we show an auxiliary lemma to relate the m ∈ N . Let d denote the Euclidean distance. Next discrete connected (k, p)-clustering problem and the consider a k-dimensional simplex with edge length 1 and continuous variant. Effectively we get that an α- the origin as center. We place all the points in A at the approximation for the continuous variant implies a origin and all the points in Bi at its own vertex of the simplex. Then d(x, y) = 0 with x, y ∈ A and d(x, y) = 0 (2p α)-approximation for the discrete variant. with x, y ∈ Bi for i ∈ [k]. Therefore it holds d(x, y) = 1 (n′ /c)(a−2)/(2a+1) 2

√ with x ∈ A, y ∈ Bi and d(x, y) = 2 for x ∈ Bi , y ∈ Bj for i ̸= j. Assume an α-approximation for the connected (k, p)-clustering on Gm and d exists. By Lemma A.8, we can construct in polynomial time a 2p α-approximation for discrete connected (k, p)-clustering on Gm and d′ , where d′ is d restricted to the vertices of Gm . Since 2p is constant and through Theorem A.7 it must be that α ∈ Ω(|V ′ |1/2−ε ). B

Application

B.1 Dataset There are two principal sources of satellite-derived SLA products: The [10], which provides climate-oriented products emphasizing long-term homogeneity, and the Copernicus Marine Service, which prioritizes accuracy regarding the sea level estimation at each time step and enhanced spatial sampling. The CMEMS dataset [13] builds upon the steady twosatellite merged constellation formed by the reference missions (TOPEX-Poseidon, Jason-1/2/3, Sentinel-6A) (also used in C3S products) by homogenizing five additional missions with respect to this reference set. Although the constellation evolves over time, the integration of more satellites improves spatial coverage and temporal accuracy, potentially at the expense of longterm stability (see Quality Information Document [13]). Given our focus on the identification and analysis of regional spatiotemporal patterns, the higher-resolution, time-accurate CMEMS product is preferred. Spatial interpolation of the altimetry data onto the uniform grid is performed by the CMEMS using the Optimal Interpolation method described by [28].

dissimilarity measure and carefully craft a distance function that takes into account both the location of a point and the time series associated with it. For Xi = (xi , ti ) and Xj = (xj , tj ) which contain both the spatial and temporal information, they use   d(xi , xj ) r(ti , tk ) D(Xi , Xj ) = 1 − exp − c0 where d(xi , xj ) describes the spatial distance between two points xi and xj , r(ti , tj ) is the correlation coefficient between the time series ti and tj , and c0 ≈ 4328 is a normalization constant. We use the Haversine distance as the spatial distance which computes the distance between two points x = (λ1 , ϕ1 ) and y = (λ2 , ϕ2 ) in the simplified model that the earth is a perfect sphere with radius R of roughly 6.370 km. The Haversine distance is defined as
√ √ R(x, y) = R · 2 · arctan 2 a, 1 − a , where a = sin2



ϕ2 − ϕ 1 2



+ cos(ϕ1 ) cos(ϕ2 ) sin2



λ2 − λ 1 2

 .

For the correlation between two time series t1 and t2 , we employ the Pearson correlation coefficient. This measure yields a value in the range [−1, 1] where 1 denotes a perfect positive linear correlation, 0 indicates no linear correlation, and –1 signifies a perfect negative linear correlation. The exponential function produces values in the interval [0, 1], assigning higher values to closer grid points and lower values to those that are B.2 Reimplementation of Thompson and more distant. Consequently, two time series that exhibit a strong Merrifield [40] develop a specific distance function for clusters that they use with an average linkage negative correlation but are geographically close are algorithm without explicit enforcement of connectivity. considered more dissimilar than two equally negatively We implemented their method following [16], resulting correlated time series that are spatially distant. In cases where two time series are uncorrelated, their resulting in Algorithm B.1. distance is 1, regardless of their spatial separation. Algorithm B.1 Hierarchical Agglomerative Clustering Overall, the greater the similarity between two time (Average Linkage), no explicit connectivity constraint series and the closer their associated grid points are, the smaller the resulting distance, approaching zero in 1: Require: P = {x1 , . . . , xn }, number of clusters k the case of high correlation and close proximity. 2: Ensure: A clustering of P into k clusters Algorithm B.1 does not enforce connectivity. Con3: Create H = {Ci = {xi } | xi ∈ P } nected clusters only emerge if the underlying data 4: Compute D(Ci , Cj ) = D(xi , xj ) for all xi , xj ∈ P matches certain criteria, such as smoothness in time and 5: while |H| > k do space. If the input data is not appropriately filtered, the ∗ 6: Add C = Ci ∪ Cj to H, remove Ci and Cj resulting segmentation lacks spatial coherence and fails 7: for all remaining clusters Cl ∈ H do to capture meaningful structures. In such cases, local ∗ 8: Compute D(C , Cl ) extrema are typically isolated into their own clusters. Figure 4.2 shows an example, and Table D.4 shows the This algorithm does not enforce connectivity di- decline in performance of the algorithm for unfiltered rectly. Thompson and Merrifield use average linkage as

Figure A.2: A point set and three one-dimensional subspaces in R2 . For each point, the distance to its closest subspace (symbolized by gray lines) is computed and squared. The sum of these squared distances is the objective that shall be minimized by the choice of the subspaces. data. B.3 EOFs and PCs for global sea level data As can be seen in Figure B.1, the first three EOFs and their corresponding PCs capture the long-term trend, the seasonal cycle, and the El Niño–Southern Oscillation (ENSO) signal, respectively. The ENSO signal is partially mixed with other modes of variability, which is a known issue in EOF analysis. B.4 Implications for sea level research Combining subspace clustering with Empirical Orthogonal Function (EOF) decomposition offers a practical way to extract regionally coherent patterns in sea level variability that may not be clearly represented in a standard global EOF analysis. Conn-Subspace can identify interesting regions with a specific sea level rise behavior. As an example for this, we discuss the dominant EOFs for the global and for two specific interesting regions, to showcase how two clusters reflect a known pattern of sea level behavior. As expected, globally, the most dominant sea level signals are the trend (see Figure B.1(a),(b)) and the seasonal cycle which is driven by solar forcing and ocean heat uptake, and land-ocean fluxes such as seasonally varying freshwater flow from ice sheets and glaciers, the world’s rivers and ocean precipitation and evapo-

ration (see Figure B.1(c),(d) in the appendix). The associated EOFs 0 and 1 shows a trend pattern and clear positive/negative separation between the northern and southern hemispheres, with the corresponding PCs showing a clear trend and annual cycle. By applying subspace clustering, the analysis can operate on more localized feature spaces, potentially revealing regional modes that are less prominent in a global context. For the subspace clustering approach, in most subregions, such as the Atlantic Ocean Figure B.2(a), the trend and annual cycle are also the dominating signals. However, in the tropical Pacific Ocean the dominant signal is driven by the so called El Niño–Southern Oscillation effect (ENSO), a climate phenomenon driven by varying winds and sea surface temperatures in the tropical Pacific that impacts climate and weather globally, e.g. leading to extensive rain in South America while causing droughts in Australia. In the global EOF decomposition, this effect is still one of the most dominant signals (see Figure B.1(e),(f)) indicating its impacts all over the globe. The subspace clustering approach is able to derive clusters, which fit well to the eastern and western Pacific regions responsible for driving the ENSO phenomenon (see Figure B.2(b), (c), see Cluster 6 and 9 in Figure 4.5). In these clusters, the dominant signal aligns

(a) EOF 0

(b) PC 0

(c) EOF 1

(d) PC 1

(e) EOF 2

(f) PC 2

Figure B.1: Empirical Orthogonal Functions (EOFs) and their corresponding Principal Components (PCs) derived via Singular Value Decomposition (SVD) of global sea level data. Panels (a), (c), and (e) show the spatial patterns (EOFs 0–2), while panels (b), (d), and (f) present their respective temporal variations (PCs 0–2). EOFs 0 and 1, along with their PCs, primarily capture the long-term trend and the seasonal cycle, respectively. The El Niño–Southern Oscillation (ENSO) signal is evident in EOF 2, although it is partially mixed with other modes of variability. closely with the ENSO signal, while the annual signal is found in the second EOF associated with these regions (see Figure B.3). Finding clear regional separations based on signal content using the subspace clustering approach enables a much more in-depth data analysis on regional scales, focusing on the dominant signals and drivers of sea level change relevant for these regions. In geodesy, EOFs are also used to reconstruct incomplete tide gauge data. To achieve this, satellite altimetry data are decomposed into EOFs and corresponding temporal components. The longer signals of tide gauge records are then fitted to the data using a least squares approach, see for example [8]. This is usually done for the entire globe. A clustering of the sea enables employing this technique regionally, e.g. for regions of good tide gauge record density. To the best of

our knowledge this has seldom been applied (e.g. by [24]) for the eastern Pacific cluster and we plan to explore this in future work. However, the regional separation may lead to inconsistencies regarding signal content at the borders to neighboring regions, which warrants further investigations. The separation of the ocean into contiguous regions additionally enables a regional reconstruction of historic sea levels prior to the altimetry era. C

Implementation

C.1 Accompanying pseudocode The following code was skipped in the main part of the paper. Notice that in Step 12 in the average linkage version of Algorithm C.1, a new computation of distances is only necessary for the clusters in the new neighborhood. For

(a) Signal from cluster 8 (red) together with the ENSO index (blue).

(a) Second EOF for cluster 6.

(b) Second EOF for cluster 9. (b) Signal from cluster 6 (red line) together with the ENSO index (blue line).

Figure B.3: Second EOFs for clusters 6 and 9. The second EOFs for both clusters show a strong annual cycle, which is the second most dominant signal in these regions.

clusters. The procedure is given in Algorithm 3.2. Notice that this method inherently uses that the time series and the locations are somewhat correlated because oth(c) Signal from cluster 9 (red) together with the erwise we would not start with small contiguous regions. ENSO index (blue). Step 3 constructs a graph where two points are Figure B.2: First EOFs for clusters 6, 8 and 9 together connected if they are neighbors in N and are assigned to with the ENSO index. The first EOFs for clusters 6 and the same cluster in the input clustering. Visually, when 9 show a strong correlation with the ENSO index, while looking at a clustering like in Figure 3.1(a), points are connected in GC if they have the same color and are the first EOF for cluster 8 shows a weak correlation. adjacent. (Notice that the set of edges E ′ is a subset of the edges E in N ). The connected components of those, it has to be recomputed (only) if a cluster was GC are the contiguous single-color regions visible in neighbor to only one of Ci , Cj since in this case the distances to the other cluster have not been computed Algorithm C.1 Agglomerative Connected Clustering yet. If a cluster C ′ was neighbor to both, then we use a ′ ′ Require: P = {x1 , . . . , xn }, graph N , k ∈ N weighted sum of d(C , Ci ) and d(C , Cj ) instead. Ensure: A clustering of P into k connected clusters C.2 Details on connectivity methods In this 1: Create H = {Ci = {xi } | xi ∈ P } section, we expand on the algorithms discussed in 2: Compute D(Ci , Cj ) = D(xi , xj ) for all xi , xj ∈ P Subsection 3.3. 3: for each cluster Ci ∈ H do 4: Set Ni ← {Cj | xi and xj are adjacent in N } C.2.1 Establishing connectivity via merg5: while |H| > k do ing The first approach converts a clustering with pos6: Identify pair (Ci , Cj ) with minimal distance sibly disconnected clusters into a connected clustering. 7: Add C ∗ = Ci ∪ Cj to H, remove Ci and Cj The idea of this algorithm is that we typically expect 8: N ∗ ← Ni ∪ Nj that, while the clusters may be disconnected, we still 9: for each Cl ∈ N ∗ do expect medium size regions that belong to the same 10: Remove Ci , Cj from Nl and add C ∗ to Nl cluster. For an example of such a situation, see Fig- 11: for each C ∈ N ∗ do j ure 3.1(a). Our aim is to iteratively remove small re- 12: Compute distance D(Ci , Cj ) gions and merge them into neighboring regions of other

Algorithm C.2 Greedy k-means++ Clustering Require: Point set P , k ∈ N,candidate parameterℓ ∈ N Ensure: centers c1 , . . . , ck and clusters C1 , . . . , Ck 1: Pick initial center c1 uniformly at random from P 2: for i = 2 to k do 3: Sample ℓ candidate points x1 , . . . , xℓ ∈ P , where each x ∈ P is chosen with probability p(x) = P 4:

min1≤j<i ∥x − cj ∥2 2 y∈P min1≤j<i ∥y − cj ∥

For each candidate xm , compute   X ϕ(xm ) = min min ∥y − cj ∥2 , ∥y − xm ∥2 y∈P

1≤j<i

number of points in Amin . Then Amin is merged into that cluster. This decreases the number of regions by one. The process stops when the number of regions reaches k. Each region now becomes a cluster. During Step 13 we make sure that if an input cluster is split by this process, one of the clusters retains the original cluster number, and if two connected components belonging to different clusters are merged, the new cluster assumes the place of one of those it was created from (mostly for picture consistency and evaluation). C.2.2 Proof for Theorem 3.1 Theorem 3.1 states that the iterative merging algorithm presented in Alg. 3.2 can be implemented in O(|E| + nk + n log n + nα(n)) time. Here we provide the accompanying proof.

Theorem 3.1. Algorithm 3.2 can be implemented in O(|E| + n(k + log n)) time. Consequently, its run5: Set ci to be xm that minimizes ϕ(xm ) ning time is O(n(k + log n)) on bounded-degree graphs, 6: repeat including the grid graphs used in our application, and 7: Assign each x ∈ P to arg min1≤i≤k ∥x − ci ∥2 and O(nk) whenever k = Ω(log n). let C1 , . . . , Ck be the resulting clusters Proof. Given the neighborhood graph N = (V, E) 8: for i = 1 to k do P 1 with |V | = n and |E| = m, a clustering C1 , . . . , Ck and 9: Compute µ(Ci ) = |Ci | x∈Ci x the subspaces associated with each cluster S1 , . . . , Sk . 10: Update center: ci ← µ(Ci ) Let z0 be the initial number of connected components 11: until convergence criterion is met and zi = |Z| in iteration i of the while-loop. For a component A ∈ Z, let n(A) be the clusters associated with the connected components adjacent to A in the Figure 3.1(a). The connected components are computed component graph GCC = (Z, E ′′ ). So, this describes in Step 4 and stored in Z. Notice that every A ∈ Z clusters a connected component can possibly be merged satisfies A ⊂ Ci for some i ∈ [k], and that A would be into. For easy lookup these neighboring clusters are a feasible connected cluster in N . However, there are stored in a bitset. We assume that for each A ∈ Z, more than k components in Z because the input clusters n(A) > 0, so each component has at least one neighbor were not connected. it can be merged into. Step 5 computes GCC which uses the connected Preprocessing First, we need to build the graph components as nodes and connects two regions A and B GC by scanning E and only keeping edges u, v, if if they contain points a ∈ A and b ∈ B that are adjacent u ∈ Ci and v ∈ Ci , for any Ci . Next, computing the in N , i. e., {a, b} ∈ E. Note that regions that belong to connected components of GC is O(n + m) by BFS/DFS. the same cluster are never connected in GCC since they Constructing GCC by scanning E and adding an edge would be one component then. We use GCC to reduce between components that have at least one adjacent pair the number of disconnected regions by merging regions of nodes in N takes O(m). into larger ones, beginning with the one of smallest Winner computation and merging In order to incardinality. For each iteration of the loop in Step 6 to crease efficiency, we can precompute a numerical scorStep 12 we compute the region with the least number of ing vector for each connected component A ∈ Z that points and then merge it into a neighboring cluster. In contains the aggregated preferences of all x ∈ A over order to choose the cluster to merge into, we compute the subspaces Si for i ∈ [k]. This is done as follows; the preference of each point x ∈ Amin in Step 9: We we precompute for each point x ∈ V a ranking of the compute the cluster Ci where the associated subspace clusters πx (1), . . . , πx (k) by increasing distance to the Si is closest to x, but only consider those Ci which have subspaces Si . Let rx [i] denote the position of cluster i a subregion B neighboring Amin . So from the clusters for point x. Fix any non-decreasing positional scoring surrounding Amin , we compute which cluster is preferred function s on ranks. The contribution of x to its comby x. Then we do a majority vote in Step 10 and choose ponents valuation of cluster [i] is s(rx [i]), which is fixed the neighboring cluster, which is preferred by the largest throughout the algorithm.

For each connected component A ∈ Z aggregate these contributions in a score vector: scoreA [i] = P s(r x [i]) for all i ∈ [k]. Given the adjacent-cluster x∈A set n(A) ⊆ [k], the chosen cluster for merging is then iw (A) = arg maxi∈n(A) scoreA [i]. This scoring rule maps the majority decision well enough for this algorithm, as updating ordinal preferences would require more time and space with limited benefits. When two components A and B merge, the new component’s score vector is updated by component-wise addition scoreA∪B = scoreA + scoreB . This leads to winner lookup being in O(|n(A)|) ⊆ O(k) and score updates are O(k). Computing all distances takes time and space O(n · k). Additional to these metadata updates, in each merging step we need to join all connected components w B i that belong to the winning cluster iw and are adjacent to the cluster that is being merged A. This information is stored in each connected component. All unions of components over the complete execution cost O((z0 − k)α(z0 )) ⊆ O(nα(n)) using union-by-rank and path compression [39]. Furthermore, we need to merge the set of neighbors n(A ∪ B) = n(A) ∪ n(B) \ {cA , cB }, where cA is the cluster that component A belongs to and cB is the one component B belongs to and adjacent components belonging to each cluster. Maintaining Amin In step 7 Amin , the component of smallest cardinality, needs to be identified. This can either be done via naive sorting in O(zi log(zi )) or improved to O(log(zi )) by using a min-heap. Note that zi decreases by at least one in each iteration, yielding O((z0 −k) log(z0 )) across the O(z0 −k) iterations, which is in O(n log(n)), since z0 ≤ n. Iterations The number of times this loop is executed depends on z0 , the initial size of Z, since each iteration decreases the number of connected components by at least one. Thus, the loop runs at most z0 − k times. Conclusion The total time is O(nk) for distance and score precomputation, O(n log n) for maintaining the heap, O(nk) for winner scans and score additions, O(nα(n)) for union–find, and O(m) for scanning the edges of the neighborhood graph. Hence the claimed bound O(m + nk + n log n + nα(n)) holds. For bounded-degree neighborhood graphs, including the grid graphs used in our application, m ∈ O(n) and this becomes O(n(k + log n)); whenever k = Ω(log n), it simplifies further to O(nk). Remark We implemented a simpler version of this algorithm without the use of the sophisticated data structures described above. Our experiments show that this variant is empirically fast enough, however, if one needs a very fast version of the merging routine, the

modifications described in Proof 6 can be used. C.2.3 Filtering For this approach, we implement Step 11 of Algorithm 3.1 as follows: We treat the grid-like graph N as an image, where each vertex corresponds to a pixel. Pixels that represent land are masked using a binary mask and ignored during the filtering process. h and for any point More specifically, set σ = 1.178 x, let B(x) denote the ball around x containing all grid points that satisfy dH (x, y) ≤ 3σ where dH (x, y) is the Haversine distance between x and y. Then assign a weight to every point y ∈ B(x) using a Gaussian kernel:   dH (x, y)2 , w(x, y) = exp − 2σ 2 Notice that the definition of σ ensures that the weight drops to half its maximum at the given half-width h. These weights are then normalized: wx,y norm wx,y =P′ . The cluster that x shall bey ∈Bx ∪{x} wx,y ′ long to is determined by weighting the clusters that the points in B(x) belong to. This is done by defining the following weighted support for a cluster Ci ∈ {C1 , . . . , Ck } at point x: X wnorm (x, y). S(x, Ci ) = y∈B(x) y∈Ci

The new cluster label assigned to x is then the one that has the highest weighted support: C ∗ (x) = arg max S(x, Ci ). Ci ∈C

As is typical for Gaussian filtering this reduces noise and creates more cohesive regions. However, this does not necessary result in exactly k connected clusters, and thus Step 15 is also executed. D

Experiments

D.1 Data Preprocessing Since the Earth is spherical and the grid is geodetic (latitude/longitude), one can not calculate the Gaussian kernel once and apply it to each data point. To our knowledge, existing tools do not provide an efficient, curvature-aware spatial smoothing pipeline; we therefore implement a parallelized spherical Gaussian filter that respects the Earth’s curvature. Given a structured spatial grid N , as defined in Subsection 4.1.1, over the Earth’s surface, where each x ∈ N corresponds to a geodetic coordinate pair (ϕx , λx ) in latitude and longitude, we define a smoothing operator based on a spherical Gaussian kernel. Let S(t, x) represent the sea level anomaly data at time t and spatial location x. We construct a binary

(a) Clustering using EnSC.

(b) Clustering using SSC-OMP.

Figure C.1: Examples for the clusterings produced by SSC-OMP [44] and EnSC [45]. validity mask V over N by identifying all x ∈ N for result is a spatially smoothed dataset, which has been which S(t, x) contains no invalid (missing) observations calculated with respecting the Earth’s curvature and across all t. For each valid x ∈ N, x ∈ / V , we define its with robustness to missing values. spherical neighborhood D.2 Subspace Clustering In addition to the experiments presented in the main text, we conducted Bx = {y ∈ N \ {x} | dH (x, y) < rc } further analyses to compare our method to other subwhere dH is the great-circle distance computed using space clustering methods and hyperspectral image segthe Haversine formula and rc = 3σ is the cutoff radius. mentation algorithms. Please note that the methods The Gaussian kernel is parameterized by a user-specified discussed here were selected to be representatives of h half-width h ∈ R, where σ = 1.178 . Each neighbor their respective fields, we performed additional experiments using other algorithms and chose these for their (x,y) ) and y ∈ Bx is assigned a weight wx,y = exp(− dH2σ 2 stability, feasibility and solution quality. the central point x has a weight wx,x = 1. All weights To compare our connected subspace clustering apw x,y norm . The filtered are normalized wx,y = P ′ y ∈Bx ∪{x} wx,y ′ proach to state-of-the-art subspace clustering (withtime series at location x is computed as out explicit connectivity constraints) we ran experiX ments using the code accompanying [45], who propose a norm e x) = S(t, wx,y · S(t, y) method called sparse subspace clustering with orthogoy∈B∪{x} nal matching pursuit (SSC-OMP), and [44], who devised where NaN-values are excluded to avoid an influence a method called elastic net subspace clustering (EnSC). These experiments were run for k = 8, 15, 20 and of the missing values from the weighted average. This 25 on unfiltered and filtered data. However, it is not operation can be performed independently for each time possible to directly specify the subspace dimensions step t and each point x. For scalability reasons the ′ m , so we are not able to include those in our results. filter is executed in parallel using multiprocessing. The We tried different configurations of the parameters and Table C.1: Number of connected components for SSC-OMP and EnSC on filtered and unfiltered input data, grouped by number of clusters. The second column indicates how many clusters should be found by the algorithm and columns 3 and 4 show how many connected components the resulting clustering has. method SSC-OMP EnSC SSC-OMP EnSC SSC-OMP EnSC SSC-OMP EnSC

# clusters wanted 8 8 15 15 20 20 25 25

filtered 225 128 796 169 1374 222 1966 219

unfiltered 52 31 90 54 114 102 137 132

(a) Clustering using EKGCSC.

(b) Clustering using EGCSC.

Figure D.1: Example of a clustering resulting from two HSI methods. ended up with γ = 50 and otherwise used the default parameters as this produced the best results. As both of these methods do not enforce connectivity, the resulting clustering is fragmented, especially when using SSCOMP (see Figure C.1). On filtered data, the clusters formed by EnSC tend to be formed by multiple coherent connected components, while SSC-OMP consistently produces less coherent clusters. This is not surprising as SSC-OMP does not consider connectivity, while it is specifically integrated into EnSC. However, even EnSC does not enforce strict connectivity, but encourages the formation of more coherent clusters by combining sparse subspace clustering with least squares regression. The number of connected components of the subgraph induced by the clusters are reported in Table C.1. The subspace clustering cost of both methods are reported in Subsection D.4.

mentation, as one explicit goal here is to produce spatially coherent clusters. We performed experiments using the code provided by [6] in their paper on Graph Convolutional Subspace Clustering (GraphConvSC). To ensure a fair comparison we modified the algorithm such that it handles spherical and incomplete input data. We ran these experiments for both efficient graph convolutional subspace clustering (EGCSC) and efficient kernel graph convolutional subspace clustering (EKGCSC), which are the two methods proposed by [6]. The parameters were set as can be seen in Table D.1. We chose these parameters as they produced the best results. It is important to note that these methods are not designed to explicitly enforce connectivity. This can be seen in Table D.2, where the number of connected components that are actually built by the algorithm proposed by [6] are listed. We additionally report D.3 Hyperspectral Imaging In addition to the summed distances of the points in a cluster to the the more traditional subspace clustering methods, we the subspace spanned by that cluster, i.e. the subspace also investigated the field of hyperspectral image seg- clustering cost in the tables shown in Subsection D.4. Table D.1: Parameters configured for GraphConvSC Category Patch Clustering parameters Clustering parameters PCA EGCSC EGCSC EGCSC EGCSC EGCSC EKGCSC EKGCSC EKGCSC EKGCSC EKGCSC EKGCSC

Parameter PATCH_SIZE N_CLASSES_LIST NB_COMPS_LIST n_components regu_coef n_neighbors ro n_clusters save_affinity regu_coef n_neighbors ro gamma n_clusters save_affinity

Value 15 {8, 15, 20, 25} {5, 10, 15, 30} nb_comps from NB_COMPS_LIST 1 × 102 15 0.8 n_classes from N_CLASSES_LIST False 1 × 103 15 0.8 0.5 n_classes from N_CLASSES_LIST False

Table D.2: Number of connected components for GraphConvSC on filtered input data, grouped by number of clusters and subspace dimension m′ . The second column indicates how many clusters should be found by the algorithm and columns 3-6 show how many connected components the resulting clustering has. (a) Filtered input data

method EKGCSC EGCSC EKGCSC EGCSC EKGCSC EGCSC EKGCSC EGCSC

# clusters wanted 8 8 15 15 20 20 25 25

m′ = 5 28 29 34 37 48 39 59 52

m′ = 10 24 27 37 36 45 43 53 52

m′ = 15 24 25 39 37 45 42 49 48

m′ = 30 25 34 38 37 44 43 50 51

m′ = 15 21 24 28 34 35 33 40 53

m′ = 30 22 18 29 31 34 33 41 39

(b) Unfiltered input data

method EKGCSC EGCSC EKGCSC EGCSC EKGCSC EGCSC EKGCSC EGCSC

# clusters wanted 8 8 15 15 20 20 25 25

m′ = 5 27 24 33 30 31 43 47 61

Figures showing examples of the resulting clusterings are given in Figure D.1. Interestingly, EKGCSC does produce less coherent clusters that EGCSC. D.4 Results In this section we present the results from our set of experiments. Each table contains the subspace clustering cost of the initial clustering methods described in Subsection 3.2 in the corresponding column. The remaining columns report how the subspace clustering cost changes when using Alg. 3.1 with the corresponding connectivity methods as described in Subsection 3.3. Additionally, there are comparisons to the subspace clustering and HSI methods discussed in Subsection D.2. Those are the last four rows in each table. As they are meant as standalone comparison methods, we did not apply the merging procedures to them. Please note that the methods SSC-OMP and EnSC do not allow for a specific m′ to be provided, which is why there is only one value in those rows.

m′ = 10 27 21 33 36 33 48 38 54

Table D.3: Distances from points to subspaces using filtered input data, with 8 clusters. Blue indicates the best value in each row. m′

Agglo-ST

Conn-Agglo-Euc

Conn-Agglo-ST

Conn-KMeans++

Conn-Ward

GraphConvSC EKGCSC

GraphConvSC EGCSC SSC-OMP EnSC

5 10 15 30 5 10 15 30 5 10 15 30 5 10 15 30 5 10 15 30 5 10 15 30 5 10 15 30

initial clustering 7692.49489 6711.79823 6058.78586 4739.49195 7577.79424 6552.16401 5892.12265 4529.30297 7665.13937 6626.54510 5925.35443 4543.34794 7672.81813 6684.77629 6019.41239 4657.46825 7818.73543 6884.08856 6250.65279 4952.76399 13124.82028 13088.99551 13084.56918 13115.11044 13128.60392 13091.91857 13250.86692 13136.15622 15276.42554 14088.99291

7465.55920 6577.24879 5885.93997 4597.55834 7516.26784 6533.40635 5880.22969 4524.68553 7480.37794 6614.75710 5891.04896 4502.91717 7618.64891 6586.35707 5930.38415 4572.23265 7547.78283 6625.74611 5948.23946 4773.37472

Smooth Merge 7478.41837 6567.92509 5916.59557 4602.77377 7521.30797 6457.07786 5867.63781 4506.63351 7563.32227 6517.96774 5712.20077 4589.32804 7536.54157 6517.60083 5903.71485 4574.05150 7519.00091 6605.33257 5941.08003 4717.35659

Integrated Conn 7606.71637 6496.16303 5678.59392 4300.00575 7472.77518 6471.56399 5793.80565 4337.03429 7553.95301 6522.68229 5841.31301 4335.89561 7412.93632 6530.95691 5725.26091 4360.19843 7662.99054 6721.36554 6059.41774 4770.03843

-

-

-

-

-

-

-

-

-

-

-

-

IterMerge

PostMerge

7449.85439 6416.84714 5681.50449 4275.62553 7407.11306 6434.20938 5656.52604 4310.62839 7471.02903 6401.78462 5716.76295 4328.07089 7456.54788 6404.43239 5741.98456 4367.77138 7559.41281 6390.55542 5728.35510 4561.76944

Table D.4: Distances from points to subspaces using unfiltered input data, with 8 clusters. Blue indicates the best value in each row. Note that the algorithm Agglo-ST performs poorly on unfiltered data (see Figure 4.2), resulting in clusters unsuitable for Conn-Subspace. This is because one large cluster dominates while the others are too small to compute subspaces for m′ dimensions. This is why the corresponding values are missing in the table. m′

Agglo-ST

Conn-Agglo-Euc

Conn-Agglo-ST

Conn-KMeans++

Conn-Ward

GraphConvSC EKGCSC

GraphConvSC EGCSC SSC-OMP EnSC

5 10 15 30 5 10 15 30 5 10 15 30 5 10 15 30 5 10 15 30 5 10 15 30 5 10 15 30

initial clustering 8357.61138 7465.22683 6934.41398 5827.25554 7731.16245 6813.48844 6189.81632 4902.32309 7587.75423 6517.61474 5821.57802 4433.68417 7787.54983 6853.38225 6208.53329 4908.12907 7731.16195 6813.49009 6189.85300 4902.18010 13529.17649 13531.45765 13720.32976 13721.31299 13534.20219 13730.57769 13717.10833 13890.14598 15122.25802 14900.37739

8269.62100 7268.01158 7613.92384 6732.86394 5974.47518 4755.52725 7616.73677 6562.13264 5802.8574 4449.47911 7576.29373 6519.65702 5861.48456 4664.09250 7534.85592 6577.85847 5997.38423 4722.56464

Smooth Merge 7581.09045 6526.79913 5977.69428 4637.71137 7554.43019 6421.3932 5873.59309 4454.67732 7580.04322 6710.61413 5941.10418 4646.05051 7569.99355 6532.33852 5991.65580 4641.04607

Integrated Conn 8363.89983 7461.79679 7607.04228 6703.52087 6001.12751 4719.07835 7477.76603 6438.73443 5739.71282 4361.23423 7511.55950 6434.85064 5744.50432 4669.39230 7607.04291 6703.51083 6000.58657 4719.93802

-

-

-

-

-

-

-

-

-

-

-

-

IterMerge

PostMerge

7544.00497 6589.54375 5877.08022 4494.57215 7459.38257 6413.93545 5680.3945 4389.54962 7402.50782 6443.82804 5756.01577 4314.02620 7544.00447 6609.47078 5866.06992 4502.06469

Table D.5: Distances from points to subspaces using filtered input data, with 15 clusters. Blue indicates the best value in each row. m′

Agglo-ST

Conn-Agglo-Euc

Conn-Agglo-ST

Conn-KMeans++

Conn-Ward

GraphConvSC EKGCSC

GraphConvSC EGCSC SSC-OMP EnSC

5 10 15 30 5 10 15 30 5 10 15 30 5 10 15 30 5 10 15 30 5 10 15 30 5 10 15 30

initial clustering 7401.75317 6373.87940 5692.49488 4324.59394 7114.02254 6072.58344 5370.74695 3956.67356 7165.77020 6043.90989 5286.28692 3830.37398 7252.59653 6185.78104 5460.38901 4035.94281 7220.75107 6267.09090 5616.07627 4308.16897 12816.98138 12696.11403 12629.34658 12607.53272 12595.44702 12671.02373 12715.13632 12833.03298 15056.41022 13462.20711

7017.48275 5936.91435 5206.89874 3940.42473 7123.69161 6099.05409 5249.62310 3831.89306 7168.83511 5971.14997 5007.97436 3829.96746 7102.93392 5955.52421 5210.23140 3765.80204 7094.40494 6066.03266 5458.88736 4210.32269

Smooth Merge 7000.30368 6085.79398 5318.16172 3949.27144 7015.46846 6014.41582 5168.60667 3816.64460 7124.80622 5908.16445 5061.30361 3678.17079 7120.01640 5835.30837 4974.27926 3617.59620 7017.41398 6039.21859 5142.28195 3879.76737

Integrated Conn 7153.61319 6261.35012 5538.33864 4027.25431 6928.5687 5828.97291 5141.67013 3724.09671 7004.66985 5825.62817 5053.81309 3708.89951 6987.27644 5815.82064 5128.15819 3668.58070 7018.86197 6036.39050 5372.64187 4000.63856

-

-

-

-

-

-

-

-

-

-

-

-

IterMerge

PostMerge

6928.48755 5866.01001 5187.14312 3737.01157 6850.53328 5847.53931 5052.96718 3655.53732 6900.80387 5738.33263 4995.24785 3633.49710 6935.95820 5794.88755 4994.50226 3659.64386 6989.61962 5934.10865 5080.41592 3905.82830

Table D.6: Distances from points to subspaces using unfiltered input data, with 15 clusters. Blue indicates the best value in each row. m′

Agglo-ST

Conn-Agglo-Euc

Conn-Agglo-ST

Conn-KMeans++

Conn-Ward

GraphConvSC EKGCSC

GraphConvSC EGCSC SSC-OMP EnSC

5 10 15 30 5 10 15 30 5 10 15 30 5 10 15 30 5 10 15 30 5 10 15 30 5 10 15 30

initial clustering 8138.18456 7204.38252 6658.85384 5568.49755 7300.91007 6357.32784 5737.88249 4460.20328 7119.17181 6005.0503 5268.05087 3823.02996 7237.21547 6140.85496 5410.36757 3960.26859 7300.91029 6357.33007 5737.91865 4460.12404 13357.03370 13437.32476 13679.37267 13684.34093 13480.88323 13266.20388 13455.86709 13666.36885 14875.43995 14671.32921

7496.84981 6598.76775 6296.71185 4973.61062 7180.58798 6187.01996 5301.06031 4005.64651 7029.95831 5883.24081 5099.97031 3792.69081 7091.07433 5993.05102 5298.03961 3944.39825 7194.19093 6177.00189 5353.26835 4001.61043

Smooth Merge 7818.48633 6779.11463 6081.50442 5276.09128 7094.45986 6074.20510 5448.38671 4225.21069 7124.18824 5977.09426 5085.76466 3616.81066 6974.68506 6030.58822 5288.48289 3956.10770 7074.05179 5944.17437 5542.68503 4223.34628

Integrated Conn 7964.90935 6955.51336 6509.29316 5469.03702 7152.81688 6191.78486 5560.49203 4274.54517 6939.76197 5838.78528 5116.52196 3736.21723 7037.43763 6007.85102 5291.30412 3817.96416 7152.80277 6188.96634 5567.75397 4059.74273

-

-

-

-

-

-

-

-

-

-

-

-

IterMerge

PostMerge

7053.37876 6334.25977 5778.67113 4440.49430 7116.94334 5987.91657 5274.70075 3862.28418 6918.34345 5778.90979 5104.49715 3643.09308 6961.75097 5821.98319 5061.52510 3739.47498 7113.16163 6045.03966 5321.37167 3793.78203

Table D.7: Distances from points to subspaces using filtered input data, with 20 clusters. Blue indicates the best value in each row. m′

Agglo-ST

Conn-Agglo-Euc

Conn-Agglo-ST

Conn-KMeans++

Conn-Ward

GraphConvSC EKGCSC

GraphConvSC EGCSC SSC-OMP EnSC

5 10 15 30 5 10 15 30 5 10 15 30 5 10 15 30 5 10 15 30 5 10 15 30 5 10 15 30

initial clustering 7127.56897 6101.62757 5406.82642 4012.58077 6865.53475 5832.31292 5126.73700 3713.09606 6900.30819 5794.35416 5049.92481 3595.09951 7064.48992 6027.69854 5329.54776 3932.44561 6911.97254 5870.80942 5154.80347 3753.07503 12585.39110 12506.02452 12733.70084 12267.75023 12215.19311 12168.04789 12060.37141 12216.14648 14910.37771 13226.72225

6810.31878 5769.48728 4912.01064 3370.35309 6750.88909 5865.71361 4919.41237 3579.53496 6927.05043 5687.2591 4768.72953 3613.7709 6935.84391 5908.07834 5119.8224 3857.21766 6851.87697 5705.9256 4893.75679 3624.77683

Smooth Merge 6819.93585 5670.97155 4889.30909 3591.39863 6831.68193 5729.15024 4946.12142 3490.69984 6786.10306 5691.85431 4764.92476 3537.25328 6826.91025 5676.37115 5046.53568 3579.53919 6816.23469 5634.00779 5002.42116 3590.02407

Integrated Conn 6971.35316 6001.03185 5076.7984 3693.48877 6651.79576 5552.82496 4862.35879 3556.92483 6714.21266 5600.75871 4840.57379 3490.14 6816.79815 5618.88795 4945.9793 3595.48081 6705.06569 5685.12404 4958.20896 3633.92843

-

-

-

-

-

-

-

-

-

-

-

-

IterMerge

PostMerge

6718.72885 5713.96147 4911.83706 3565.78342 6584.89575 5551.75766 4801.07043 3446.30952 6654.56156 5520.3482 4802.62151 3373.72568 6633.74908 5529.78331 4789.39302 3542.17183 6685.90989 5607.57718 4880.22641 3507.30482

Table D.8: Distances from points to subspaces using unfiltered input data, with 20 clusters. Blue indicates the best value in each row. Note that the algorithm Agglo-ST performs poorly on unfiltered data (see Figure 4.2), resulting in clusters unsuitable for Conn-Subspace. This is because one large cluster dominates while the others are too small to compute subspaces for m′ dimensions. m′

Agglo-ST

Conn-Agglo-Euc

Conn-Agglo-ST

Conn-KMeans++

Conn-Ward

GraphConvSC EKGCSC

GraphConvSC EGCSC SSC-OMP EnSC

5 10 15 30 5 10 15 30 5 10 15 30 5 10 15 30 5 10 15 30 5 10 15 30 5 10 15 30

initial clustering 7783.69078 6825.71836 6245.5336 5098.51104 6802.40018 5661.22009 4892.76619 3413.25815 6798.59783 5620.91742 4848.81668 3369.86502 6891.16075 5762.35931 5007.12359 3553.72154 7142.41807 6185.16338 5540.21058 4224.92307 13597.25370 13509.21213 13604.91292 13604.26299 13165.23867 12989.97725 13504.44223 13620.71771 14689.48839 14458.67851

7242.28414 6224.32097 5301.69381 4093.73849 6711.08028 5681.3181 4901.48325 3405.64826 6784.49403 5920.02513 4691.61664 3458.22741 6977.26662 5638.61196 4901.24059 3460.08787 6920.79741 5931.42334 5293.17786 3988.19091

Smooth Merge 7179.5592 6005.97634 3934.63734 6917.901 5495.87935 4963.01987 3589.80958 6790.3771 5797.54924 4916.4418 3332.83429 6800.17259 5569.52283 4874.53135 3564.00204 6852.91056 5984.3858 5328.62205 3755.66247

Integrated Conn 7528.35562 6318.75352 5771.76803 4865.15531 6703.29191 5482.47043 4821.95063 3386.6587 6689.41092 5544.80481 4768.92015 3285.58812 6713.21372 5503.7572 4761.15466 3452.07467 6997.92878 5944.67975 5142.70403 3655.8976

-

-

-

-

-

-

-

-

-

-

-

-

IterMerge

PostMerge

6650.41124 5924.41911 4936.45871 4195.20874 6560.96462 5430.55569 4647.6343 3237.95838 6605.36159 5448.75877 4648.76742 3255.90307 6596.3755 5508.66299 4702.8914 3346.93706 6702.55058 5776.47192 4936.7126 3621.75654

Table D.9: Distances from points to subspaces using filtered input data, with 25 clusters. Blue indicates the best value in each row. m′

Agglo-ST

Conn-Agglo-Euc

Conn-Agglo-ST

Conn-KMeans++

Conn-Ward

GraphConvSC EKGCSC

GraphConvSC EGCSC SSC-OMP EnSC

5 10 15 30 5 10 15 30 5 10 15 30 5 10 15 30 5 10 15 30 5 10 15 30 5 10 15 30

initial clustering 6935.05471 5912.08858 5211.14154 3802.41450 6622.15663 5538.29201 4799.06170 3340.89806 6765.17821 5649.72727 4884.40787 3395.02750 6846.60590 5777.17176 5045.16003 3609.66688 6804.97968 5748.56173 5018.93476 3602.40923 12231.26304 12240.14617 11942.24078 11934.06822 11972.89186 11747.07094 11735.43384 11825.74653 14677.86061 12982.57873

6604.73022 5439.0021 4618.39123 3352.27087 6667.53307 5462.91872 4590.01205 3107.19511 6548.14257 5325.49801 4608.88795 3285.17893 6546.40152 5859.8821 4944.36625 3385.55874 6714.11371 5517.56823 4739.56481 3312.15312

Smooth Merge 6664.46905 5413.75331 4528.27139 3522.85521 6530.66165 5352.08113 4454.23534 3230.85019 6532.95274 5358.70727 4606.41312 3224.51674 6404.29306 5395.78375 4616.90696 3247.91839 6704.62957 5495.12792 4805.99741 3343.36821

Integrated Conn 6684.20784 5599.36789 4714.31267 3447.12802 6385.11851 5276.7694 4549.961 3195.57073 6498.09127 5354.7229 4574.74659 3256.04637 6645.10238 5476.8037 4645.33832 3230.89165 6562.28748 5499.27583 4809.59483 3397.42704

-

-

-

-

-

-

-

-

-

-

-

-

IterMerge

PostMerge

6477.34752 5433.27462 4679.71563 3325.62566 6350.81598 5229.74159 4521.96939 3126.14016 6399.39178 5303.26748 4527.33796 3200.7043 6412.88878 5327.48569 4581.94485 3239.08517 6544.75466 5387.06249 4653.33822 3373.04519

Table D.10: Distances from points to subspaces using unfiltered input data, with 25 clusters. Blue indicates the best value in each row. Note that the algorithm Agglo-ST performs poorly on unfiltered data (see Figure 4.2), resulting in clusters unsuitable for Conn-Subspace. This is because one large cluster dominates while the others are too small to compute subspaces for m′ dimensions. m′

Agglo-ST

Conn-Agglo-Euc

Conn-Agglo-ST

Conn-KMeans++

Conn-Ward

GraphConvSC EKGCSC

GraphConvSC EGCSC SSC-OMP EnSC

5 10 15 30 5 10 15 30 5 10 15 30 5 10 15 30 5 10 15 30 5 10 15 30 5 10 15 30

initial clustering 7611.33551 6658.96184 6109.23328 4992.56736 6648.23889 5488.59572 4708.85936 3218.62384 6585.85218 5375.51503 4585.05939 3080.66145 6791.59278 5689.86451 4942.47421 3501.83054 6964.9354 6010.5925 5355.69749 4044.59545 13250.46585 13528.28216 13596.14924 13576.80709 12476.49544 12998.51352 13064.42761 13583.19723 14473.09811 14351.16981

6985.54103 5882.39635 3934.40832 6603.82734 5564.22328 4760.16073 3301.88318 6500.81257 5347.90321 4523.11896 2992.33895 6651.56208 5688.4807 4604.54066 3264.31383 6725.90101 5508.17897 4895.42358 3385.23974

Smooth Merge 6780.14297 6533.52799 5428.61674 4606.55632 3157.41492 6552.52463 5268.33666 4528.61575 2938.29087 6653.9628 5548.82075 4541.47083 3192.40442 6791.40374 5585.96961 4868.43329 3423.68021

Integrated Conn 7378.82363 6242.66174 5520.47426 4761.62842 6444.52568 5306.91418 4680.51634 3200.42369 6328.12887 5283.01416 4452.95515 3028.98765 6434.88511 5281.69382 4715.21944 3302.94252 6687.39939 5530.5792 4645.91934 3335.34227

-

-

-

-

-

-

-

-

-

-

-

-

IterMerge

PostMerge

6431.27588 5484.5064 4963.57033 3122.85764 6401.5226 5239.18674 4414.99156 2947.91908 6375.37258 5267.75215 4467.76037 3016.05834 6459.77151 5285.76723 4501.4064 3122.85764 6586.38741 5504.71253 4757.14492 3361.74465

Related documents

• SCAN: Enhance Time Series Anomaly Detection via Multi-Scale Neighborhood-Centered Clustering#287123
• How to evaluate clustering with ground truth?#310858
• Integrating local spatial outliers into context-aware spatial prediction models#398547
• El estudio científico del Turismo a través de la dimensión geográfica#645046
• Spatiotemporal subspace feature tracking by mining discriminatory characteristics#957608