--- title: "Toward Convex Manifolds: A Geometric Perspective for Deep Graph Clustering of Single-cell RNA-seq Data" authors: - "Nairouz Mrabah" - "Mohamed Mahmoud Amar" - "Mohamed Bouguessa" - "Abdoulaye Banire Diallo" doi: "10.24963/ijcai.2023/540" url: "https://ijcai.org/proceedings/2023/540" pdf: "https://ijcai.org/proceedings/2023/0540.pdf" venue: "Proceedings of the Thirty-Second International Joint Conference on Artificial Intelligence" year: 2023 pages: "4855-4863" track: "Main Track" keywords: - "Machine Learning: ML: Clustering" - "Multidisciplinary Topics and Applications: MDA: Bioinformatics" - "Machine Learning: ML: Unsupervised learning" --- # Toward Convex Manifolds: A Geometric Perspective for Deep Graph Clustering of Single-cell RNA-seq Data Nairouz Mrabah, Mohamed Mahmoud Amar, Mohamed Bouguessa, Abdoulaye Banire Diallo *Proceedings of the Thirty-Second International Joint Conference on Artificial Intelligence*, pages 4855-4863. DOI: [10.24963/ijcai.2023/540](https://doi.org/10.24963/ijcai.2023/540) ## Abstract The deep clustering paradigm has shown great potential for discovering complex patterns that can reveal cell heterogeneity in single-cell RNA sequencing data. This paradigm involves two training phases: pretraining based on a pretext task and fine-tuning using pseudo-labels. Although current models yield promising results, they overlook the geometric distortions that regularly occur during the training process. More precisely, the transition between the two phases results in a coarse flattening of the latent structures, which can deteriorate the clustering performance. In this context, existing methods perform euclidean-based embedding clustering without ensuring the flatness and convexity of the latent manifolds. To address this problem, we incorporate two mechanisms. First, we introduce an overclustering loss to flatten the local curves. Second, we propose an adversarial mechanism to adjust the global geometric configuration. The second mechanism gradually transforms the latent structures into convex ones. Empirical results on a variety of gene expression datasets show that our model outperforms state-of-the-art methods. ## Full text ## 1 Introduction Single-cell RNA sequencing (scRNA-seq) facilitates the study of individual cells, thus providing an in-depth understanding of the basic unit of biology. The assessment of transcriptional heterogeneity within a cell population can provide answers to a wide range of biological research questions. For example, scRNA-seq datasets have been commonly used in the discovery of new cell types [Saviano et al., 2020], the identification of highly pathogenic cells in tumor tissue [Tirosh et al., 2016], the detection of hyper-responsive immune cells [Shalek et al., 2014], the analysis of resistance to treatments [Ocasio et al., 2019], and the study of developmental processes in specific cellular conditions such as cancer [Brady et al., 2017] and pulmonary epithelium differentiation [Treutlein et al., 2014]. Since manual cell annotation is labor intensive and can be prone to mislabeling, unsupervised cell identification remains a critical research direction. ∗Contact Author The analysis of scRNA-seq data heavily relies on cell type identification. To accomplish this task, various clustering algorithms have been used. In particular, researchers have applied traditional methods such as k-means, spectral clustering [Park and Zhao, 2018], density-based clustering [Ester et al., 1996], and hierarchical clustering [Johnson, 1967] to detect cell populations. Furthermore, dimensionality reduction techniques including PCA, TSNE, and UMAP have been applied before clustering to overcome the curse of dimensionality. Unfortunately, the traditional clustering algorithms often produce suboptimal results. This is because the transcriptomic data have complex relationships between cells within a cluster, significant sparsity caused by the dropout events, and strong technical variability of the gene expression levels. To address the shortcomings of traditional clustering approaches, CIDR [Lin et al., 2017], MAGIC [van Dijk et al., 2017], and SAVER [Huang et al., 2018], first impute the missing values, also known as dropouts, and then cluster the imputed data. Despite the positive effects of imputation, CIDR, MAGIC, and SAVER struggle to capture the complex structure inherent to scRNA-seq data. Some other approaches, such as SIMLR [Wang et al., 2017], and MPSSC [Park and Zhao, 2018], leverage multi-kernel spectral clustering to learn robust similarity measures. However, computing the Laplacian matrix is time-consuming, which in turn precludes its applicability to large-scale datasets. Moreover, these methods overlook important characteristics of the transcriptional data, such as zero inflation and over-dispersion. Recently, the deep clustering paradigm has found success in many biological applications. This paradigm relies on a two-fold strategy consisting of self-supervised learning followed by a pseudo-supervision task. During the initial phase, a pretext task is solved as a proxy to acquire high-level representations. The second phase simultaneously learns clustering assignments and clustering-oriented embeddings by training with pseudo-labels. Existing deep clustering methods for scRNA-seq data include scziDesk [Chen et al., 2020], scDC [Tian et al., 2019], and scDCC [Tian et al., 2021]. These models operationalize auto-encoding architectures. Furthermore, all of them integrate the Zero-Inflated Negative Binomial distribution (ZINB) to model the key properties of the gene expression data (over-dispersion, zero inflation, and discreteness). However, they fail to consider cell-cell relationships, which can make the clustering task more challenging. As a result, more recent models, such as scTAG [Yu et al., 2022] and scGAE [Luo et al., 2021], leverage graph neural networks (GNNs) to preserve the neighborhood relationships. Although deep clustering has achieved promising results for scRNA-seq data, existing studies fail to consider two important issues: 1) the suitability of the geometric configuration obtained at the end of the pretraining phase to the clustering task, and 2) the geometric distortions that regularly occur during the training process. Several measures can be used to analyze the geometry of a manifold, such as curvature and torsion. In this work, we focus on two measures, namely Intrinsic Dimension (ID) and Linear Intrinsic Dimension (LID). Intuitively, ID computes the number of independent directions a point can move in without leaving the manifold. The LID computes the smallest dimension of a linear subspace that best approximates the manifold in a least-squares sense at each point. An increasing discrepancy between ID and LID during the training process indicates that the ongoing geometric transformation makes the manifold more curved. The authors in [Mrabah et al., 2022] have conducted a geometric investigation to examine the behavior of ID and LID under the deep clustering paradigm. They have shown that self-supervised learning (i.e., pretraining phase) leads to the emergence of curved latent manifolds with low intrinsic dimensions. Moreover, introducing pseudo-supervised learning after the first phase flattens the latent structures coarsely, which in turn degrades the clustering performance by twisting the curvatures [Mrabah et al., 2022]. The main goal of our approach is to avoid this sharp geometric transition caused by minimizing a euclidean-based clustering loss after the first phase. Without prior knowledge, there is no systematic way to identify the optimal non-euclidean metric that can capture the latent similarities. Therefore, it is crucial to gradually adjust the geometric configuration of the latent manifolds, so that it becomes appropriate for euclidean-based clustering. We propose a novel single-cell graph auto-encoder model that follows the deep clustering paradigm. Unlike previous methods, our approach tackles the coarse geometric transition between pretraining and clustering. We incorporate two mechanisms to address this problem. We argue that enforcing local flatness and global convexity in a gradual manner is favorable for euclidean-based embedding clustering. Our first mechanism targets the local geometric configuration. In particular, we supply the pseudo-supervision module with an overclustering loss. This loss is used to flatten the strong local curves without destroying the global structure. The second mechanism adjusts the global geometric configuration. More precisely, we introduce an adversarial loss to gradually transform the latent manifolds into convex ones. We build a discriminator for each cluster. Our discriminators are trained to make the auto-encoder push convex combinations of samples with high-confidence clustering assignments inside their corresponding manifold. By ensuring convexity, the clusters become linearly separable [Boyd et al., 2004], which in turn improves the clustering performance. Contributions. (1) We establish the first single-cell deep clustering model that tackles the coarse geometric transition between the first and second training phases. (2) We propose a mechanism that performs a local geometric transformation of the latent manifolds. In particular, we introduce an overclustering loss to flatten the local curves without destroying the global structures. (3) We design a second mechanism that adjusts the global geometric configuration. More precisely, we introduce an adversarial loss to gradually transform the latent manifolds into convex ones. (4) We conduct several experiments on real-world scRNA-seq datasets. The obtained results confirm the validity of our contributions and show that the proposed model outperforms state-of-the-art methods. ## 2 RelatedWork We discuss two main strategies: 1) the deep embedding methods, which perform clustering and feature learning separately, 2) the deep clustering models, which perform clustering and feature learning jointly. Other relevant methods outside the field of single-cell research are discussed in Appendix A (all appendices are provided in the Supplementary Material † ). Deep embedding methods. In [Yu et al., 2021], scGMAI leverages a deep auto-encoder as an imputation technique, followed by the Fast Independent Component Analysis (FastICA) for dimensionality reduction, and a Gaussian Mixture Model (GMM) for clustering. Instead of performing vanilla reconstruction, scDCA [Eraslan et al., 2019] operationalizes a denoising auto-encoder. The decoding task of scDCA is modeled by the ZINB distribution to capture the characteristics of the gene expression data. In another work, scGNN [Wang et al., 2021] achieves multi-modal reconstruction using several auto-encoders. Moreover, this model relies on a GNN to aggregate cell-cell relationships. A left-truncated GMM is applied to uncover regulatory signals within the gene expression count matrix. Unfortunately, the existing deep embedding methods overlook the suitability of the geometric configuration for the clustering task. More precisely, selfsupervised learning leads to curved latent manifolds, thus we do not expect euclidean-based approaches such as GMM or k-means to identify the clustering structures effectively. Deep clustering methods. Compared with scDCA, scDC [Tian et al., 2019] combines the denoising task with an embedding clustering process. Both tasks are optimized jointly. More precisely, scDC adopts the clustering loss of DEC [Xie et al., 2016]. To improve scDC, scDCC [Tian et al., 2021] expresses prior information and domain-specific knowledge as soft pairwise constraints, then integrate these constraints as an additional term in the final loss function. In another work, scziDesk [Chen et al., 2020] uses a different loss function compared with scDC. It employs a weighted k-means loss to strengthen the initial pairwise similarities after encoding. Recently, several methods leverage a graph convolutional encoder to maintain the proximity of similar data points in the latent space. For instance, scGAE [Luo et al., 2021] replaces the vanilla auto-encoding framework of scDC with a graph convolutional auto-encoder. In [Yu et al., 2022], scTAG adopts a topology adaptive graph convolutional auto-encoder. Unlike scGAE and scTAG, scDSC [Gan et al., 2022] encodes the gene expression data and the associated structural similarities using two encoders, and then maximizes the mutual information between their outputs. However, all these methods overlook the coarse geometric transformation caused by the sharp transition between pretraining and clustering. [† https://github.com/MMAMAR/scTConvexMan](https://github.com/MMAMAR/scTConvexMan) ## 3 Proposed Approach We define the gene expression count matrix X = (xij) ∈ Rns×ng , where ng represents the number of genes, ns is the number of cells, and the xij element is the count of the j th gene in the ith cell. We consider the problem of clustering the scRNA-seq data X and we denote by nc the corresponding number of clusters. Let C = (c1, ..., cns ) be the sequence of clustering indices for each cell such that ci ∈ {1, ..., nc}. The cell representations for scRNA-seq data can have pronounced similarities based on the gene expression count. To this end, we build a non-directed graph to capture the local structure characterizing these pronounced similarities. We define the graph structure by considering the relationship between each cell and its nearest neighbors. Formally, let G = (V, E, X) be an attributed graph, where V = {v1, v2, ..., vns } is a set of ns nodes associated with different cells; eij ∈ E specifies the existence of an edge between the ith and j th nodes; and X is used as the attribute matrix of G. Two nodes vi and vj are connected if and only if xi (xj, respectively) is located within the k nearest neighbors of xj (xi, respectively) based on the Euclidean distance. Let A = (aij) ∈ Rns×ns be the resultant adjacency matrix of G, where aij = 1 if vi and vj are connected, otherwise aij is set equal to zero. We build a GNN architecture to extract high-level representations from the graph G. In particular, we devise a graph auto-encoder [Kipf and Welling, 2016] to learn the node embeddings. The first component is an encoder fE. The encoding process relies on LE graph convolutional layers [Kipf and Welling, 2017] that project the graph G into a lowdimensional latent space. We denote by Z = fE(X,A) ∈ Rns×d the matrix of embedded codes obtained by the encoding process, where d is the dimension of the latent space. The graph convolutional operation of the lth layer is described: where W (l) E represents the training weights of the lth layer; ReLU is the corresponding activation function; Z (l) is the output of this layer, Z (0) =X, and Z (LE) = Z; ˜∆- 12 ˜ A ˜∆- 1 2 is the normalized adjacency matrix such that ˜ A = A + I, I ∈ Rns×ns is the identity matrix, and ˜∆ = Diag( ˜ A 1ns ). The second component of the auto-encoding model is a two-head decoder. The decoding process aims to reconstruct the input graph G. We denote by GD the attributed graph generated by the decoding process. The first head f (1) D is devoted to computing the structure of GD using the inner product operation. The output of the first head is the adjacency matrix of GD denoted ˆ A = f (1) D (Z). The second head f (2) D is devoted to computing the attributes of GD . It has LD fully-connected layers with ReLU activations. The output of the second head is the feature matrix ˆX = f (2) D (Z). We denote by W (l) D the training weights of the lth decoding layer and all the training weights of our auto-encoder are denoted by the set W. In addition to the auto-encoding framework, our approach includes an adversarial training strategy. We use the graph auto-encoder constrained to the second decoding head (that is, f (2) D ◦ fE) as the generator network. Moreover, we devise nc fully connected discriminators f (k) adv : Rd → [0, 1] with Ladv layers, such that k ∈ {1, ..., nc}. We denote by W (k,l) adv the training weights of the lth layer in the kth discriminator. ## 3.1 Training Losses According to the deep clustering paradigm, we train our model based on self-supervision for the pretraining phase, then we introduce a pseudo-supervision loss for the clustering phase. To avoid the coarse geometric transition from pretraining to clustering, we propose two mechanisms. The first mechanism targets the local geometric configuration of the latent manifolds. In particular, we supply the pseudosupervision module with an overclustering loss to flatten the strong local curves without destroying the globally-curved shape. The second mechanism adjusts the global geometric configuration. More precisely, we introduce an adversarial loss to gradually transform the manifolds into convex ones. ## Self-supervision Module It consists of tackling a secondary task that requires a certain level of data understanding to be solved. The self-supervision loss is optimized for the first and second phases. The pretext task of our approach consists of reconstructing the input graph G. We model the decoding process for the reconstruction task by the probability distribution p( ˆ A, ˆX |Z). We develop a two-task decoding design to perform both structure reconstruction and feature reconstruction. Accordingly, we factorize the joint distribution p( ˆ A, ˆX |Z) into two independent distributions p( ˆ A|Z) and p( ˆX |Z). We parameterize the probability distribution p( ˆ A | Z) by the first decoding head. We opt for the Bernoulli distribution to capture the binary nature of each generated edge as described by: where the probability of generating each edge B(ˆaij |βij) is a Bernoulli distribution parametrized by βij = sigmoid(zT i zj). For feature reconstruction, we parameterize the probability distribution p( ˆX |Z) by the second decoding head. We operationalize the ZINB distribution [Wang et al., 2021] to capture the key aspects of the gene expression data: 1) zeroinflation (high-sparsity) caused by the true and dropout zeros, 2) discreteness, and 3) over-dispersion (variance greater than the mean). We model the distribution p( ˆX | Z) as follows: where Z(ˆxi|πi, µi, σ2 i ) is the ZINB distribution parameterized by πi, µi, σ2 i ∈ Rng ; πi = sigmoid(Wπf (2) D (zi)) is the zero rate vector of this distribution and Wπ denotes the weight matrix of π; µi = exp(Wµf (2) D (zi)) is the mean of the associated negative binomial and Wµ is the weight matrix of µ; σ2 = exp(Wσf (2) D (zi)) is the variance of the negative binomial and Wσ is the weight matrix of σ2 . The complete expression of the ZINB mass function with respect to its parameters πi, µi, σ2 i is provided in Appendix B to save space. The self-supervision module reconstructs the input graph by minimizing the loss function LSS. More precisely, we optimize the log-likelihood of the input graph based on the distribution p( ˆX, ˆ A| Z) as described by: where LX(LA, respectively) is the log-likelihood of the gene count matrix X (adjacency matrix A, respectively). ## Pseudo-supervision Module It is trained using pseudo-labels during the second phase to achieve the primary task (i.e., clustering). We design a pseudo-supervision strategy that puts into action two loss functions: a clustering loss function Lc and an overclustering loss function Lo. The overclustering procedure refers to clustering the data into a number of categories no greater than the real number of clusters nc. The goal of the clustering loss is to construct clusteringfriendly latent structures that reduce the within-cluster variance and maximizes the between-cluster variance. We model the clustering process using two distributions p(C|Z) and q(C|Z). Each cluster is characterized by its center. Let {Φcj }nc j=1 be the set of trainable clustering centers. We initialize the points of this set, after the pretraining phase, by encoding the centers obtained from applying spectral clustering on A. The first distribution p(C|Z) measures the soft clustering assignments and factorizes as p(C|Z) = ∏ns i=1 p(ci | zi). We compute pcij , which is the probability of assigning the latent code zi to the j th cluster (i.e., p(ci = j|zi)), according to the Student's t-distribution as described by: The second distribution q(C|Z) models the target clustering assignments (i.e., pseudo-labels). The factorization of this distribution is q(C|Z) = ∏nc i=1 q(ci | zi). We compute the target assignments q(ci|zi) of the ith sample iteratively based on the high-confidence assignments of p(ci|zi). We denote by qc ij the probability value q(ci = j|zi). Let τ 1 i (τ 2 i , respectively) be the highest (second highest, respectively) assignment score of the set {pcij }nc j=1 . We denote by Ωc = {xi ∈ Rng | τ 1 i - τ 2 i ⩾ βc} the set of reliable samples for the clustering task, and βc is a fixed threshold. To emphasize the high-confidence assignments, we compute qc ij using: Our clustering loss Lc is the KL (Kullback-Leibler) divergence between the soft assignment distribution p(C|Z) and the target distribution q(C|Z) as described by: In addition to the clustering loss, we design an overclustering task to flatten the local structures without twisting the latent manifolds. More precisely, the goal of the overclustering process is to smooth the local curves while preserving the globally-curved geometric structures as shown by [Mrabah et al., 2022]. Let no be the number of clusters for the overclustering task. We model the overclustering process using two distributions p(O|Z) and q(O|Z), where O = (o1, ..., ons ) is a sequence of random variables representing the overclustering indices for each cell, such that oi ∈ {1, ..., no}. Let {Φoj }no j=1 be the set of trainable overclustering centers. We initialize these points, after the pretraining phase, by applying k-means to Z. The first distribution p(O|Z) measures the soft overclustering assignments and factorizes as p(O|Z) = ∏ns i=1 p(oi | zi). We denote by poij the probability p(oi = j|zi), which is computed as described by: The second distribution q(O|Z) models the target overclustering assignments. The factorization of this distribution is q(O|Z) = ∏ns i=1 q(oi | zi). We denote by qo ij the probability value q(oi = j|zi). Let λ1i (λ2i , respectively) be the highest (second highest, respectively) assignment score of the set {poij }no j=1 . We denote by Ωo = {xi ∈ Rng | λ1i -λ2i ⩾ βo} the set of reliable samples for the overclustering task and βo is a fixed threshold. We compute qo ij as described by: The overclustering loss Lo is the KL (Kullback-Leibler) divergence between the soft assignment distribution p(O|Z) and the target distribution q(O|Z) as described by: ## Adversarial Module It handles the global configuration of the latent space. Our adversarial mechanism is responsible for gradually transforming the globally-curved structures into convex ones. During the training process, the number of reliable samples for clustering increases progressively. We train the adversarial module using samples from the set Ωc. As a result, the adversarial module gradually enforces the convexity constraint from the core to the outer parts of the manifolds. Definition 1. A convex combination x of k points x1, ..., xk is a linear combination of these points: x = θ1x1 + ... + θkxk, where θ1, ..., θk ∈ R+ and ∑i θi = 1. Let Ck be the set of samples belonging to the kth cluster. We compute the clustering index of the sample xi using the formula maxj(pij). Given the representations of m cells from the kth cluster, x1, ..., xm ∈ Ωc ∩ Ck, we randomly select m coefficients θ1, ..., θm ∈ R+, such that ∑i θi = 1. We construct the point ˆx(k) cvx as a convex combination of the points x1, ..., xm in the latent space, as articulated by: For each cluster Ck, we generate ncvx points according to Eq. (11). We obtain different samples by randomly selecting the initial points x1, ..., xm. We denote by Γk the set of samples generated for the kth cluster. Each discriminator f (k) adv is trained to distinguish between real samples reconstructed from the set Ωc ∩ Ck and generated samples belonging to the set Γk. The loss function L(k) d of the kth discriminator is: We use the graph auto-encoder constrained to the second decoding head as the generator for the adversarial module. We train the generator to fool the discriminator into considering the generated points ˆxcvx as real reconstruction samples. The loss function of the generator Lg is described by: The generator is trained to incorporate convex combinations of the latent codes into the embedded manifolds. The gradual inclusion of each cluster's convex hull leads to smoother latent structures and a higher level of convexity. Clustering convex sets is relatively straightforward because they have well-defined boundaries. Thus, it is easy to separate the points inside the set from the points outside it. Theorem 1 shows that disjoint convex sets are linearly separable. Theorem 1. Hyperplane separation [Boyd et al., 2004]: if A and B are disjoint nonempty convex sets from Rd , then there exists a vector u ∈ Rd - {0} and a number c ∈ R, such that: and the hyperplane {x | uT x ⩾ c} separates A and B. ## 3.2 Training Strategy Our model undergoes a two-phase training process. We pretrain for T1 iterations. For the first phase, we update the autoencoder parameters to minimize the self-supervision function LSS. For the second phase, we alternate between three training steps. The first step consists of minimizing the loss function LSS + γcLc. In the second step, we minimize the loss function LSS + γoLo. In the third step, we minimize the loss function LSS + γgLg. The hyperparameters γc, γo, and Figure 1: Illustration of the clustering phase of scTCM. The matrix Zcvx ∈ Rncvx×d represents the convex combinations of the latent codes. It worth to highlight that we use a discriminator for each cluster. We present a single discriminator to simplify the illustration. ![Image](https://ijcai.org/proceedings/2023/fulltext/540/image_000000_330e8b455295e52d6bcd47734c250d051df1d142e72b6238ba79b47b40d9111e.png) γg control the trade-off between self-supervision and the remaining tasks (i.e., clustering, overclustering, and adversarial training). Moreover, we train the discriminators at each step by minimizing the loss functions L(k) d . Adam optimizer is used for all training phases and steps. The second phase ends when the average ID of all manifolds becomes equal to the LID. We illustrate the framework of the second phase in Figure 1. We provide our algorithm and its computational complexity analysis in Appendix C and Appendix D, respectively. ## 4 Experiments We carry out multiple experiments on eight real-world scRNA-seq datasets of various species. These datasets are collected using popular sequencing platforms. Information on the data and preprocessing steps can be found in Appendix E. We compare the performance of our method, scTCM (deep clustering of single-cell RNA-seq data: Toward Convex Manifolds), against state-of-the-art clustering techniques for scRNA-seq. Our comparison includes three deep embedding methods (scDCA [Eraslan et al., 2019], scGMAI [Yu et al., 2021], and scGNN [Wang et al., 2021]), three deep clustering methods (scziDesk [Chen et al., 2020], scDC [Tian et al., 2019], and scDCC [Tian et al., 2021]), and three deep graph clustering methods (scTAG [Yu et al., 2022], scDSC [Gan et al., 2022], and scGAE [Luo et al., 2021]). All of these approaches are discussed in the related work section. We use the official code provided by the authors for each baseline and we tune the hyperparameters if no specific instructions are given. To ensure reproducibility, our code is available on https://github.com/MMAMAR/scTConvexMan. The geometric configuration is evaluated by measuring the ID and LID of the latent manifolds. A comprehensive explanation of the strategies applied to estimate ID and LID can be found in Appendix F. We compare the efficiency of the different models based on their execution time in seconds. All of our experiments are performed under consistent hardware and software setups, as described in Appendix G. In Appendix H, we thoroughly outline the design and configuration of our approach, including the architecture, the learning rates, and any other hyperparameters. In particular, our approach has three data-dependent hyperparameters (βc, βo, and no), which are fixed using grid search as explained in Appendix H. We hold the remaining design choices constant across all the datasets. A qualitative evaluation is included in Appendix I. Table 1: Clustering results on eight scRNA-seq datasets. Best methods in red and second best in blue. | Metric | Dataset | Deep Graph Clustering - scTCM | Deep Graph Clustering - scTAG | Deep Graph Clustering - scGAE | Deep Graph Clustering - scDSC | Deep Clustering - scziDesk scDC | Deep Clustering - scziDesk scDC | Deep Clustering - scDCC | Deep Embedding - scDCA | Deep Embedding - scGMAI | Deep Embedding - scGNN | |----------|-----------------------------------|---------------------------------|---------------------------------|---------------------------------|---------------------------------|-----------------------------------|-----------------------------------|---------------------------|--------------------------|---------------------------|--------------------------| | | Muraro | 96.41 | 89.35 | 89.75 | 86.52 | 73.66 | 75.12 | 84.35 | 73.61 | 61.40 | 66.30 | | | Plasschaert | 97.46 | 89.01 | 83.83 | 91.76 | 93.56 | 52.60 | 61.34 | 62.66 | 48.07 | 53.56 | | | QX LM | 99.77 | 97.77 | 75.23 | 78.49 | 96.67 | 75.06 | 82.04 | 74.52 | 60.88 | 63.54 | | | QS Diaph | 99.08 | 98.39 | 55.51 | 94.94 | 97.82 | 71.84 | 71.26 | 72.30 | 47.93 | 68.04 | | | ACC QS Heart | 98.95 | 96.13 | 89.46 | 92.05 | 94.46 | 55.53 | 68.23 | 70.72 | 60.09 | 63.50 | | | QS LM | 99.63 | 99.17 | 62.56 | 86.06 | 97.80 | 60.64 | 74.68 | 58.53 | 57.98 | 85.41 | | | Wang Lung | 99.30 | 95.14 | 77.11 | 96.02 | 97.59 | 98.97 | 90.51 | 98.62 | 38.56 | 97.97 | | | Young | 85.30 | 81.35 | 70.72 | 67.88 | 79.81 | 58.75 | 67.41 | 64.57 | 42.00 | 50.30 | | | Average | 96.99 | 93.28 | 75.52 | 86.71 | 91.42 | 68.56 | 74.97 | 71.94 | 52.11 | 68.57 | | | Muraro Plasschaert QX LM QS Diaph | 89.80 | 82.92 | 83.53 | 85.00 | 77.76 | 75.49 | 83.84 | 76.21 | 71.68 | 64.10 | | | | 90.54 | 73.79 | 68.44 | 80.02 | 79.80 | 61.22 | 64.97 | 65.52 | 57.11 | 63.17 | | | | 99.01 | 93.64 | 79.39 | 80.61 | 91.28 | 84.05 | 87.98 | 83.69 | 76.39 | 74.43 | | | | 95.93 | 93.28 | 68.04 | 91.67 | 91.40 | 78.07 | 81.04 | 78.38 | 68.36 | 71.34 | | | NMI QS Heart | 95.90 | 89.43 | 78.39 | 89.58 | 87.75 | 65.31 | 72.75 | 70.61 | 69.41 | 73.20 | | | QS LM | 98.40 | 96.44 | 72.71 | 82.73 | 92.15 | 70.48 | 79.94 | 71.71 | 71.98 | 80.35 | | | Wang Lung | 92.30 | 71.50 | 58.62 | 74.54 | 81.46 | 90.11 | 57.99 | 87.38 | 34.32 | 84.61 | | | Young | 84.86 | 79.16 | 65.27 | 72.18 | 76.72 | 61.78 | 67.66 | 63.66 | 49.54 | 40.03 | | | Average | 93.34 | 85.02 | 71.79 | 82.04 | 84.79 | 73.31 | 74.52 | 74.64 | 62.34 | 68.90 | | | Muraro Plasschaert | 93.72 | 87.02 | 88.21 | 89.51 | 66.30 | 66.09 | 74.10 | 64.59 | 51.32 | 49.19 | | | | 94.70 | 76.97 | 70.64 | 86.99 | 85.88 | 40.70 | 49.76 | 50.69 | 57.11 | 46.40 | | | QX LM | 99.61 | 95.45 | 63.31 | 74.79 | 91.28 | 75.00 | 80.23 | 74.43 | 50.61 | 58.33 | | | QS Diaph | 97.88 | 96.37 | 42.90 | 96.15 | 93.69 | 64.79 | 66.03 | 64.98 | 41.11 | 56.09 | | | ARI QS Heart | 98.15 | 94.09 | 84.89 | 94.78 | 92.08 | 46.73 | 58.18 | 59.87 | 43.68 | 56.74 | | | QS LM | 99.32 | 98.11 | 53.79 | 90.89 | 94.29 | 53.84 | 65.77 | 54.44 | 48.99 | 79.68 | | | Wang Lung | 96.90 | 80.14 | 59.98 | 83.48 | 89.74 | 95.50 | 64.10 | 94.02 | 13.25 | 91.31 | | | Young | 79.74 | 70.61 | 58.87 | 54.61 | 66.81 | 44.69 | 52.85 | 50.06 | 32.87 | 28.98 | | | Average Muraro | 95.00 | 87.34 | 65.32 61 | 83.90 88 | 85.00 106 | 60.91 74 | 63.87 83 | 64.13 240 | 42.36 152 | 58.34 733 | | | Plasschaert | 330 766 | 54 159 | 490 | 242 | 378 | 165 | 216 | 546 | 589 | 1725 | | | QX LM | 387 | 70 | 189 | 120 | 278 | 73 | 82 | 198 | 321 | 816 | | | QS Diaph | 227 | 50 | 46 | 64 | 98 | 58 | 73 | 52 | 78 | 341 | | | Time QS Heart | | | | | | | | | | | | | QS LM | 576 | 87 | 221 | 134 | 406 | 164 | 273 | 261 | 300 | 1475 | | | Wang Lung | 240 | 63 | 46 | 52 | 115 | 56 | 78 | 66 | 86 | 413 | | | Young | 956 647 | 519 213 | 718 351 | 150 202 | 632 601 | 307 223 | 190 203 | 202 432 | 572 367 | 1867 1243 | Table 2: Impact of the proposed adversarial and overclustering mechanisms on the clustering performance of scTCM. Best results in bold. | Method - Adverserial | Method - Overclustering | Muraro - ACC | Muraro - NMI | Muraro - ARI | Muraro - Time | Plasschaert - ACC | Plasschaert - NMI | Plasschaert - ARI | Plasschaert - Time | Young - ACC | Young - NMI | Young - ARI | Young - Time | |------------------------|---------------------------|----------------|----------------|----------------|-----------------|---------------------|---------------------|---------------------|----------------------|---------------|---------------|---------------|----------------| | ✗ | ✗ | 93.49 | 86.80 | 90.75 | 101 | 96.60 | 88.51 | 93.06 | 646 | 82.11 | 82.41 | 74.33 | 548 | | ✗ | ✓ | 95.75 | 88.38 | 92.63 | 120 | 97.18 | 89.26 | 94.04 | 649 | 84.57 | 84.08 | 75.71 | 577 | | ✓ | ✗ | 93.87 | 86.93 | 91.54 | 325 | 97.00 | 89.29 | 93.44 | 732 | 84.22 | 84.01 | 77.59 | 632 | | ✓ | ✓ | 96.41 | 89.80 | 93.72 | 330 | 97.46 | 90.54 | 94.70 | 766 | 85.30 | 84.86 | 79.74 | 647 | Clustering results. Table 1 shows the clustering performance of our method against multiple state-of-the-art techniques. For each method, we report the average results of ten executions. We observe that the deep clustering approaches generally outperform the deep embedding methods. This finding underscores the importance of combining embedded learning and clustering. Additionally, we can see that the deep graph clustering methods tend to perform better than the deep clustering approaches. This finding supports the significance of extracting the graph from the gene expression count matrix to use it in the encoding process. Last but not least, ing transforms the global configuration into a convex one. In addition, we can see from Table 2 that our mechanisms do not induce significant run-time overhead. These results align with the computational complexity provided in Appendix D. Figure 3: Sensitivity of scTCM to the hyperparameters no and βo in terms of ACC and NMI. ![Image](https://ijcai.org/proceedings/2023/fulltext/540/image_000001_bb64250f968a43e192e09e474ef62978a75ab3066e18c32abec90b6fe0212131.png) Sensitivity analysis. We analyze the sensitivity of our model to the data-dependent hyperparameters (βc, βo, and no). We hold the remaining design decisions constant across all datasets. In Fig. 3, we study the sensitivity to the hyperparameters of the overclustering loss (i.e., βo and no). In Appendix J, we study the sensitivity to the hyperparameter of the clustering loss (i.e., βc). We found that the proposed model produces consistent results for a wide range of values. Especially, fixing no, scTCM yields stable clustering results as βo varies. As we can see, the ACC and NMI spike when no is equal to 30 and 10 on Muraro and Young, respectively. ## 5 Conclusion In this research, we propose a single-cell deep clustering model that addresses the coarse geometric transition between the pretraining and clustering phases. Our approach utilizes two mechanisms to make the latent representations suitable for euclidean-based clustering. More precisely, we focus on two geometric aspects: local flatness and global convexity. The first mechanism is a local geometric transformation that flattens the local curves without destroying the global structures. It consists of minimizing an overclustering loss. The second mechanism adjusts the global geometric configuration by transforming the latent manifolds into convex ones through an adversarial loss. Our method shows higher clustering performance against state-of-the-art approaches for scRNA-seq data. Empirical results provide strong evidence that this performance is imputed to the proposed mechanisms and particularly their ability to tackle the coarse geometric transition between pretraining and clustering. Figure 2: ID and LID of scTCM(-adv-overclus) and scTCM on Muraro (first row) and Plasschaert (second row). P1 stands for the first training phase (self-supervision) and P2 stands for the second training phase (pseudo-supervision). the proposed method yields the best performance in terms of ACC, NMI, and ARI compared to all the other methods, across all the datasets. On average, the NMI and ARI scores of scTCM are more than 7% higher than the scores of the second-best method (that is, scTAG). Unlike previous methods, scTCM has two mechanisms to avoid the coarse geometric transition between the pretraining and clustering phases. Geometric study. In Fig. 2, we investigate the evolution of the latent manifolds during the training process. To this end, we build a modified version of scTCM. We obtain this version, denoted scTCM(-adv-overclus), by ablating the adversarial module and the overclustering loss. For scTCM(-adv-overclus), we observe that the difference between ID and LID decreases rapidly after the pretraining phase. Within a few iterations, the latent manifolds undergo a coarse geometric transformation: from curved (evidenced by a significant difference between ID and LID) and lowdimensional structures to flattened and higher-dimensional ones. In contrast to scTCM(-adv-overclus), scTCMhas a less pronounced decrease in LID and a more moderate increase in ID even after a considerable number of iterations. The proposed mechanisms gradually flatten the latent manifolds. Unlike scTCM(-adv-overclus), the clustering-oriented representations of scTCM (Euclidean geometry perspective) are slowly constructed while protecting the curves from twisting. Ablation study. In this experiment, we analyze the impact of our contributions on clustering performance. First, we can see from Table 2 that the overclustering mechanism improves the clustering results consistently. Second, performing adversarial training also brings consistent improvement compared with scTCM(-adv-overclus). It is worth noting that the full model regularly outperforms all the other variants. These results validate the synergy of the proposed mechanisms. While overclustering flattens the local structures, adversarial train- ## Acknowledgments This work has been supported by Research Grants from the Natural Sciences and Engineering Research Council of Canada (NSERC). ## References 1. [Boyd et al., 2004] Stephen Boyd, Stephen P Boyd, and Lieven Vandenberghe. Convex optimization. Cambridge university press, 2004. 2. [Brady et al., 2017] Samuel W Brady, Jasmine A McQuerry, Yi Qiao, Stephen R Piccolo, Gajendra Shrestha, David F Jenkins, Ryan M Layer, Brent S Pedersen, Ryan H Miller, Amanda Esch, et al. Combating subclonal evolution of resistant cancer phenotypes. Nature communications, 8(1):1-15, 2017. 3. [Chen et al., 2020] Liang Chen, Weinan Wang, Yuyao Zhai, and Minghua Deng. Deep soft k-means clustering with self-training for single-cell rna sequence data. NAR ge- nomics and bioinformatics, 2(2):lqaa039, 2020. 4. [Eraslan et al., 2019] Gökcen Eraslan, Lukas M Simon, Maria Mircea, Nikola S Mueller, and Fabian J Theis. Single-cell rna-seq denoising using a deep count autoen- coder. Nature communications, 10(1):1-14, 2019. 5. [Ester et al., 1996] Martin Ester, Hans-Peter Kriegel, Jörg Sander, Xiaowei Xu, et al. A density-based algorithm for discovering clusters in large spatial databases with noise. In KDD, volume 96, pages 226-231, 1996. 6. [Gan et al., 2022] Yanglan Gan, Xingyu Huang, Guobing Zou, Shuigeng Zhou, and Jihong Guan. Deep structural clustering for single-cell rna-seq data jointly through au- toencoder and graph neural network. Briefings in Bioin- formatics, 23(2):bbac018, 2022. 7. [Huang et al., 2018] Mo Huang, Jingshu Wang, Eduardo Torre, Hannah Dueck, Sydney Shaffer, Roberto Bonasio, John I Murray, Arjun Raj, Mingyao Li, and Nancy R Zhang. Saver: gene expression recovery for single-cell rna sequencing. Nature methods, 15(7):539-542, 2018. 8. [Johnson, 1967] Stephen C Johnson. Hierarchical clustering schemes. Psychometrika, 32(3):241-254, 1967. 9. [Kipf and Welling, 2016] Thomas N Kipf and Max Welling. Variational graph auto-encoders. In NeurIPS workshop, volume 32, pages 1-3, 2016. 10. [Kipf and Welling, 2017] Thomas N Kipf and Max Welling. Semi-supervised classification with graph convolutional networks. In ICLR, 2017. 11. [Lin et al., 2017] Peijie Lin, Michael Troup, and Joshua WK Ho. Cidr: Ultrafast and accurate clustering through im- putation for single-cell rna-seq data. Genome biology, 18(1):1-11, 2017. 12. [Luo et al., 2021] Zixiang Luo, Chenyu Xu, Zhen Zhang, and Wenfei Jin. A topology-preserving dimensionality re- duction method for single-cell rna-seq data using graph au- toencoder. Scientific reports, 11(1):1-8, 2021. 13. [Mrabah et al., 2022] Nairouz Mrabah, Mohamed Bouguessa, and Riadh Ksantini. Escaping feature twist: A variational graph auto-encoder for node clustering. In IJCAI, pages 3351-3357, 2022. 14. [Ocasio et al., 2019] Jennifer Karin Ocasio, Benjamin Bab- cock, Daniel Malawsky, Seth J Weir, Lipin Loo, Jeremy M Simon, Mark J Zylka, Duhyeong Hwang, Taylor Dismuke, Marina Sokolsky, et al. scrna-seq in medulloblastoma shows cellular heterogeneity and lineage expansion sup- port resistance to shh inhibitor therapy. Nature communi- cations, 10(1):1-17, 2019. 15. [Park and Zhao, 2018] Seyoung Park and Hongyu Zhao. Spectral clustering based on learning similarity matrix. Bioinformatics, 34(12):2069-2076, 2018. 16. [Saviano et al., 2020] Antonio Saviano, Neil C Henderson, and Thomas F Baumert. Single-cell genomics and spatial transcriptomics: Discovery of novel cell states and cellular interactions in liver physiology and disease biology. Jour- nal of hepatology, 73(5):1219-1230, 2020. 17. [Shalek et al., 2014] Alex K Shalek, Rahul Satija, Joe Shuga, John J Trombetta, Dave Gennert, Diana Lu, Peilin Chen, Rona S Gertner, Jellert T Gaublomme, Nir Yosef, et al. Single-cell rna-seq reveals dynamic paracrine control of cellular variation. Nature, 510(7505):363-369, 2014. 18. [Tian et al., 2019] Tian Tian, Ji Wan, Qi Song, and Zhi Wei. Clustering single-cell rna-seq data with a model-based deep learning approach. Nature Machine Intelligence, 1(4):191-198, 2019. 19. [Tian et al., 2021] Tian Tian, Jie Zhang, Xiang Lin, Zhi Wei, and Hakon Hakonarson. Model-based deep embedding for constrained clustering analysis of single cell rna-seq data. Nature communications, 12(1):1-12, 2021. 20. [Tirosh et al., 2016] Itay Tirosh, Benjamin Izar, Sanjay M Prakadan, Marc H Wadsworth, Daniel Treacy, John J Trombetta, Asaf Rotem, Christopher Rodman, Christine Lian, George Murphy, et al. Dissecting the multicellular ecosystem of metastatic melanoma by single-cell rna-seq. Science, 352(6282):189-196, 2016. 21. [Treutlein et al., 2014] Barbara Treutlein, Doug G Brown- field, Angela R Wu, Norma F Neff, Gary L Mantalas, F Hernan Espinoza, Tushar J Desai, Mark A Krasnow, and Stephen R Quake. Reconstructing lineage hierarchies of the distal lung epithelium using single-cell rna-seq. Na- ture, 509(7500):371-375, 2014. 22. [van Dijk et al., 2017] David van Dijk, Juozas Nainys, Roshan Sharma, Pooja Kaithail, Ambrose J Carr, Kevin R Moon, Linas Mazutis, Guy Wolf, Smita Krishnaswamy, and Dana Pe'er. Magic: A diffusion-based imputation method reveals gene-gene interactions in single-cell rna- sequencing data. BioRxiv, page 111591, 2017. 23. [Wang et al., 2017] Bo Wang, Junjie Zhu, Emma Pierson, Daniele Ramazzotti, and Serafim Batzoglou. Visualization and analysis of single-cell rna-seq data by kernel-based similarity learning. Nature methods, 14(4):414-416, 2017. 24. [Wang et al., 2021] Juexin Wang, Anjun Ma, Yuzhou Chang, Jianting Gong, Yuexu Jiang, Ren Qi, Cankun Wang, Hongjun Fu, Qin Ma, and Dong Xu. scgnn is a novel graph neural network framework for single-cell rna-seq analyses. Nature communications, 12(1):1-11, 2021. 25. [Xie et al., 2016] Junyuan Xie, Ross Girshick, and Ali Farhadi. Unsupervised deep embedding for clustering analysis. In ICML, pages 478-487, 2016. 26. [Yu et al., 2021] Bin Yu, Chen Chen, Ren Qi, Ruiqing Zheng, Patrick J Skillman-Lawrence, Xiaolin Wang, An- jun Ma, and Haiming Gu. scgmai: a gaussian mixture model for clustering single-cell rna-seq data based on deep autoencoder. Briefings in Bioinformatics, 22(4):316, 2021. 27. [Yu et al., 2022] Zhuohan Yu, Yifu Lu, Yunhe Wang, Fan Tang, Ka-Chun Wong, and Xiangtao Li. Zinb-based graph embedding autoencoder for single-cell rna-seq interpreta- tions. AAAI, 36(4):4671-4679, 2022.