SNEkhorn: Dimension Reduction with Symmetric Entropic Affinities
Authors: Hugues Van Assel
Titouan Vayer
Rémi Flamary
Nicolas Courty
Venue: NeurIPS 2023 poster (2023)
Keywords: Dimension Reduction
Optimal Transport
Affinities

Abstract: Many approaches in machine learning rely on a weighted graph to encode the
similarities between samples in a dataset. Entropic affinities (EAs), which are notably used in the popular Dimensionality Reduction (DR) algorithm t-SNE, are particular instances of such graphs. To ensure robustness to heterogeneous sampling densities, EAs assign a kernel bandwidth parameter to every sample in such a way that the entropy of each row in the affinity matrix is kept constant at a specific value, whose exponential is known as perplexity. EAs are inherently asymmetric and row-wise stochastic, but they are used in DR approaches after undergoing heuristic symmetrization methods that violate both the row-wise constant entropy and stochasticity properties. In this work, we uncover a novel characterization of EA as an optimal transport problem, allowing a natural symmetrization that can be computed efficiently using dual ascent. 
The corresponding novel affinity matrix derives advantages from symmetric doubly stochastic normalization in terms of clustering performance, while also effectively controlling the entropy of each row thus making it particularly robust to varying noise levels. Following, we present a new DR algorithm, SNEkhorn, that leverages this new affinity matrix. We show its clear superiority to state-of-the-art approaches with several indicators on both synthetic and real-world datasets.

SNEkhorn: Dimension Reduction with Symmetric Entropic Affinities
Hugues Van Assel
Titouan Vayer
ENS de Lyon, CNRS UMPA UMR 5669
Univ. Lyon, ENS de Lyon, UCBL, CNRS, Inria
LIP UMR 5668
hugues. van, www.eeensesseseeeeee..ee.aeeeeem dess-yon.f
titouan. vayer
Rémi Flamary
Nicolas Courty
École polytechnique, IP Paris, CNRS
Université Bretagne Sud, CNRS IRISA UMR 6074
CMAP UMR 7641
remi -llmaryyyplllytcciiiuuuueeeeee edu
nicolas, nicolas.cuurt// couttyirria.
Abstract
Many approaches in machine learning rely on weighted graph to encode the similarities between samples a dataset. Entropic affinities (EAs), which are notably used in the popular Dimensionality Reduction (DR) algorithm t-SNE, particular instances of such graphs. To ensure robustness to heterogeneous sampling densities, EAs assign kernel bandwidth parameter to every sample in such way that the entropy of each row in the affinity matrix is kept constant at a specific value, whose exponential is known perplexity. EAs are inherently asymmetric and row-wise stochastic, but they are used in DR approaches after undergoing heuristic symmetrization methods that violate both the row-wise constant entropy and stochasticity properties. In this work, we uncover a novel characterization of EA an optimal transport problem, allowing natural symmetrization that can be computed efficiently using dual ascent. The corresponding novel affinity matrix derives advantages from symmetric doubly stochastic normalization terms clustering performance, while also effectively controlling the entropy of each row thus making particularly robust to varying noise levels. Following, we present a new DR algorithm, SNEkhorn, that leverages this new affinity matrix. We show its clear superiority to existing approaches with several indicators on both synthetic and real-world datasets.
Introduction
Exploring and analyzing high-dimensional data a core problem of data science that requires building low-dimensional and interpretable representations of the data through dimensionality reduction (DR). Ideally, these representations should preserve the data structure by mimicking, in the reduced representation space (called latent space), notion of similarity between samples. We call affinity the weight matrix of a graph that encodes this similarity. has positive entries and the higher the weight in position (i. the higher the similarity or proximity between samples i and j. Seminal approaches relying on affinities include Laplacian eigenmaps 22). spectral clustering [50] and semi-supervised learning (55). Numerous methods be employed to construct such affinities. A common choice is to use kernel (e.g., Gaussian) derived from a distance matrix normalized by a bandwidth parameter that usually has a large influence on the outcome of the algorithm. Indeed, excessively small kernel
37th Conference on Neural Information Processing Systems (NcurIPS 2023).

---

bandwidth can result in solely capturing the positions of closest neighbors, at the expense of large- scale dependencies. Inversely, setting too large bandwidth blurs information about close-range pairwise relations. Ideally, one should select a different bandwidth for each point to accommodate varying sampling densities and noise levels. One approach is to compute the bandwidth of a point based on the distance from its k-th nearest neighbor [53]. However, this method fails to consider the entire distribution of distances. In general, selecting appropriate kernel bandwidths can be a laborious task, and many practitioners resort to greedy search methods. This can be limiting in some settings, particularly when dealing with large sample sizes.
Entropic Affinities and SNE/t-SNE. Entropic affinities (EAs) were first introduced in the seminal paper Stochastic Neighbor Embedding (SNE) [16]. It consists in normalizing each row of distance matrix by bandwidth parameter E, such that the distribution associated with each row of the corresponding stochastic (i.e, row-normalized) Gaussian affinity has a fixed entropy. The value of this entropy, whose exponential is called the perplexity, is then the only hyperparameter left to tune and has an intuitive interpretation as the number of effective neighbors of each point [49]. EAs are notoriously used to encode pairwise relations in high-dimensional space for the DR algorithm I-SNE [46), among other DR methods including [6]. t-SNE is increasingly popular in many applied fields 220 32] mostly due to its ability to represent clusters in the data p7 Nonetheless, one major flaw of EAs is that they are inherently directed and often require post- prose-poesssing symmetrization.
Doubly Stochastic Affinities. Doubly stochastic (DS) affinities are negative matrices whose rows and columns have unit E norm. In many applications, it has been demonstrated that DS affinity normalization (i.e., determining the nearest DS matrix to given affinity matrix) offers numerous benefits. First, it can seen as a relaxation of k-means [51] and it is well-established that it enhances spectral clustering performances [10 [22 D. Additionally, DS matrices present the benefit of being invariant to the various Laplacian normalizations (50). Recent observations indicate that the DS projection of the Gaussian kernel under the KL geometry is more resilient to heteroscedastic noise compared to its stochastic counterpart (23) It also offers a more natural analog to the heat kernel 330 These properties have led to a growing interest in DS affinities, with their use expanding to various applications such as smoothing filters B3 subspace clustering (25 and transformers [41].
Contributions. In this work, we study the missing link between EAs, which are easy to tune and adaptable to data with heterogeneous density, and DS affinities which have interesting properties in practical applications as aforementioned. Our main contributions are as follows. We uncover the convex optimization problem that underpins classical entropic affinities, exhibiting novel links with regulae-re Optimal Transport (OT) (Section We then propose in Section 3.2 principled symmetrization of entropic affinities. The latter enables controlling the entropy in each point, unlike t-SNE's post-processing symmetrization, and produces genuinely doubly stochastic affinity. We show how to compute this new affinity efficiently using a dual ascent algorithm. In Section
4 we introduce SNEkhorn: DR algorithm that couples this new symmetric entropic affinity with a doubly stochastic kernel in the low-dimensional embedding space, without sphere concentration issue [29]. We finally showcase the benefits of symmetric entropic affinities on variety of applications in Section 5 including spectral clustering and DR experiments on datasets ranging from images genomies data.
Notations. [r:]] denotes the set (1,..., m). exp and log applied to vectors/matrices are taken element-
wise. 1 (1, ..., ..TT..)TTTT the vector of 1. (. is the standard inner product for matrices/vectors. S is the space of n x symmetric matrices. P.. denotes the i-th row of a matrix P. @ (resp. ( stands for clement-wise multiplication (resp. division) between vectors/matrices. For a, B E E" a b B e Rax= is (o, + 35)ay. The entropy of p E TR i H()) = - E, Pt (log((); - 1) = - -(p, log p - 1). The Kullback-Leibles divergence between two matrices P, Q with nonnegative entries such that (log( Pi Q
Qi 0 = P.j is KL(P)Q) Eig Py
11) log (P 2 Q) 100117
2 Entropic Affinities, Dimensionality Reduction and Optimal Transport
Given a dataset X E R nf samples in dimension p. most DR algorithms compute represen- tation of X in a diver-dimensional latent space Z E xx with q < p that faithfully captures and represents pairwise dependencies between the samples (or rows) in X. This is generally achieved
With the convention 0 log 01 0.

---

by optimizing Z such that the corresponding affinity matrix matches another affinity matrix defined from X. These affinities are constructed from matrix C R that encodes a notion of "distance" between the samples, e.g., the squared Euclidean distance Ct. = |X -x;:12 2 or more generally any cost matrix C E D == {C E R : C = C' and Cij = 0 - i= j}. A commonly used option is the Gaussian affinity that is obtained by performing row-wise normalization of the kernel exp(-C/e), where e > 0 is the bandwidth parameter.
Entropic Affinities (EAs). Another frequently used approach to generate affinities from C C D is to employ entropic affinities [16). The main idea is to consider adaptive kernel bandwidths (et) [r(r) to capture finer structures in the data compared to constant bandwidths [47). Indeed, EAs rescale distances to account for the varying density across regions of the dataset. Given & E In - II] the goal
of EAs is to build Gaussian Markov chain transition matrix Pe with prescribed entropy as
exp (xC(g/g//E))
Vi, Vj. P. E,
(EA)
with of E RR s.t. = =ogg)
The hyperparameter S, which is also known as perplexity, can be interpreted as the effective number of neighbors for each data point [49] Indeed, perplexity of Ç means that each row of P (which is discrete probability since P is row-wise stochastic) has the same entropy as a uniform distribution over & neighbors. Therefore, it provides the practitioner with interpretable parameter specifying which scale of dependencies the affinity matrix should faithfully capture. In practice, a root-finding al- gorithm is used to find the bandwidth parameters uef) that satisfy the constraints (49) Hereafter, with a slight abuse of language, we call H(P.))-1 the perplexity of the point i.
Dimension Reduction with SNE/t-SNE. One of the main applications of EAs is the DR algorithm SNE [6] We denote by Cx = (I X,:113) ij and Cz (IIZ: Zz:11B);; the cost matrices derived from the rows (i.e., the samples) of X and Z respectively. SNE focuses on minimizing in the latent coordinates Z R the objective
where Pe solves (EA) with cost Cx and Qz)is = exp( exp(-(Czi) //( exp(- - (epp--CC]at) [Cz]aa)) In the seminal paper [46]. a newer proposal for a symmetric version was presented, which has since replaced SNE in practical applications. Given sym- metric normalization for the similarities in latent space [Qzlis it consists in solving
Σe, exp(-ICz)e))
min ZER
KL PL Qz) where
P
(P* + PP
(Symmettic-SNE)
In other words, the affinity matrix P is the Euclidean projection of P on the space of symmetric matrices S: Pe = Projs (Pa) = arg millpes P P' (see Appendix A.I). Instead of the
Gaussian kernel, the popular extension t-SNE [46] considers different distribution in the latent space [åz] = (1 + (C + [Cz]ut) In this formulation, Qz is a joint Student t-distribution that accounts for crowding effects: relatively small distance in high-dimensional space can be accurately represented by a significantly greater distance in the low-dimensional space.
Considering symmetric similarities is appeal- ing since the proximity between two points is inherently symmetric. Nonetheless, the Eu- clidean projection in (Symmetric-SNe) does not preserve the consiruction of entropic affinities. In particular, pe is not stochastic in general and H(P:) s: * (log + 1) thus the entropy associated with each point is no longer controlled after symmetrization (see the bottom left plot of Figure Ib. This is ar-
I-SNE
I-SNEkborn
pst
guably one of the main drawbacks of the ap- proach. By contrast, the pse affinity that will be introduced in Section 3 can accurately set the entropy in each point to the desired value log& + 1. As shown in Figure this leads to more faithful embeddings with better separa- tion of the classes when combined with the t-SNEkhom algorithm (Section 4).
Sample
Sample
embeddings with silhou- Figure 1: Top: COIL 34
ette scores produced by t-SNE and t-SNEkhorn (our method introduced in Section 4 for & = 30. Bottom: e H(P..)-1 (perplexity) for each point i.
3

---

Symmetric Entropy-Constrained Optimal Transport. Entropy-regularizn OT [37 and its con- nection to affinity matrices are crucial components in our solution. the special case of uniform marginals, and for >> 0, entropic OT computes the minimum of P - (P, C) E. H(P..) over the space of doubly stochastic matrices {P E _ P1 = P 1]. The optimal solution is the unique doubly stochastic matrix pds of the form Pdc diagí u K diag(v) where K exp( (ppC-c/)) is the Gibbs energy derived from C and u, V are positive vectors that can be found with the celebrated Sinkhorn- Knopp''' algorithm [8 42 Interestingly, when the cost C is symmetric (e.g., C E D) can take u = v 17 Section 5.2] so that the unique optimal solution is itself symmetric and writes
p = exp ((f + f - C)vv) where f E R
(DS)
In this case, by relying on convex duality as detailed in Appendix A.2 an equivalent formulation for the symmetric entropic OT problem is
min PER
(P.C) s.t.
and Σ H(P)))
(EOT)
where 0 S n (logn 1) is constraint the global entropy Σ, H.)).) of the OT plan P which happens to be saturated at optimum (Appendix [A.)). This constrained formulation of symmetric entropic OT will provide new insights into entropic affinities, detailed in the next sections.
3 Symmetric Entropic Affinities
In this section, we present our first major contribution: symmetric entropic affinities. We begin by providing new perspective on EAs through the introduction of an equivalent convex problem.
3.1
Entropic Affinities Entropic Optimal Transport
We introduce the following set of matrices with row-wise stochasticity and entropy constraints:
Hgg = [P E g s.t. P1 = and Vi, H(P))) ) 2 log E + 1}
(I)
This space is convex since p E R. + H()) is concave, thus its superlevel set is convex. In contrast to the entropic constraints utilized in standard entropic optimal transport which set a lower-bound on the global entropy, as demonstrated in the formulation (EOT). Hg imposes a constraint on the entropy of each row of the matrix P. Our first contribution is to prove that EAs can be computed by solving specific problem involving Hę (see Appendix A the proof).
Proposition 1. Let C E Raxx without constant rows. Then Pe solves the entropic affinity problem EA with cost C if and only if p is the unique solution of the convex problem
min PEH
(..)).
(EA as OT)
Interestingly, this result shows that EAs boil down to minimizing a transport objective with cost C and row-wise entropy constraints He where & is the desired perplexity. As such, (EA as OT) can be seen as a specific semi-relaxed OT problem (39 14) (i.e., without the second constraint on the marginal P PT1 1) but with entropic constraints on the rows of P. We also show that the optimal solution P of (EA as OT) has satarated entropy ie., Vi, H(P:) Iog In other words, relaxing the equality constraint in (EA) as an inequality constraint in P e He does not affect the solution while it allows reformulating entropic affinity as a convex optimization problem. To the best of our knowledge, this connection between OT and entropic affinities is novel and is an essential key to the method proposed in the next section.
The kernel bandwidth parameter e from the original formulation of entropic affinities (EA) is the Lagrange dual variable associated with the entropy constraint in (EA as OT). Hence computing =* in (EA) exactly corresponds to solving the dual problem of (EA as OT).
Remark 2.
Remark 3. Let K, then P Proj KT. (K,) arg minpe KL(P|K,). Therefore P can be seen KL Bregman projection 33 of Gaussian kernel onto He. Hence the input matrix in (Symmetric-SNE) is P KL Projšs & (Proj (K,)) which corresponds to surprising mixture of KL and orthogonal projections.
CCoo). As shown in Appendix A.5 if e* solves (EA and < min(s"),

---

pil
Perpheity
P
Perpleiity
Laplacian Eigenvalues
p
0170
p
15
15
0180
10
10
16
16
Figure 2: Samples from a mixture of three Gaussians with varying standard deviations. The edges' strength is proportional to the weights in the affinities pds (DS) and pse (SEA) computed with & = 5 (for pdb, is the average perplexity such that Σ, H(P... = E, H(P:)) }). Points' color represents the perplexity (H(P.)-1, Right plot: smallest eigenvalues of the Laplacian for the two affinities.
3.2 Symmetric Entropic Affinity Formulation
Based the previous formulation version of EAs that enables kceping the entropy associated with each row (or equivalently column) to the desired value of log & + 1 while producing symmetric doubly stochastic affinity matrix. Our strategy to enforce symmetry through an additional constraint in (EA as OT), in a similar fashion (OOT. More precisely we consider the convex optimization problem
propose symmetric entropic affinities: symmetric
min (P.C).
(SEA)
PER
where we recall that S is the set of n x n symmetric matrices. Note that for any & S n 1, 111 Hg S hence the He n S is non-empty and convex set. We first detail some important properties of problem (SEA) (the proofs of the following results can be found in Appendix A.4).
Proposition (Saturation of the entropies). Let C with zero diagonal, then SEA) with cost C has a unique solution that denote by p If moreover C CE? D, then for at least indices i E [n] the solution satisfies H(P:.) = log { +
saturated entropies i.e., the corresponding n 1 points have exactly a perplexity of ç. In practice, with the algorithmic solution detailed below, have observed that all n entropies are saturated. Therefore, we believe that this proposition can be extended with a few more assumptions on C. Accordingly, problem (SEA) allows accurate control over the wose entropies while providing a ssmmmetric doubly stochastic matrix, unlike P defined in Syym as summarized in Table 1 In the sequel, denote by (P(P) (H(P::)) the vector of rowe entropies of P. We rely on the following result to compute p
In other words, the unique solution se has at least
Proposition (Solving for SEA). Let C E (A, 11-P1 P1) and g(>, X) minl PER".""S C(P,y,X) Strong duality holds for (SEA). Moreover, let Y X* E argmax,20.x 4(>, X) be the optimal dual variables respectively associated with the entropy and marginal constraints. Then, for at least n 1 indices E [m], > O. When Vi E [m]. Y 0 then (Pee) (log 1)1 and pse has the form
7, (P,, + (r, (log & + 1)1 -OT-HHA H,
a exp ((X* + X- - (1")
By defining the symmetric matrix P(7, A) exp (A A 2C) 2 (y ① 7)), we prove that, when 0. minpes L(P, Y, A) has unique solution given by P(y,A) which implies g(v, A) C(P(r, A), 7, A). Thus the proposition shows that when Y > 0, p P(y", X+) where Y X
solve the following convex problem (as maximization of concave objective)
max ((P(Y,A),T. 7,0
(Dual-SEA)
5

---

Consequently, to find a'd we solve the problem (Dual-SEA). Although the form of pe pre- sented in Proposition{Sis only valid when is
Table 1: Properties of Pe, Pe, pds and pae
AFFINITY MATRIX
P
P
pds
psc
positive and we have only proved it for n indices, we emphasize that if (Dual-SEA) has finite solution, then equal to aae Indeed in this case the solution satisfies the KKT system associated with (SEA).
REFERENCE
[6
146
29
SEAA
P PT
P1 P 1= H. P
(log 5 + 1)1
Numerical optimization. The dual problem (Dall-EEA) antees through dual ascent approach with closed-form gradients (using e.g., SGD, BFGS [28] or ADAM (8]. At each gradient step, one can compute the current estimate P X) while the gradients of the loss w.t.. " and À are given respectively by the constraints (log & 1)1 (P 2)) and 1 - P(y, A)1 (see e.g., [4] Proposition 6.1.1J). Concerning time complexity, each step can be performed with (n²)) algebraic operations. From practical perspective, found that using change of variable y f- x- and optimize 7 E I leads to enhanced numerical stability.
concave and can solved with guar-
Remark 6. In the same spirit as Remark one can express p as a KL projection of K, exp(-C/on Indeed, we show in Appendix A.5 that if 0 p ProjHers (Ko,)
min; 7 then
Comparison between pds and pse. In Figure 2 we illustrate the ability of our proposed affinity p= to adapt varying noise levels. In the OT problem that we consider, each sample is given a mass of one that is distributed over its neighbors (including itself since self- are allowed). For each sample, we refer to the entropy of the distribution over its neighbors as the spreading its mass. One can notice that for pds (DS) (OT problem with global entropy constraint (EOT)) the samples do not spread their mass evenly depending on the density around them. On the contrary, the per-row entropy constraints pae force equal spreading among samples. This can have benefits, particularly for clustering, illustrated in the rightmost plot, which shows eigenvalues of the associated Laplacian matrices (recall that the number of connected components equals the dimension the null space of its Laplacian (7]). As can be seen, pds results in many unwanted clusters, unlike Pde, which is robust to varying noise levels (its Laplacian matrix has only vanishing eigenvalues). We further illustrate this phenomenon on Figure with varying noise levels.
Adjusted rand index (ARI)
6.
P
do
Noise level
Figure 3: ARI spectral clustering the example of three Gaussian clusters with variances: o², 20 and 30² (as in Figure 2).
4 Optimal Transport for Dimension Reduction with SNEkhorn
In this section, we build upon symmetric entropic affinities to introduce SNEkhom,, a new DR algorithm that fully benefits from the advantages of doubly stochastic affinities.
SNEkhorn's objective. Our proposed method relies on doubly stochastic affinity matrices to capture the dependencies among the samples in both input and latent spaces. The KT. divergence, which is the central criterion in most popular DR methods (44, used to measure the discrepancy between the
two affinities. As detailed in sections 2 and 3 ead computed using the cost [Cx]i = || X,:11
corrects for heterogeneity in the input data density by imposing point-wise entropy constraints. As we do not need such correction for embedding coordinates Z since they must be optimized, opt for the standard affinity (DS) built as an OT transport plan with global entropy constraint (EOT). This OT plan be efficiently computed using Sinkhorn's 's algorithm. More precisely, propose the optimization problem
min zER
(P*
(SNEkhom)
where Q exp (fz + fz Cz) stands for the (DS) affinity computed with cost [Czlip IZ, Zj 11 and fz is the optimal dual variable found by Sinkhorn's algorithm. We set the bandwidth to 1 in Qd similarly to [46] the bandwidth in the low dimensional space only affects the scales

---

P
Symmetrie
p
SNEXboo
04
04
1-01
Po
Figure 5: From left to right: entries Pe (Symmetric- Syymetric-SNE) SNE and associated embeddings generated
using P Then p (EAA) matrix and associated SNEkhorn embeddings. Perplexity & = 30.
of the embeddings and not their shape. Keeping only the terms that depend Z and relying on the double stochasticity of pse, the objective in (SNEkhorn) can be expressed as (P Cz) 2ffz, 1).
Heavy-tailed kernel in latent space. Since it is well known that heavy-tailed kernels can be beneficial DR (1I, we propose an extension called t-SNEkhorn that simply amounts to computing a doubly stochastic student-t kernel in the low-dimensional space. With our construction, it corresponds to choosing the cost [Czlis log(1 + IZ Zj: I12) instead of IZ. Zz:112.
Inference. This new DR objective involves computing doubly stochastic normalization for each update of Z. Interestingly, compute the optimal dual variable fz in Q we leverage well- conditioned Sinkhorn fixed point iteration setting:
[19 13,, which converges extremely fast in the symmetric
1 2
Vi, [[z]
[[z] log Σ exp ([fz]= [Cz]xa)
(Sinkhorn)
the right side of Figure we plot ||Q 1 20 function of (Sinkhorn) iterations for toy example presented Section 5 In most practical cases, we found that about 10 iterations were enough to reach sufficiently small error. Z is updated through gradient descent with gradients obtained by performing back- propagation through the Sinkhorn iterations. These iterations be further accelerated with warm start strategy by plugging the last lz to initialize the current one.
minz KT. Qz)
Q
1 1| S tol
Her.. needed
10
tol. level
Figure 4: Left: SNEkhorn embedding on the simulated data of Section 5 using Qz instead Q with & = 30. Right: number of iterations needed to achieve ||Qu21 || tol with (Sinkhorn).
Related work. Using doubly stochastic affinities for SNE has been proposed in (29). with two key differ- ences from our work. First, they do not consider EAs and resort to d (SS) This affinity, unlike pae, is not adaptive to the data heterogeneous density (as illusrated in Figure 2) Second, they use the affinity Qz in the low-dimensional space and illustrate empirically that matching the latter with doubly stochastic matrix (e.g., pds psc) can sometimes impose spherical constraints on the embedding Z. This is detrimental for projections onto 2D flat space (typical use case of DR) where embeddings tend to form circles. This can verified on the left side of Figure 4 In contrast, in SNEkhorn, the latent affinity is also doubly stochastic so that latent coordinates Z are not subject to spherical constraints anymore. The corresponding SNEkhorn embedding is shown in Figure (bottom right).
Numerical experiments
This section aims to illustrate the performances of the proposed affinity matrix pee (SEA) and DR method SNEkhorn faithfully representing dependencies and clusters in low dimensions. First, we showcase the relevance of our approach on simple synthetic dataset with heteroscedastic noise. Then,
7

---

we evaluate the spectral clustering performances of symmetric entropic affinities before benchmarking t-SNEkhomn with t-SNE and UMAP [1] on real-world images and genomics datasets H
Simulated data. We consider the toy dataset with heteroscedastic noise from 223 It consists of sampling uniformly two vectors P1 and P2 in the 10 -omedimensional probability simplex. n. 10
samples are then generated as X, = / (2 Tig)
where
M(1000, P1), 1 15 < 500
x;
M(1000, p2),
501 S i S 750
M2000,, p2),
751 S 1000.
where M stands for the multinomial distribution. The goal of the task is test the robustness to heteroscedastic noise. Indeed, points generated using P2 exhibit different levels of noise due to various numbers of multinomial trials to form an estimation of P2. This typically occurs in real- world scenarios when the same entity is measured using different experimental setups thus creating heterogeneous technical noise levels (e.g., in single-cell sequencing 20)). This phenomenon is known as batch effect (3). In Figure 5] we show that, unlike P (Symmetric-SNE). pse (SEA) manages to properly filter the noise (top row) to discriminate between samples generated by P1 and P2, and represent these two clusters separately in the embedding space (bottom row). In contrast, Pe and SNE are misled by the batch effect. This shows that P doesn' 't fully benefit from the adaptivity of EAs due to poor normalization and symmetrization. This phenomenon partly explains the superiority of SNEkhorn and t-SNEkhorn over current approaches on real-world datasets as illustrated below.
Real-world datasets. We then experi- ment with various labeled classification datasets including images and genomic data. For images, we use COIL 20 34 OLIVETTI faces 112. UMNIST 155
Table 2: ART (x clustering scores on genomies.
DATA SET
P
pda
pst
Pse
LIVER (14520)
75.8
75.8
84.9
80.8
85.9
BREAST (20947)
30.0
30.0
26.5
23.5
28.5
and CIFAR 10 (22. For CIFAR, we ex- periment with features obtained from the last hidden layer of pre-trained ResNet [38 while for the other three datasets, we take as input the raw pixel data. Re- garding genomies data, consider the Curated Microarray Database (CuMiDa)
LEUKEMIA (28497)
43.7
44.1
19.7
12.5
50.6
COLORECTAL (44076)
95.9
95.9
93.9
95.9
95.9
LIVER (76427)
76.7
76.7
83.3
81.1
81.1
BREAST (45827)
43.6
53.8
74.7
71.5
77.0
COLORECTAL (21510)
57.6
57.6
54.7
94.0
79.3
RENAL (53757)
47.6
47.6
49.5
49.5
49.5
PROSTATE (6919)
12.0
13.0
13.2
16.3
17.4
[ made of microarray datasets for var- ious types of cancer, as well as the pre- processed SNAREseq (chromatin acces- sibility) and scGEM (gene expression) datasets used in 9 For CuMiDa, we retain the datasets with most samples. For all the datasets, when the data dimension exceeds 50 we apply pre-processing step of PCA in dimension 50, as usually done in practice [46). In the following experiments, when specified the hyperparameters are set to the value leading to the best average score on five different seeds with grid-search. For perplexity parameters, we test all multiples of 10 in the interval [10, min(n, 300] where n is the number of samples in the dataset. We use the same grid for the & of the self-tuning affinity pst (53) and for the n_ neineibbors parameter of UMAP. For scalar bandwidths, we consider powers of 10 such that the corresponding affinities' average perplexity belongs to the perplexity range.
THROAT (42743)
9.29
9.29
11.4
11.8
44.2
SCGEM SNARESEQ
57.3
58.5
74.8
69.9
71.6
8.89
9.95
46.3
55.4
96.6
Spectral Clustering. Building the strong connections between spectral clustering mechanisms and t-SNE 44 27] we first consider spectral clustering tasks to evaluate the affinity matrix p (EAA) and compare against P (Symmetric-SNe SNE). We also consider two versions of the Gaussian affinity with scalar bandwidth K
the symmetrized row-stochastic p = Proj's (Pe)
where P is K normalized by row and dds (DS). We also consider the adaptive Self-Tuning pst affinity from (53) which relies
adaptive bandwidth corresponding the distance from the k-th nearest neighbor of each point. We use the spectral clustering implementation of scikit- -learn (36) with default parameters which uses the unnormalized graph Laplacian. We measure the quality of clustering using the Adjusted Rand Index (ARI). Looking at both Table and Figure X one can notice that, in general, symmetric entropic affinities yield better results than usual entropic affinities with significant improvements in some datasets (e.g., throat microarray and SNAREseq). Overall
Our code is available at http:::///hthaaa...a.o/otott//t/ wow./PPyytoooOO////tttrro.

---

0.70
0.33
Adputed Bost Sooe
0.00
0.45
0:29
0.91
0023
0.00
D.71
0.35
0.20
0.30
0.50
0.25
100
210
100)
300
1.00
.000
2001
Pirpiestty
Peeplecity
Porplecty
P
p
p
pos
p
Figure 6: ARI spectral clustering score function of the perplexity parameter for image datasets.
Table 3: Scores for the UMAP, t-SNE and t-SNEkhorn embeddings.
Silhouette (x100)
Trustworthiness (x (x100)
UMAP
(-SNE
L-SNEkhorn
UMAP
I-SNE
L-SNEkhorn
COIL
20.4 3.3
30.7 6.9
52.3
99.6 1
99.6 0.1
99.9
OLIVETTI
6.4 = 4.2
4.1
15.7 1 2.2
96.5 1 1.3
96.2
98.0 0.4
UMNIST
2.7
1.5
2.94
93.0
9.2
99.8 1 0.1
CIFAR
13.6
18.3
1.3
90.2 0.8
90.1
92.4 0.3
Liver (14520)
49.7 1.3
61.1 0.7
89.2 0.7
9..4
92.3 + 0.3
Breast (1997)
0.8
29.0
31.2
9.3
93.2 0.4
Leukemia (28477)
0.7
26.2 2.7
1.1
94.3 9.8
0.5 0.5
Colorectal (44m7)
67.6 2.2
6.5
+
93.2
93.7
94.3 + 0.6
Liver (76427)
4.3
38.
2.5
0.4
8..4
1.0
Breast (48827)
3.3
1.9
4.5
9.4
0.2
94.7 + 0.3
Colorectal (21510)
1.3
4.6
35.1
0.9 0.7
88.2 08.2 0.7
Renal (55757)
1.5
4.3
47.8
9.2
94.0 0.2
0.2
Prostate (6919)
2.4
0.2
9.1
8.2 1.8
73.1 13.1000 0.5
Throat (43743)
2.4
2.30
0.1
1.3
0.8
86.8 1.0
scGEM
26.9
33.0
0.7
95.0 1.3
96.2
96.8 0.3
SNAREseq
6.8 6.8=64 6.0
35.8 + 5.2
67.9
93.1 2.8
9.11
99.2 0.1
pse outperforms all the other affinities in is crucial. Figure (ouure/]alss also shows that this superiority is verified for the whole range of perplexities. This can be attributed doubly stochastic normalization in terms of clustering and of EAs in terms of adaptivity. In the next experiment, we show that these advantages translate into better elustering and neighborhood retrieval at the embedding level when running SNEkhor..
12 datasets. This shows that the adaptivity of EAs
the fact that symmetric entropic affinities combine the advantages of
Dimension Reduction. To guarantee fair parison, we implemented not only SNEkhorn, but also t-SNE and UMAP in PyTorch (35] Note that UMAP also relies on adaptive affini- ties but sets the degree of each node (related to the hyperparameter Deneernnniiiiiiie i_neighbors which plays similar role to the perplexity) rather than the entropy. All models were optimized using ADAM 18 with default parameters and the same stop- ping criterion: the algorithm stops whenever the relative variation of the loss becomes smaller than 10-5 For each run, we draw independent (0, 1) coordinates and use this same matrix initialize all the methods that we wish to compare. To evaluate the embeddings' quality, make use of the silhouette [40] and trustworthiness [48
+-NNE
-33
-3
-0
50
-35
Figure SNAREseq embeddings produced by t-SNE and t-SNEkhorn with 50.

---

scores from scikit-lears [36] with default parameters. While the former relies on class labels, the latter measures the agreement between the neighborhoods in imput and output spaces, thus giving two complementary metrics to properly evaluate the embeddings. The results, presented in Table 3 demonstrate the notable superiority of t-SNEkhorn compared to the commonly used t-SNE and UMAP algorithms. A sensitivity analysis on perplexity can also be found in Appendix B Across the 16 datasets examined, t-SNEkhorn almost consistently outperformed the others, achieving the highest silhouette score 15 datasets and the highest trustworthiness score on 12 datasets. To visually assess the quality of the embeddings, we provide SNAREseq embeddings in Figure 7 Notably, one can notice that the use of t-SNEkhorm results in improved class separation compared to t-SNE.
6 Conclusion
We have introduced a new principled and efficient method for constructing symmetric entropic affinities. Unlike the current formulation that enforces symmetry through orthogonal projection, our approach allows control over the entropy in each point thus achieving entropic affinities' primary goal. Additionally, it produces DS-normalized affinity and thus benefits from the well- advantages of this normalization. Our affinity takes as input the same perplexity parameter as EAs and can thus be used with little hassle for practitioners. We demonstrate experimentally that both our affinity and DR algorithm (SNEkhorn), leveraging a doubly stochastic kernel in the latent space, achieve substantial improvements over existing approaches.
Note that in the present work, we do not address the issue of large-scale dependencies that are not faithfully represented in the low-dimensional space [44). The latter shall be treated in future works. Among other promising research directions, one could focus on building multi-scale versions of
symmetric entropic affinities [24] as well as fast approximations for SNEkhorn forces by adapting methods (26) to the doubly stochastic setting. could
egg., Barnes-Hul (45] or interpolation-
also be interesting to use SEAs in order to study the training dynamics of transformers [54].
Acknowledgments
The authors are grateful to Mathurin Massias, Jean Feydy and Aurélien Garivier for insightful discussions. This project was supported in part by the ANR projects AllegroAssai ANR-19-CHIA- 0009, SingleStatOmmies ANR-18-CE45-00 also supported by the ACADEMICS grant of the IDEXLYON, project of the Université de Lyon, PIA operated by ANR-16-IDEX-0005.
and OTTOPIA ANR-20-CHIA--00000 This work was
References
[1] Mario Beauchemin. On affinity matrix normalization for graph cuts and spectral clustering. Puttern Recognition Letters, 68:90-96, 2015.
[2] Mikhail Belkin and Partha Niyogi. Laplacian eigenmaps for dimensionality reduction and data representation. Neural computation, 15(6): 1373-1637 2003.
[3] Jean-David Benamou, Guillaume Carlier, Marco Cuturi, Luca Nenna, and Gabriel Peyré. Itera- tive bregman projections for regularized transportation problems. SIAM Journal on Scientific Computing, 37()::1111-A11 1138,,, 2015.
[4] Dimitri P Bertsekas. Nonlinear programming. Journal of the Operational Research Society, 48(3):334-334, 1997.
[5] T. Tony Cai and Rong Ma. Theoretical foundations of t-sne for visualizing high-dimensional clustered data. Journal of Machine Learning Research (JMLR), 23(301): 2022.
[6] Miguel Carreira-Perpinán. The elastic embedding algorithm for dimensionality reduction. In International Conference on Machine Learning (ICML), volume 10, pages 167-174, 2010.
[7] Fan RK Chung. Spectral graph theory, volume 92. American Mathematical Soc., 1997.
[8] Marco Cuturi. Sinkhorn distances: Lightspeed computation of optimal transport. Neural Information Processing Systems (NeurIPS), 26, 2013.
10

---

[9] Pinar Demetci, Rebecca Santorella, Bjöm Sandstede, William Stafford Noble, and Ritambhara Singh. wrssev-waseerstein optimal transport align single-cell multi-omic data. bioRxiv. 2020.
[10] Tianjiao Ding, Derek Lim, Rene Vidal, and Benjamin D Haeffele. Understanding doubly stochastic clustering. In International Conference on Machine Learning (ICML), 2022.
[11] Bruno César Feltes, Eduardo Bassani Chandelier, Bruno lochins Grisci, and Márcio Dorn. Cumida: An extensively curated microarray database for benchmarking and testing of machine learning approaches in cancer research. Journal of Computational Biology, 2664):376-386, 2019. PMID: 30789283.
[12] Samaria Ferdinando and Harter Andy. Parameterisation of stochastic model for human face identification, 1994.
[13] Jean Feydy, Thibault Séjourné, François-Xave Xavier Vialard, Shun-ichi Amari, Alain Trouvé, and Gabriel Peyré. Interpolating between optimal transport and mmd using sinkhorn divergences. In International Conference on Artificial Intelligence and Statistics (AISTATS), pages 2681-2690. PMLR, 2019.
[14] Rémi Flamary, Cédric Févotte, Nicolas Courty, and Valentin Emiya. Optimal spectral transporta- tion with application music transcription. Neural Information Processing Systems (NeurIPS), 29, 2016.
[15] Daniel B Graham and Nigel M Allinson. Characterising virtual cigensignatures for general purpose face recognition. Face recognition: from theory applications, pages 446-456, 1998.
[16] Geoffrey E Hinton and Sam Roweis. Stochastic neighbor embedding. Neural Information Processing Systems (NeurIPS), 15, 2002.
[17] Martin Idel. A review of matrix scaling and sinkhorn's normal form for matrices and positive maps. arXiv preprint arXiv: 169iv:1609.06349 2016.
[18] Diederik P Kingma and Jimmy Ba. Adam: A method for stochastic optimization. arXiv preprint 1rXiv:1412.6980, 2014.
[19] Philip Knight, Daniel Ruiz, and Bora Uçar. symmetry preserving algorithm for matrix scaling. SIAM journal on Matrix Analysis and Applications, 35(3):931-955. 2014.
[20] Dmitry Kobak and Philipp Berens. The art of using t-sne for single-cell transcriptomics. Nature communications, 10(1):1-14, 2019.
[21] Dmitry Kobak, George Linderman, Stefan Steinerberger, Yuval Kluger, and Philipp Berens. Heavy-tailed kernels reveal finer cluster structure in t-sne visualisations. In Joint European Conference on Machine Learning and Knowledge Discovery in Databases, pages 124-139. Springer, 2020.
[22] Alex Krizhevsky, Geoffrey Hinton, et al. Learning multiple layers of features from tiny images. 2009.
[23] Boris Landa, Ronald R Coifman, and Yuval Kluger. Doubly stochastic normalization of the gaussian kernel robust to heteroskedastic noise. SIAM journal on mathematics of data science, 3(1):388-413, 2021.
[24] John A Lee, Diego H Peluffo-Ordóñezz and Michel Verleysen. Multi-scale similarities in stochastic neighbour embedding: Reducing dimensionality while preserving both local and global structure. Neurocomputing, 169:246-261 2015.
[25] Derek Lim. René Vidal, and Benjamin D Haeffele. Doubly stochastic subspace clustering. arXiv preprint arXiv:2011. 148599, 2020.
[26] George C Linderman, Manas Rachh, Jeremy G Hoskins, Stefan Steinerberger, and Yuval Kluger. Fast interpolation t-sne for improved visualization of single-cell rna-seq data. Nature methods, 16(3):243-245, 2019.

---

[27] George C Linderman and Stefan Steinerberger. Clustering with t-sne, provably. SIAM Journal on Mathematics of Data Science, 1(2):313-332, 2019.
[28] Dong C Liu and Jorge Nocedal. On the limited memory bfgs method for large scale optimization. Mathematical progrumming, 45(1-3):503-528, 1989.
[29] Yao Lu, Jukka Corander, and Zhirong Yang. Doubly stochastic neighbor embedding on spheres. Puttern Recognition Letters, 128: 128:100-106, 2019.
[30] Nicholas F Marshall and Ronald R Coifman. Manifold learning with bi-stochastic kernels. IMA Journal of Applied Mathematics, 84(3):455 08255-44,,22 2019.
[31] Leland Mclnnes, John Healy, and James Melville. Umap: Uniform manifold approximation and projection for dimension reduction. arXiv preprint arXiv:1802.03440 2018.
[32] Binu Melit Devassy, Sony George, and Peter Nussbaum. Unsupervised clustering of hyperspec- tral paper data using t-sne. Journal of Imaging, 6(5):29, 2020.
[33] Peyman Milanfar. Symmetrizing smoothing filters. SIAM Journal on Imaging Sciences,
6(1 6(3)2263,, 2013.
[34] Sameer A Nene, Shree K Nayar, Hiroshi Murase, et al. Columbia object image library (coil-20). 1996.
[35] Adam Paszke, Sam Gross, Soumith Chintala, Gregory Chanan, Edward Yang, Zachary De Vito, Zeming Lin, Alban Desmaison, Luca Antiga, and Adam Lerer. Automatic differentiation in pytorch. 2017.
[36] Pedregosa, G. Varoquaux, A. Gramfort, V. Michel, B. Thirion, O. Grisel, M. Blondel, P. Prettenhofer, R. Weiss, V. Dubourg, J. Vanderplas, A. Passos, D. Cournapeau, M. Brucher, M. Perrot, and Duchesnay. Scikit-learn: Machine learning in Python. Journal of Machine Learning Research, 12:2825 2830, 2011.
[37] Gabriel Peyré, Marco Cuturi, et al. Computational optimal transport: With applications to data science. Foundations and Trends@ in Machine Learning, 11(--6)355-600, 2019.
[38] Huy Phan. Pytorch models trained on cifar-10 dataset. PyTorch_CIFARI 2021.
[39] Julien Rabin, Sira Ferradans, and Nicolas Papadakis. Adaptive color transfer with relaxed optimal transport. In 2014 IEEE international conference on image processing (ICIP), pages 4852 4856. IEEE, 2014.
[40] Peter Rousseeuw. Silhouettes: a graphical aid the interpretation and validation of cluster analysis. Journal of computational and applied mathematics, 20:53-65, 1987.
[41] Michael E Sander, Pierre Ablin, Mathieu Blondel, and Gabriel Peyré. Sinkformers: Transform- with doubly stochastic attention. In International Conference on Artificial Intelligence and Statistics, pages 3515-3530. PMLR, 2022.
[42] Richard Sinkhorn. A relationship between arbitrary positive matrices and doubly stochastic matrices. The annals of mathematical statistics, 35(2):876-87, 1964.
[43] Hoa Thi Nhu Tran, Kok Siong Ang, Marion Chevrier, Xiaomeng Zhang, Nicole Yee Shin Lee, Michelle Goh, and Jinmiao Chen. A benchmark of batch-effect correction methods for single-cell ma sequencing data. Genome biology, 21:1-32, 2020.
[44] Hugues Van Assel, Thibault Espinasse, Julien Chiquet, and Franck Picard. A probabilistic graph coupling view of dimension reduction. Neural Information Processing Systems (NeurIPS), 2022.
[45] Laurens Van Der Maaten. Barnes-hut-sne, arXiv preprint arXiv: 2013.
[46] Laurens Van der Maaten and Geoffrey Hinton. Visualizing data using t-sne. Journal of Machine Learning Research (JMLR), 9(11), 2008
12

---

[47] David Van Dijk, Roshan Sharma, Juozas Nainys, Kristina Yim, Pooja Kathail, Ambrose J Carr, Cassandra Burdziak, Kevin R Moon, Christine L Chaffer, Diwakar Pattabiraman, al. Recovering gene interactions from single-cell data using data diffusion. Cell, 174(3):716-729, 2018.
[48] Jarkko Venna and Samuel Kaski. Neighborhood preservation in nonlinear projection methods: An experimental study. In Artificial Neural Networks -ICNN 2001: International Conference Vienna, Austria, August 21-25, 2001 Proceedings 11, pages 485-49) Springer, 2001.
[49] Max Vladymyrov and Miguel Carreira-Perpinan. Entropic affinities: Properties and efficient numerical computation. In International Conference Machine Learning (ICML), pages 477-485. PMLR, 2013.
[50] Ulrike Von Luxburg. A tutorial on spectral clustering, Statistics and computing, 17(4):395 416, 2007.
[51] Ron Zass and Amnon Shashua. A unifying approach to hard and probabilistic clustering. In Tenth IEEE International Conference Computer Vision (ICCV'0S) Volume 1, volume 1, pages 294-301, IEEE, 2005.
[52] Ron Zass and Amnon Shashua. Doubly stochastic normalization for spectral clustering. In B. Schölkopf, J. Platt, and T. Hoffman, editors, Neural Information Processing Systems (NeurIPS). MIT Press, 2006.
[53] Lihi Zenni and Pietro Perona. Self-tuning spectral clustering. Advances in neural information processing systems, 17, 2004.
[54] Shuangfei Zhai, Tatiana Likhomanenko, Etai Littwin, Jason Ramapuram, Dan Busbridge, Yizhe Zhang, Jiatao Gu, and Joshua M. Susskind. $isigmaSreparam: Stable transformer training with spectral reparametrization. 2023.
[55] Dengyong Zhou, Olivier Bousquet, Thomas Lal, Jason Weston, and Bernhard Schölkopf. Learning with local and global consistency. Neural Information Processing Systems (NeurIPS), 16, 2003.
13

---

A Proofs
A.1 Euclidean Projection onto S
For the problem arg minpes || - K[ the Lagrangian takes the form, with W e R
L(P, W) |P P-
KI (W, P -PT).
Cancelling the gradient of & with respect P gives 2(P' K) + W W Thus P
K + (WT- (-)-) -W). Using the symmetry constraint on P yields P 1 (K + +K] Hence have:
arg min | KI
K [K+K])
PES
A.2 From Symmetric Consoay- OT Sinkhorn Iterations
In this section, we derive Sinkhorn iterations from the problem EOT) Let C E D. We start by making the constraints explicit.
min PER
(P,C)
(5)
s.t.
Σ H(P)))
(6)
(E[²]
P1=1, P-P
For the above convex problem the Lagrangian writes, where IR,, f R and r Raxn,
((P,f,u,r== = (P,
(u. (un
_H(P.) + 2(f, 1 - P1) + (T,P-PT).
(8)
iE[m]
Strong duality holds and the first order KKT condition gives for the optimal primal p and dual (v", fr I*) variables:
VpL(P" .f N
(U,"",,,,""" r
v/' log P 1 +I* - I+T = -.
Since p", C E S we have 1* P = f*1 18 Hence C + v² log P - f e f 0. Suppose that n² 0 then the previous reasoning implies that V(i,j), Cij S + fj Using that C E 'D we have C₁ = Cjj thus Vi, ft = 0 and thus this would imply that C = which is not
allowed by hypothesis. Therefore v* * 0 and the entropy constraint saturated at the optimum by complementary slackness. Isolating p then yields:
P = exp ((f @ f C)/v")
(10)
P must be primal feasible in particular P*1 1. This constraint gives us the Sinkhorn fixed point relation for f*:
[r], - LSE ((f C+))v²), where for vector a, use the notation LSE(a) log Σσ exp(ax).
Vi [n].
(11)
A.3 Proof of Proposition I
We recall the result
Proposition 1. Ler C E R without constant rows. Then P solves the entropic affinity problem
EA with cost C if and only if is the unique solution of the convex problem
min
(P.C).
(EA as OT)
PEN
Proof. We begin by rewriting the above problem to make the constraints more explicit.
min
(P,C)
PER Xm
s.t. Vi, H(P,)) 2 log {
11
14

---

By concavity of entropy, one has that the entropy constraint is convex thus the above primal problem convex optimization problem. Moreover, the latter is strictly feasible for any & [n 1]. Therefore Slater's condition is satisfied and strong duality holds.
Introducing the dual variables X E I and R the Lagrangian of the above problem writes:
C(P.A. = (P. C) + (E, (logÉ + 1)1 HP))) P1).
(12)
where we recall that (P) = (H:P::)) Note that we will deal with the constraint P E R directly, hence there is no associated dual variable. Since strong duality holds, for any solution P to the primal problem and any solution (e* X²) to the dual problem, the pair P (e", must satisfy the Karush- Kaai--Kuunn-uhhhee (KKT) conditions. The first-order optimality condition gives:
VpC(P EE" ;², C + diag(=") log P -'1 0.
(first-order)
Assume that there exists é [r:]] such that e O. Then (first-orde)) gives that the 8 row of C is constant which is not allowed by hypothesis. Therefore ex > 0 (i.e., ex has positive entries). Thus isolating P the first order condition results in:
P diag(u) (u) e- diag(E")
(13)
where u exp (A+ 2 E*). This matrix must satisfy the stochasticity constraint P1 1. Hence one has u = 1 Ø (exp (diag(E*) -C)) and p has the form
expl (-Cg//E;)) v(i,)) E [n] P 3 Exex (Ce//E)))
(14)
As consequence E* > 0, complementary slackness the KKT conditions gives us that for i the entropy constraint is saturated i.e., H(P)) ) = log { 1. Therefore P solves the problem (EA).
Conversely any solution of EAA P
E. exp (=C,)/)) with (st) such that H(P;)) = log gives
exp C_
an admissible matrix for minpeN, (P, C) and the associated variables satisfy the KKT conditions
which are sufficient conditions for optimality since the problem is convex.
A.4 Proof of Proposition and Proposition
5
The goal of this section is to prove the following results:
Proposition 4 (Saturation of the entropies). Let C E S with zero diagonal, then SEA with cost C has a unique solution that we denote by pae_ moreover C DD then for at least -1 indices E [n] the solution satisfies H(P:)) logg 1.
Proposition (Solving for SEA). Let C E E D,D,D,(PPP... (P. (A.1- - P1) and g(>,\) mill PERT x(), Strong duality holds for (SEA). Moreover, let Y" X E argmax, 4(7, A) be the optimal dual variables respectively associated with the entropy and marginal constraints. Then, for at least n 1 indices E [m], 78 > O. When Vi E [m]] 2; > 0 then H,(P = (log & 1)1 and pse has the form
(log + 11111H HP)
a exp ((X* - 2 (r +7*)).
(2)
The unicity of the solution Proposition is a consequence the following lemma
Lemma 7. Ler C + e with zero diagonal. Then the problem minpeH,ns (P, C) has a unique solution.
Proof. Making the constraints explicit, the primal problem of symmetric entropic affinity takes the following form
min
(P.C)
PER Km
Vi, H(P,) > log { +
(SEA)
PI=1. P PT
Suppose that the solution is not unique i.e., there exists a couple of optimal solutions (P1, P2) that satisfy the constraints of (SEA) and such that (Pi, C) (P2, C). For i [n]. denote the
15

---

function fi P + (log + 1) - H.).). Then fi is continuous, strictly convex and the entropy
conditions of (SEA) can be written vi [n]] (( S
Now consider Q
1 (Pi + P2). Then clearly Q1 1,Q QT. Since ft is strictly convex
we have [{(Q) ( + P2) < ( S.. Thus [+(Q) < 0 for any i [m].
Take any E > 0 and / E [n]. By continuity of ʃi there exists 5; > 0 such that, for any H with
||||| S 6g we have [i(Q H) <1Q + 8. Take > O such that Vi E [m],0 < - (i())
(this is possible since for any i E [=], < 0) and H with ||||| S min 6,- Then for any
i E [n]. fi(Q + H) < O. In other words, we have proven that there exists 1j > 0 such that for any H
such that [H]] F S 17, it holds: Vi E [m], fQQ H) < O.
Now let us take H as the Laplacian matrix associated to C.e. for any (i.j) [m] Hij - Cyp if i A j and Σ; Ca otherwise. Then we have (H, C) = - -in C2jy +0 = - 0 < 0 since
Eig
C has zero diagonal (and is nonzero). Moreover, H Н since C is symmetric and H1
by
construction. Consider for 0 < B S
the matrix H2 = 3H. Then |H> B H | S 1).
By the previous reasoning one has: Vi E [m], ( + H3) < O. Moreover, (Q + (Q+ = Q + H and (Q+ + Hg)1 = 1. For & small enough we have Q + H E R and thus there is a 3 (that depends on P1 and P2) such that Q + H3 is admissible i.e. satisfies the constraints of (SEA). Then, for such 3,
(C, Q - (1. 1 (C + + (C, (C.PP)
(15)
(C, = BH, <.
Thus (C, Q + H2) < (C,Pi) which leads to a contradiction.
We can now prove the rest of the claims of Proposition and Proposition S
Proof. Let C E D. We first prove Proposition 4 The unicity is a consequence of Lemma 7 For the saturation of the entropies we consider the Lagrangian of the problem (SEA) that writes
= (P, C) + (7, (log E + 1)1 - (P)) + (A, 1 - P1) + (T,P-
for dual variables 7 E TR E I. and r E Rax,, Strong duality holds by Slater's conditions
because "11 T is stricly feasible for & S n 20-11 1. Since strong duality holds, for any solution P to the primal problem and any solution (y*,, >², ,T" to the dual problem, the pair P (r X², I*) must satisfy the KKT conditions. They can be stated as follows:
C+ (-adga((")) log P - X'1 T +T* T =0
P'1=1, H, (log + 1)1, p'=p''
(KKT-SEA)
Y' 20 0
i (H(::) - (log & 1))= 0.
Let us denote = =f E [n| s.t. 7ē =0}. For I, using the first-order condition, one has for i E [m], CR = - - it. Since C E D. we have C2 = 0 thus X; 0 and Cs = Di - T For (L, E 12, one has Cie = - 1 = -(Tr - r = -Crrg C is symmetric thus Cie = O. Since C only has null entries on the diagonal, this shows that E = f" and therefore I has at most one element, By complementary slackness condition (last row of the KKT-SHA conditions) it holds that V' + f, (P ; = log { 1. Since the solution of (SEA) is unique P = s and thus Vi + f, H(P))) = log f + I which proves Proposition 4 but also that for at least n indices 7, > 0.
Moreover, from the KKT conditions we have
v(i.j) E [n]] rj Tig = Cipy + 7/ log Pij X;
(16)
Now take (i,j) E [n] fixed. From the previous equality 1 Li = Cij + of log Pr
x; but also rú - rí = Cs + 7 log Pj - X3. Using that P = (P*)T T and C E S we get Nij rj =
Ci + 2 log P;; - Xj. But r - -T; = -TT) Ti which gives
Og + 73 P.g - -15
>j== = - + 25 log P.j X7)
(17)
16

---

This implies
5( E [m]², 2C,j + (7) + +25 log P
(y) +
(18)
Consequently, y > 0 we have the desired form from the above equation and by complementary slackness (P (log& Ç + 1)1 which proves Proposition 5 Note that otherwise, it holds
+ 2011
v(i,) * (e,e), P.j exp
(19)
7/ + 73
A.5
EA and SEA as KL projection
We prove the characterization as a projection Lemma Let C D and K. exp(-C/o) Then for any min; eq, it holds Pe ProjHs (K,))) arg minpeRs KL(P)K,).
EA in Lemma 8 and of (SEA) in Lemma 9
Proof. The KL projection K onto He reads
min PER
KL(P|K)
(20)
s.t.
Vi, H(P,)) 2 log
(21)
P1
(22)
Introducing the dual variables X E R and & E RR, the Lagrangian of this problem reads:
C(P.A, KK(P (KS (log E 1)1 H())) + (A.1-P1) P1)
(23)
Strong duality holds hence for any solution P to the above primal problem and any solution (k". X*)
the dual problem, the pair P* (x=* X²) must satisfy the KKT conditions. The first-order optimality condition gives:
PpC(P', K", X*) = log (P* @ K) + diag(n' " log P - =0. -.
(24)
Solving for X given the stochasticity constraint and isolating P gives
exp (loog (1,)/(1++5])
Vii,j) E [n] P.
(25)
Σ. exp (llog (1t)
We now consider P* function of s.. Plugging this expression back in yields the dual function + (Kx). The latter is concave as any dual function and its gradient reads:
(xG(K) (log + 1)1 (k))
(26)
Denoting by p = 1 + K. and taking the dual feasibility constraint к > 0 into account gives the
solution: for any i, P max(e;, 1) where e* solves (EA) with cost C = loog K. Moreover
have that a < min(e" where ex E (IR-)" solves EA Therefore for any i E [no]] one has haase//a et/o 2 1. Thus there exists K R, such that a KT) e..
This R* cancels the above gradient i.e., (log + 1)1 = H(P*(K))) thus solves the dual problem. Therefore given the expression of P we have that Proj (K) x pe
KL
Lemma 9. Let C D 0 and K, exp(-C/o). Suppose that the optimal dual variable 7t associated with the entropy constraint of (SEA) is positive. Then for any a S. min, Y. holds p Projkkens (.,))
Proof. Let a > 0. The KL projection of K onto Hg S boils down to the following optimization problem.
min
(L(P|
PER P
s.t.
Vi, H(P.) ..) > log
(SEA-Proj)
P1 P1=1, 1, P P.
17

---

By strong convexity of P - KL(P)Kg) and convexity of the constraints the problem (SEA-Proj) admits unique solution. Moreover, the Lagrangian of this problem takes the following form, where E RR μ E R and r R
14, w. r) KL(P)K.) + (w, (log 1)1 - H,()) (14, - PT) + P
Strong duality holds by Slater's conditions thus the KKT conditions necessary and sufficient. In particular if P and (w-,"² ,u, (²) satisfy
VpC(P" 10, w", I log (P* 2 K) + diag(w") 1²11 8 B3
P-1 =1 (P*) (log + 1)1, P p*
(KKT-Proj)
w² 20 Vi,wi (H(P,) (log { 1)) 0.
then P* solution (SEA-Proj and (w u=, ,8*) are optimal dual variables. The first condition rewrites
1 Cip + wt log(Pi)) - 14 + B3g GJ. O,
v(i,)), log P.t) 0
(27)
which is equivalent to
v(,,), (11 + w/) log( Pú Cú an a ( (1) = 0.
(28)
Now take pse the optimal solution of SEA As written in the proof Proposition 5 of p and the optimal dual variables (7*, ""," satisfy the KKT conditions:
(i,,), Ci + 78 log Pij X + rij r
pº = 1, H, (P*) 2 (log + 1)1, aSr = (P=O)
(KKT-SEA)
?" 2 Vi, 7, (H(P:)) - (log 5 + 1)) = 0.
By hypothesis y* > 0 which gives Vi, H(P - (log & + 1) = 0. Now take 0 < a S min, of
and define Wi zi -. Using the hypothesis on have Vi,wi 0 and w satisfies Vi, a(1 + w = r. Moreover for any E [m]
oa (H(P)) (log E + 1))= =0
(29)
Define also Vi, 14 = / and V(I.). = rit/a. Since pse, A conditions (KTTSEAA) then by the previous reasoning p (w p 3r) (KKT-Proj) and in particular pse is optimal solution of Proj-roo since KKT conditions are sufficient. Thus we have proven that pse E arg minpsegns the solution this is in fact an equality.
X r") satisfies the KKT
satisfy the KKT conditions
KL and by the uniqueness of
B Sensitivity Analysis for Dimensionality Reduction Experiments
In Figure 8
we extend the sensitivity analysis performed for spectral clustering (Figure 5) to DR scores. One can notice that (SNEkhorn outperforms (SNE on wide range of perplexity values.
Silhonetse COIL
Trustworthinees COIL
Silhouctte UMNIST
Trustworthiness UMNIST
0.97
03882
45-
0.08
022
01
0.915
0.35
0.1
LSM
0.063
-02
20
101
131
001
60
L'impasttyy
Perplectyy
Exepirsity
Pieplccity
Figure 8: Dimensionality reduction scores as a function of the perplexity parameter.
18


Reviews:
[{"review_id": "7kDQMfvgwV", "reviewer": "Reviewer_WhX7", "summary": "Existing DR approaches employ EAs after heuristic symmetrisation (symmetric-SNE). This leads to less faithful embeddings with low silhouette scores. This paper avoids such heuristic symmetrisation by enforcing symmetrisation in a related, new OT based formulation.\n\nTowards this goal, firstly, the EA problem is equivalently written as a (semi-relaxed) OT problem, where the transport plan recovers the EAs. Then symmetrisation is explicitly enforced as constraints on the transport map in this OT problem. Using optimality conditions it is shown that such EAs, maintain the crucial row-wise entropic equalities (prop4). As a result EAs are symmetric, doubly stochastic, and satisfy the entropic equalities, making the organically superior to the existing heuristics.\n\nUsing the proposed EAs, a DR formulation is proposed (leading to so-called SNEkhorn). Details of solving this problem are presented. The proposed embeddings are empirically compared to the state-of-the-art wrt. spectral clustering and dimension reduction. The improvements are significant.", "questions": "1. It seems a primary reason for having Q_Z^{ds} in SNEkhorn is to enable fast objective-gradient computation via sinkhorn. Is this true?\n2. Will it be more appropriate to consider other affinities like \\tilde{Q}_Z etc. in place of Q_Z^{ds} ? If \\tilde{Q}_Z leads to spherical embeddings, can there be applications where this is desirable? Any discussion on this and alternatives may help the reader.\n\n", "limitations": "Limitations were discussed in the concluding section.", "rating": 6, "confidence": 4, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "1. The organisation and write-up are well polished and makes it an easy read.\n2. Connections made between EAs and OT are interesting, especially in the light of the symmetrisation issue with existing DR approach.\n3. Empirical improvements over baselines are impressive.\n", "weaknesses": "..", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We would like to thank the reviewer for her/his insightful comments and questions. Our answers to the questions raised can be found below.\n\n> It seems a primary reason for having Q_Z^{ds} in SNEkhorn is to enable fast objective-gradient computation via sinkhorn. Is this true?\n\nIt is true that it enables fast objective-gradient computation via the Sinkhorn algorithm. However the primary reason for using such Q_Z^{ds} is to counter the sphere concentration phenomenon for the embedding occuring when matching a non DS affinity (for instance \\tilde{Q}_Z or Q_Z) with a DS affinity for input data [1]. This geometrical effect disapears when the embeddings affinity is itself DS. This is very useful as soon as one wants to embed onto a flat space.\n\n> Will it be more appropriate to consider other affinities like \\tilde{Q}_Z etc. in place of Q_Z^{ds} ? If \\tilde{Q}_Z leads to spherical embeddings, can there be applications where this is desirable? Any discussion on this and alternatives may help the reader.\n\nThank you for this very interesting comment. Spherical embeddings are desirable if the embedding space is indeed a sphere. This is what is done in [1]. In our work, we rather focused on Euclidean flat spaces as it is what practionners are usually interested in when performing dimensionality reduction (for instance with scRNA seq data). Note that this also suggests as future works to study SNEkhorn for embedding on spherical manifolds. This discussion will be added in the paper.\n\n[1] Yao Lu, Jukka Corander, and Zhirong Yang. Doubly stochastic neighbor embedding on spheres. Pattern Recognition Letters, 128:100–106, 2019."}, {"review_id": "FfgSzpgIx1", "reviewer": "Reviewer_s51e", "summary": "The paper presents a novel approach to dealing with entropic affinities (EAs) used in machine learning for dimensionality reduction tasks, specifically in the popular t-SNE algorithm. It addresses the limitations of current symmetrization methods applied to EAs, which can compromise the entropy and stochasticity properties of the affinity matrix. The proposed method uses optimal transport to achieve a natural symmetrization, leading to a new affinity matrix. This new matrix is then leveraged in a new dimensionality reduction algorithm called SNEkhorn, which is demonstrated to outperform state-of-the-art methods on both synthetic and real-world datasets.", "questions": "1. What sensitivity analysis was performed on the hyperparameters introduced by the new symmetrization approach? How robust is the SNEkhorn algorithm to changes in these parameters?\n\n2. Can you provide additional insights into the interpretability of the reduced representations obtained by SNEkhorn? How well does it preserve meaningful structure in the data compared to other methods? Is the visualization in Figure 6 a special case?\n\n3.Could you elaborate on the limitations or potential failure cases of the proposed approach? What are the scenarios where SNEkhorn might not perform as well or could face challenges?", "limitations": "While the paper proposes a novel approach to symmetrizing entropic affinities and introduces the SNEkhorn algorithm with promising results, it is essential to consider the potential limitations and shortcomings of the work. Some possible limitations include:\n\nScalability: Dimensionality reduction algorithms often face challenges with scalability when dealing with large datasets. The paper should address the computational efficiency and scalability of the proposed SNEkhorn algorithm, especially when applied to high-dimensional and massive datasets commonly encountered in real-world applications.\n\nParameter Sensitivity: It is essential to investigate the sensitivity of hyper parameters and assess how they impact the performance of the SNEkhorn algorithm, e.g. perplexity.\n", "rating": 7, "confidence": 3, "soundness": 4, "presentation": 3, "contribution": 3, "strengths": "1. The paper introduces a novel approach to symmetrizing entropic affinities by formulating them as an optimal transport problem. This brings a fresh perspective to the problem and could open up new avenues for research and application in machine learning.\n2. The proposed symmetrization method claims to maintain the constant entropy and stochasticity properties of the affinity matrix while being computationally efficient through dual ascent. This combination of robustness and efficiency is essential for practical applications.\n3. The development of the SNEkhorn algorithm based on the new affinity matrix offers a practical implementation of the proposed approach. The demonstrated superiority of SNEkhorn over other baseline methods on various datasets further strengthens the paper's claims.\n4. The paper evaluates the proposed SNEkhorn algorithm on both synthetic and real-world datasets, providing a comprehensive assessment of its performance and general applicability.", "weaknesses": "Clarity of Presentation: While the abstract provides a high-level overview, some concepts, such as entropic affinities, optimal transport, and dual ascent, might be challenging for readers not already familiar with the domain. The paper should provide a clear and concise introduction to these concepts, ideally with intuitive explanations and illustrative examples. For example, what is P and C on paper 3 \"Symmetric Entropy-Constrained Optimal Transport\"\n\nRobustness Analysis: Since the proposed method claims to be robust to varying noise levels, it would be valuable to include a thorough analysis of its performance under different levels of noise in the datasets instead of two level. This could strengthen the claim of robustness and highlight the algorithm's practicality in real-world scenarios.\n", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We would like to thank the reviewer for the appreciation of our work, insightful comments and questions. \n\n### Answers to weaknesses:\n\n> Clarity of Presentation: While the abstract provides a high-level overview, some concepts, such as entropic affinities, optimal transport, and dual ascent, might be challenging for readers not already familiar with the domain. The paper should provide a clear and concise introduction to these concepts, ideally with intuitive explanations and illustrative examples. For example, what is P and C on paper 3 \"Symmetric Entropy-Constrained Optimal Transport\"\n\nWe will add details about these concepts in the paper as it may indeed appear challenging for some readers. Note that in our article C is the transportation cost matrix (C_{ij} is the cost of transporting a unit of mass from sample i to sample j) while P is the transport plan (P_{ij} is the amount of mass transported from sample i to sample j) that we propose to use as affinity matrix. We will make it clearer in the revised version.\n\n> Robustness Analysis: Since the proposed method claims to be robust to varying noise levels, it would be valuable to include a thorough analysis of its performance under different levels of noise in the datasets instead of two level. This could strengthen the claim of robustness and highlight the algorithm's practicality in real-world scenarios.\n\nThis is a very interesting question. We provide some new results studying the robustness to noise in the global reply (Figure 1 in the pdf).\n\n### Answers to questions\n\n> What sensitivity analysis was performed on the hyperparameters introduced by the new symmetrization approach? How robust is the SNEkhorn algorithm to changes in these parameters?\n\nWe provide a sensitivity analysis in Fig. 5 in the paper for spectral clustering applications. We also added in the global reply (see Figure 2 in the adjoining PDF) a new sensitivity analysis for the end task of dimensionality reduction. This shows the superiority of t-SNEkhorn over t-SNE over a wide range of perplexity values.\n\n> Can you provide additional insights into the interpretability of the reduced representations obtained by SNEkhorn? How well does it preserve meaningful structure in the data compared to other methods? Is the visualization in Figure 6 a special case?\n\nSimilarly to SNE (or UMAP), SNEkhorn's axis cannot be directly interpreted due to the non-linearity of the model.\nAs shown in section 5, using symmetric entropic affinities and SNEkhorn lead to better robustness to heteroscedastic noise as well as enhanced clustering abilities compared to previous approaches. As such, the partitioning structure as well as close range pairwise relations are better preserved in SNEkhorn. \nWe displayed figure 1 and 6 to give a visual interpretation of the figures given in table 3. It is indeed an example where SNEkhorn performs much better than other methods.\n\n> 3.Could you elaborate on the limitations or potential failure cases of the proposed approach? What are the scenarios where SNEkhorn might not perform as well or could face challenges?\n\nFirst note that while our proposed symmetric EAs are much more robust than the l2 symmetrization carried out in t-SNE, we will share many failure cases with tSNE, i.e. when the perplexity is too large (loss of local subtleties) or too small (fake structures appearing). In addition, our method comes with a more involved computational cost since we need to optimize twice the number of parameters (vectors $\\boldsymbol\\gamma$ and $\\boldsymbol\\lambda$) for affinity matrix estimation (instead of only $\\boldsymbol\\epsilon$ for the directed EAs). For runtime values, we refer to table 1 of the pdf in the global answer above."}, {"review_id": "Fk0Ue2Q73U", "reviewer": "Reviewer_W2ip", "summary": "This is a very interesting paper about an application of the Sinkhorn algorithm to symmetrize the matrix of entropic affinities in methods of DR like SNE, t-SNE, etc.\nStrong theoretical contribution.\nSome experiments to illustrate.", "questions": "The \"crowding effect\" is just an intuitive explanation found in LvdM t-SNE paper; would it not be better to refer to something more conventional (empty space phenomenon, norm concentration, or just the reduction in volume from HD to LD?)?\nIn Figure 1, SNEkhorn or t-SNEkhorn ?", "limitations": "No possible negative societal impact here.", "rating": 8, "confidence": 5, "soundness": 4, "presentation": 3, "contribution": 4, "strengths": "The paper is very interesting for its vision, state of the art across multiple domains, strong theoretical contribution, revisiting of existing methods and proposal of a new method.\nThe paper is very well written, not rushed.\nThe symmetrization with Sinkhorn is well done and this a nice \"tour de force\" like only Brittons can do.\nThe experiments include a grid search on the perplexity, and multiple runs for each perplexity.", "weaknesses": "The paper is well written but it is also very/too dense: much information is packed, that is nice, but it takes quite some time to digest and to decode.\nThe notation is sometimes unconventional and difficult to follow, probably to keep it as compact as possible.\nWhile the theoretical part is very convincing and interesting for itself, it will only appeal to a very restricted audience of people specialised in the methodological design of advanced methods of DR. Also, the experiments are not sufficient to convince practitioners that the Sinkhorn symmetrization is from now on something necessary to get the best results. The experiments lack broader comparisons with other methods that are in the state of the art (here it looks like the authors have reimplemented t-SNE and UMAP; also, there are methods that do not look for symmetry at all, although this can be detrimental in the presence of many isolated outliers). The quantitative assessment might be partly questionable: the trustworthiness is not the best DR QA indicator; what was the neighbourhood size, by the way?\nThe lack of accelerated version is a possible weakness for dissemination, while NeurIPS would be precisely the ideal launching pad.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We would like to thank the reviewer for his/her careful reading of the manuscript and insightful comments.\n\n> The quantitative assessment might be partly questionable: the trustworthiness is not the best DR QA indicator; what was the neighbourhood size, by the way?\n\nThis is an interesting point. We wanted to fairly evaluate performances of the DR methods methods and to avoid ad-hoc performance measures that are used in the optimization problem of specific methods. This is why we chose classical (if a bit unorthodox) performance measures that evaluate the preservation of global (silhouette score) and local (trustworthiness) structures from the data. We are interested if the reviewer suggests an alternative that we might add in supplementary. For trustworthiness, the neighbourhood size was set to the default value of 5.\n\n> The \"crowding effect\" is just an intuitive explanation found in LvdM t-SNE paper; would it not be better to refer to something more conventional (empty space phenomenon, norm concentration, or just the reduction in volume from HD to LD?)? In Figure 1, SNEkhorn or t-SNEkhorn ?\n\nThe reviewer is raising a good point. We mentionned the \"crowding effect\" to be faithful to the story of t-SNE and explain how heavy-tailed kernel first appeared in these methods. We agree indeed that a more conventional explanation will be better and will modify our article accordingly for the sake of clarity. \n\nFor figure 1, we indeed used SNEkhorn in the embedding space as we found that it gave descent results on the COIL dataset and we lacked space to provide both SNEkhorn and t-SNEkhorn embeddings. We will plot both methods in the supplementary in the final version of the paper.\n\nNote that quantitative results for t-SNE and t-SNEkhorn on COIL are provided in table 3."}, {"review_id": "nKT5g39VeS", "reviewer": "Reviewer_8Btp", "summary": "This paper presents a new dimensionality reduction algorithm, named SNEkhorn. By uncovering the novel links between Entropic affinities (EAs) and Optimal Transport (OT), the authors derive EAs with symmetric doubly stochastic normalization and the fixed row-wise entropy, which is the key to SNEkhorn. Besides, the authors show a dual ascent algorithm to compute this new affinity efficiently. Beyond this theoretical contribution, the authors showcase the benefits of SNEkhorn through numerical experiments on simulated data and real data.", "questions": "- How does the computational complexity for the proposed method compare to t-SNE? Does the proposed method spend the same amount of time as t-SNE while achieve better results?\n- How to choose the dimensionality of PCA? In the reference the authors mentioned, the dimensionality of PCA is usually set to be 30 instead of 50. Does the dimensionality affects the performance of the algorithm?", "limitations": "Yes", "rating": 6, "confidence": 3, "soundness": 4, "presentation": 3, "contribution": 3, "strengths": "- Entropic affinities (EAs) needs to be symmetrized when used in popular Dimensionality Reduction (DR) algorithm t-SNE, but the symmetrization destroys the row-wise constant entropy and stochasticity properties of EAs. To derive EAs that can maintain both symmetry and its own properties is a very natural. The proposed method is novel and simple to describe.\n- The authors discover novel connections between EAs and OT, and show a dual ascent algorithm to compute EAs from OT. The proposed algorithm provides a new perspective on EAs.\n- For conventional DR algorithm, the latent affinity is not doubly stochastic, which imposes spherical constraints on the embedding space. The proposed methods avoid this problem.\n- The simulation study and real data analysis are very sufficient. The authors illustrate the effectiveness from many aspects.", "weaknesses": "- There is not any discussion about the computational complexity of proposed method. Although the proposed algorithm performs better than t-SNE and UMAP, how does the computational complexity of proposed method compare with t-SNE or UMAP. Since the computational cost is also important in practice. \n- Some of the formulas and references in the article are incorrectly hyperlinked.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We thank the reviewer for the careful reading of the manuscript, her/his assessment and relevant remarks.\n\n> How does the computational complexity for the proposed method compare to t-SNE? Does the proposed method spend the same amount of time as t-SNE while achieve better results?\n\nThis is a very good point. We provide some additional insights about computational aspects in the general answers as well as some runtime experiments in the pdf. Although both methods have a complexity that scales quadratically with the number of samples, one can notice in these experiments that the current implementation of t-SNEkhorn is slower than t-SNE. We will add a discussion and computational times in the final version of the paper.\n\nImportantly, note that in our paper we did not use the most recent acceleration procedures that have been proposed for t-SNE [2,3] but our approach will directly benefit from many of those computational tricks.\n\n> How to choose the dimensionality of PCA? In the reference the authors mentioned, the dimensionality of PCA is usually set to be 30 instead of 50. Does the dimensionality affects the performance of the algorithm?\n\nThis point indeed deserves interest as PCA is often applied as a pre-processing step before computing the distance matrix. Choosing the dimensionality of this PCA step seems to depend on the nature of the data (see e.g. [4] for scRNA data that recommends using 50 axes). In his t-SNE website (unlike in the paper as the reviewer rightfully pointed out), Laurens Van der Maaten also recommends using 50 axes. In our experiments, as shown for instance in table 2 of the pdf above provided with the global response, we did not find any significant difference between using 30 axes and 50 axes (although this might not be true for all datasets).\n \n[2] Laurens Van Der Maaten. Barnes-hut-sne. arXiv preprint arXiv:1301.3342, 2013.\n\n[3] Linderman, G. C., Rachh, M., Hoskins, J. G., Steinerberger, S. and Kluger, Y. (2019). Fast interpolation-based t-SNE for improved vi- sualization of single-cell RNA-seq data. Nature Methods 16 243-245.\n\n[4] Dmitry Kobak and Philipp Berens. The art of using t-SNE for single-cell transcriptomics.\nNature Communication, 10:5416, 2019."}]

Author rebuttal:
We first would like to thank all the reviewers for their remarks and questions.

You may find attached a pdf with some new results to answer the various points raised. The new results are as follows.

## New results

### Figure 1 : Robustness to noise

In Figure 1, we focus on reviewer s51e's question about the robustness analysis. We plot spectral clustering scores on the clustering of three Gaussian clusters with variance $\sigma$, $2\sigma$ and $3\sigma$. This figure shows that symmetric entropic affinities are more robust than the Sinkhorn kernel when $\sigma$ increases. Note that the later is known for its robustness to heteroscedastic noise as shown in [5].

### Figure 2 : Sensitivity analysis

In figure 2, we extend the sensitivity analysis perfomed for spectral clustering (Figure 5 in the submitted paper) to DR scores. This figure shows that tSNEkhorn outperforms tSNE on a wide range of perplexity values. We hope this partly answers reviewer s51e's request about sensitivity analysis.

### Table 1 : Runtime experiments

In table 1, as requested by reviewer 8Btp, we display the runtimes of tSNE and tSNEkhorn on a few datasets. Even though both algorithms have a quadratic complexity with respect to the number of samples, tSNEkhorn's optimization is a bit more involved and requires more time than tSNE.

### Table 2 : Sensibility to PCA dimensions

Finally, to answer the question of reviewer 8Btp, we test different values of dimensions for the pre-processing PCA step. We find that the performance of both tSNE and tSNEkhorn, measured by the silhouette score on COIL/Olivetti/UMNIST, is stable accross these values.


[5] Boris Landa, Ronald R Coifman, and Yuval Kluger. Doubly stochastic normalization of the
gaussian kernel is robust to heteroskedastic noise. SIAM journal on mathematics of data science, 3(1):388–413, 2021.

Meta-review: The paper brings new insights by developing a connection between Entropic Affinities and Optimal Transport. This new connections yields a novel Dimensionality Reduction algorithm whose name is the title of the paper. There was consensus that this is a paper which will be of interest to Neuritis community.
Decision: Accept (poster)