A TRAINABLE OPTIMAL TRANSPORT EMBEDDING
Year: 2021

Abstract:
We address the problem of learning on sets of features, motivated by the need of performing pooling operations in long biological sequences of varying sizes, with long-range dependencies, and possibly few labeled data. To address this challenging task, we introduce a parametrized representation of fixed size, which embeds and then aggregates elements from a given input set according to the optimal transport plan between the set and a trainable reference. Our approach scales to large datasets and allows end-to-end training of the reference, while also providing a simple unsupervised learning mechanism with small computational cost. Our aggregation technique admits two useful interpretations: it may be seen as a mechanism related to attention layers in neural networks, or it may be seen as a scalable surrogate of a classical optimal transport-based kernel. We experimentally demonstrate the effectiveness of our approach on biological sequences, achieving state-of-the-art results for protein fold recognition and detection of chromatin profiles tasks, and, as a proof of concept, we show promising results for processing natural language sequences. We provide an open-source implementation of our embedding that can be used alone or as a module in larger learning models at https://github.com/claying/OTK.

heading: 1 INTRODUCTION
text: Many scientific fields such as bioinformatics or natural language processing (NLP) require processing sets of features with positional information (biological sequences, or sentences represented by a set of local features). These objects are delicate to manipulate due to varying lengths and potentially long-range dependencies between their elements. For many tasks, the difficulty is even greater since the sets can be arbitrarily large, or only provided with few labels, or both.
Deep learning architectures specifically designed for sets have recently been proposed (Lee et al., 2019; Skianis et al., 2020). Our experiments show that these architectures perform well for NLP tasks, but achieve mixed performance for long biological sequences of varying size with few labeled data. Some of these models use attention (Bahdanau et al., 2015), a classical mechanism for aggregating features. Its typical implementation is the transformer (Vaswani et al., 2017), which has shown to achieve state-of-the-art results for many sequence modeling tasks, e.g, in NLP (Devlin et al., 2019) or in bioinformatics (Rives et al., 2019), when trained with self supervision on large-scale data. Beyond sequence modeling, we are interested in this paper in finding a good representation for sets of features of potentially diverse sizes, with or without positional information, when the amount of training data may be scarce. To this end, we introduce a trainable embedding, which can operate directly on the feature set or be combined with existing deep approaches.
More precisely, our embedding marries ideas from optimal transport (OT) theory (Peyré & Cuturi, 2019) and kernel methods (Schölkopf & Smola, 2001). We call this embedding OTKE (Optimal
∗Equal contribution. †Univ. Grenoble Alpes, Inria, CNRS, Grenoble INP, LJK, 38000 Grenoble, France. ‡D.I., UMR 8548, École normale supérieure, Paris, France.
Transport Kernel Embedding). Concretely, we embed feature vectors of a given set to a reproducing kernel Hilbert space (RKHS) and then perform a weighted pooling operation, with weights given by the transport plan between the set and a trainable reference. To gain scalability, we then obtain a finite-dimensional embedding by using kernel approximation techniques (Williams & Seeger, 2001). The motivation for using kernels is to provide a non-linear transformation of the input features before pooling, whereas optimal transport allows to align the features on a trainable reference with fast algorithms (Cuturi, 2013). Such combination provides us with a theoretically grounded, fixed-size embedding that can be learned either without any label, or with supervision. Our embedding can indeed become adaptive to the problem at hand, by optimizing the reference with respect to a given task. It can operate on large sets with varying size, model long-range dependencies when positional information is present, and scales gracefully to large datasets. We demonstrate its effectiveness on biological sequence classification tasks, including protein fold recognition and detection of chromatin profiles where we achieve state-of-the-art results. We also show promising results in natural language processing tasks, where our method outperforms strong baselines. Contributions. In summary, our contribution is three-fold. We propose a new method to embed sets of features of varying sizes to fixed size representations that are well adapted to downstream machine learning tasks, and whose parameters can be learned in either unsupervised or supervised fashion. We demonstrate the scalability and effectiveness of our approach on biological and natural language sequences. We provide an open-source implementation of our embedding that can be used alone or as a module in larger learning models.

heading: 2 RELATED WORK
text: Kernel methods for sets and OT-based kernels. The kernel associated with our embedding belongs to the family of match kernels (Lyu, 2004; Tolias et al., 2013), which compare all pairs of features between two sets via a similarity function. Another line of research builds kernels by matching features through the Wasserstein distance. A few of them are shown to be positive definite (Gardner et al., 2018) and/or fast to compute (Rabin et al., 2011; Kolouri et al., 2016). Except for few hyperparameters, these kernels yet cannot be trained end-to-end, as opposed to our embedding that relies on a trainable reference. Efficient and trainable kernel embeddings for biological sequences have also been proposed by Chen et al. (2019a;b). Our work can be seen as an extension of these earlier approaches by using optimal transport rather than mean pooling for aggregating local features, which performs significantly better for long sequences in practice. Deep learning for sets. Deep Sets (Zaheer et al., 2017) feed each element of an input set into a feed-forward neural network. The outputs are aggregated following a simple pooling operation before further processing. Lee et al. (2019) propose a Transformer inspired encoder-decoder architecture for sets which also uses latent variables. Skianis et al. (2020) compute some comparison costs between an input set and reference sets. These costs are then used as features in a subsequent neural network. The reference sets are learned end-to-end. Unlike our approach, such models do not allow unsupervised learning. We will use the last two approaches as baselines in our experiments. Interpretations of attention. Using the transport plan as an ad-hoc attention score was proposed by Chen et al. (2019c) in the context of network embedding to align data modalities. Our paper goes beyond and uses the transport plan as a principle for pooling a set in a model, with trainable parameters. Tsai et al. (2019) provide a view of Transformer’s attention via kernel methods, yet in a very different fashion where attention is cast as kernel smoothing and not as a kernel embedding.

heading: 3 PROPOSED EMBEDDING

heading: 3.1 PRELIMINARIES
text: We handle sets of features in Rd and consider sets x living in X = { x | x = {x1, . . . ,xn} such that x1, . . . ,xn ∈ Rd for some n ≥ 1 } .
Elements of X are typically vector representations of local data structures, such as k-mers for sequences, patches for natural images, or words for sentences. The size of x denoted by n may vary, which is not an issue since the methods we introduce may take a sequence of any size as input, while providing a fixed-size embedding. We now revisit important results on optimal transport and kernel methods, which will be useful to describe our embedding and its computation algorithms.
Optimal transport. Our pooling mechanism will be based on the transport plan between x and x′ seen as weighted point clouds or discrete measures, which is a by-product of the optimal transport problem (Villani, 2008; Peyré & Cuturi, 2019). OT has indeed been widely used in alignment problems (Grave et al., 2019). Throughout the paper, we will refer to the Kantorovich relaxation of OT with entropic regularization, detailed for example in (Peyré & Cuturi, 2019). Let a in ∆n (probability simplex) and b in ∆n ′ be the weights of the discrete measures ∑ i aiδxi and ∑ j bjδx′j with respective locations x and x′, where δx is the Dirac at position x. Let C in Rn×n ′
be a matrix representing the pairwise costs for aligning the elements of x and x′. The entropic regularized Kantorovich relaxation of OT from x to x′ is
min P∈U(a,b) ∑ ij CijPij − εH(P), (1)
where H(P) = − ∑ ij Pij(log(Pij) − 1) is the entropic regularization with parameter ε, which controls the sparsity of P, and U is the space of admissible couplings between a and b:
U(a,b) = {P ∈ Rn×n ′ + : P1n = a and P >1n′ = b}.
The problem is typically solved by using a matrix scaling procedure known as Sinkhorn’s algorithm (Sinkhorn & Knopp, 1967; Cuturi, 2013). In practice, a and b are uniform measures since we consider the mass to be evenly distributed between the points. P is called the transport plan, which carries the information on how to distribute the mass of x in x′ with minimal cost. Our method uses optimal transport to align features of a given set to a learned reference.
Kernel methods. Kernel methods (Schölkopf & Smola, 2001) map data living in a space X to a reproducing kernel Hilbert space H, associated to a positive definite kernel K through a mapping function ϕ : X → H, such that K(x,x′) = 〈ϕ(x), ϕ(x′)〉H. Even though ϕ(x) may be infinitedimensional, classical kernel approximation techniques (Williams & Seeger, 2001) provide finitedimensional embeddings ψ(x) in Rk such that K(x,x′) ≈ 〈ψ(x), ψ(x′)〉. Our embedding for sets relies in part on kernel method principles and on such a finite-dimensional approximation.

heading: 3.2 OPTIMAL TRANSPORT EMBEDDING AND ASSOCIATED KERNEL
text: We now present the OTKE, an embedding and pooling layer which aggregates a variable-size set or sequence of features into a fixed-size embedding. We start with an infinite-dimensional variant living in a RKHS, before introducing the finite-dimensional embedding that we use in practice.
Infinite-dimensional embedding in RKHS. Given a set x and a (learned) reference z in X with p elements, we consider an embedding Φz(x) which performs the following operations: (i) initial embedding of the elements of x and z to a RKHS H; (ii) alignment of the elements of x to the elements of z via optimal transport; (iii) weighted linear pooling of the elements x into p bins, producing an embedding Φz(x) inHp, which is illustrated in Figure 1. Before introducing more formal details, we note that our embedding relies on two main ideas:
• Global similarity-based pooling using references. Learning on large sets with long-range interactions may benefit from pooling to reduce the number of feature vectors. Our pooling rule follows an inductive bias akin to that of self-attention: elements that are relevant to each other for the task at hand should be pooled together. To this end, each element in the reference set corresponds to a pooling cell, where the elements of the input set are aggregated through a weighted sum. The weights simply reflect the similarity between the vectors of the input set and the current vector in the reference. Importantly, using a reference set enables to reduce the size of the “attention matrix” from quadratic to linear in the length of the input sequence. • Computing similarity weights via optimal transport. A computationally efficient similarity score between two elements is their dot-product (Vaswani et al., 2017). In this paper, we rather consider that elements of the input set should be pooled together if they align well with the same part of the reference. Alignment scores can efficiently be obtained by computing the transport plan between the input and the reference sets: Sinkhorn’s algorithm indeed enjoys fast solvers (Cuturi, 2013).
We are now in shape to give a formal definition. Definition 3.1 (The optimal transport kernel embedding). Let x = (x1, . . . ,xn) in X be an input set of feature vectors and z = (z1, . . . , zp) in X be a reference set with p elements. Let κ be a
positive definite kernel, e.g., Gaussian kernel, with RKHSH and ϕ : Rd → H, its associated kernel embedding. Let κ be the n× p matrix which carries the comparisons κ(xi, zj), before alignment. Then, the transport plan between x and z, denoted by the n × p matrix P(x, z), is defined as the unique solution of (1) when choosing the cost C = −κ, and our embedding is defined as
Φz(x) := √ p× ( n∑ i=1 P(x, z)i1ϕ(xi), . . . , n∑ i=1 P(x, z)ipϕ(xi) ) = √ p×P(x, z)>ϕ(x),
where ϕ(x) := [ϕ(x1), . . . , ϕ(xn)]>.
Interestingly, it is easy to show that our embedding Φz(x) is associated to the positive definite kernel
Kz(x,x ′) := n∑ i,i′=1 Pz(x,x ′)ii′κ(xi,x ′ i′) = 〈Φz(x),Φz(x′)〉, (2)
with Pz(x,x′) := p × P(x, z)P(x′, z)>. This is a weighted match kernel, with weights given by optimal transport in H. The notion of pooling in the RKHS H of κ arises naturally if p ≤ n. The elements of x are non-linearly embedded and then aggregated in “buckets”, one for each element in the reference z, given the values of P(x, z). This process is illustrated in Figure 1. We acknowledge here the concurrent work by Kolouri et al. (2021), where a similar embedding is used for graph representation. We now expose the benefits of this kernel formulation, and its relation to classical non-positive definite kernel. Kernel interpretation. Thanks to the gluing lemma (see, e.g., Peyré & Cuturi, 2019), Pz(x,x′) is a valid transport plan and, empirically, a rough approximation of P(x,x′). Kz can therefore be seen as a surrogate of a well-known kernel (Rubner et al., 2000), defined as
KOT(x,x ′) := n∑ i,i′=1 P(x,x′)ii′κ(xi,x ′ i′). (3)
When the entropic regularization ε is equal to 0, KOT is equivalent to the 2-Wasserstein distance W2(x,x
′) with ground metric dκ induced by kernel κ. KOT is generally not positive definite (see Peyré & Cuturi (2019), Chapter 8.3) and computationally costly (the number of transport plans to compute is quadratic in the number of sets to process whereas it is linear for Kz). Now, we show the relationship between this kernel and our kernel Kz, which is proved in Appendix B.1. Lemma 3.1 (Relation between P(x,x′) and Pz(x,x′) when ε = 0). For any x, x′ and z in X with lengths n, n′ and p, by denoting W z2 (x,x ′) := 〈Pz(x,x′), d2κ(x,x′)〉 1/2 we have
|W2(x,x′)−W z2 (x,x′)| ≤ 2 min(W2(x, z),W2(x′, z)). (4)
This lemma shows that the distanceW z2 resulting fromKz is related to the Wasserstein distanceW2; yet, this relation should not be interpreted as an approximation error as our goal is not to approximate W2, but rather to derive a trainable embedding Φz(x) with good computational properties. Lemma 3.1 roots our features and to some extent self-attention in a rich optimal transport literature. In fact, W z2 is equivalent to a distance introduced by Wang et al. (2013), whose properties are further studied by Moosmüller & Cloninger (2020). A major difference is that W z2 crucially relies on Sinkhorn’s algorithm so that the references can be learned end-to-end, as explained below.

heading: 3.3 FROM INFINITE-DIMENSIONAL KERNEL EMBEDDING TO FINITE DIMENSION
text: In some cases, ϕ(x) is already finite-dimensional, which allows to compute the embedding Φz(x) explicitly. This is particularly useful when dealing with large-scale data, as it enables us to use our method for supervised learning tasks without computing the Gram matrix, which grows quadratically in size with the number of samples. When ϕ is infinite or high-dimensional, it is nevertheless possible to use an approximation based on the Nyström method (Williams & Seeger, 2001), which provides an embedding ψ : Rd → Rk such that
〈ψ(xi), ψ(x′j)〉Rk ≈ κ(xi,x′j).
Concretely, the Nyström method consists in projecting points from the RKHS H onto a linear subspace F , which is parametrized by k anchor points F = Span(ϕ(w1), . . . , ϕ(wk)). The corresponding embedding admits an explicit form ψ(xi) = κ(w,w)−1/2κ(w,xi), where κ(w,w) is the k × k Gram matrix of κ computed on the set w = {w1, . . . ,wk} of anchor points and κ(w,xi) is in Rk. Then, there are several ways to learn the anchor points: (a) they can be chosen as random points from data; (b) they can be defined as centroids obtained by K-means, see Zhang et al. (2008); (c) they can be learned by back-propagation for a supervised task, see Mairal (2016).
Using such an approximation within our framework can be simply achieved by (i) replacing κ by a linear kernel and (ii) replacing each element xi by its embedding ψ(xi) in Rk and considering a reference set with elements in Rk. By abuse of notation, we still use z for the new parametrization. The embedding, which we use in practice in all our experiments, becomes simply
Φz(x) = √ p× ( n∑ i=1 P(ψ(x), z)i1ψ(xi), . . . , n∑ i=1 P(ψ(x), z)ipψ(xi) ) = √ p×P(ψ(x), z)>ψ(x) ∈ Rp×k, (5)
where p is the number of elements in z. Next, we discuss how to learn the reference set z.
3.4 UNSUPERVISED AND SUPERVISED LEARNING OF PARAMETER z
Unsupervised learning. In the fashion of the Nyström approximation, the p elements of z can be defined as the centroids obtained by K-means applied to all features from training sets in X . A corollary of Lemma 3.1 suggests another algorithm: a bound on the deviation term between W2 and W z2 for m samples (x 1, . . . ,xm) is indeed
E2 := 1 m2 m∑ i,j=1 |W2(xi,xj)−W z2 (xi,xj)|2 ≤ 4 m m∑ i=1 W 22 (x i, z). (6)
The right-hand term corresponds to the objective of the Wasserstein barycenter problem (Cuturi & Doucet, 2013), which yields the mean of a set of empirical measures (here the x’s) under the OT metric. The Wasserstein barycenter is therefore an attractive candidate for choosing z. K-means can be seen as a particular case of Wasserstein barycenter when m = 1 (Cuturi & Doucet, 2013; Ho et al., 2017) and is faster to compute. In practice, both methods yield similar results, see Appendix C, and we thus chose K-means to learn z in unsupervised settings throughout the experiments. The anchor points w and the references z may be then computed using similar algorithms; however, their mathematical interpretation differs as exposed above. The task of representing features (learning w in Rd for a specific κ) is decoupled from the task of aggregating (learning the reference z in Rk).
Supervised learning. As mentioned in Section 3.1, P(ψ(x), z) is computed using Sinkhorn’s algorithm, recalled in Appendix A, which can be easily adapted to batches of samples x, with possibly varying lengths, leading to GPU-friendly forward computations of the embedding Φz. More important, all Sinkhorn’s operations are differentiable, which enables z to be optimized with stochastic gradient descent through back-propagation (Genevay et al., 2018), e.g., for minimizing a classification or regression loss function when labels are available. In practice, a small number of Sinkhorn iterations (e.g., 10) is sufficient to compute P(ψ(x), z). Since the anchors w in the embedding layer below can also be learned end-to-end (Mairal, 2016), our embedding can thus be used as a module injected into any model, e.g, a deep network, as demonstrated in our experiments.

heading: 3.5 EXTENSIONS
text: Integrating positional information into the embedding. The discussed embedding and kernel do not take the position of the features into account, which may be problematic when dealing with structured data such as images or sentences. To this end, we borrow the idea of convolutional kernel networks, or CKN (Mairal, 2016; Mairal et al., 2014), and penalize the similarities exponentially with the positional distance between a pair of elements in the input and reference sequences. More precisely, we multiply P(ψ(x), z) element-wise by a distance matrix S defined as:
Sij = e − 1 σ2pos (i/n−j/p)2 ,
and replace it in the embedding. With similarity weights based both on content and position, the kernel associated to our embedding can be viewed as a generalization of the CKNs (whose similarity weights are based on position only), with feature alignment based on optimal transport. When dealing with multi-dimensional objects such as images, we just replace the index scalar i with an index vector of the same spatial dimension as the object, representing the positions of each dimension.
Using multiple references. A naive reconstruction using different references z1, . . . , zq in X may yield a better approximation of the transport plan. In this case, the embedding of x becomes
Φz1,...,zq (x) = 1/ √ q (Φz1(x), . . . ,Φzq (x)) , (7)
with q the number of references (the factor 1/√q comes from the mean). Using Eq. (4), we can obtain a bound similar to (6) for a data set of m samples (x1, . . . ,xm) and q references (see Appendix B.2 for details). To choose multiple references, we tried a K-means algorithm with 2-Wasserstein distance for assigning clusters, and we updated the centroids as in the single-reference case. Using multiple references appears to be useful when optimizing z with supervision (see Appendix C).

heading: 4 RELATION BETWEEN OUR EMBEDDING AND SELF-ATTENTION
text: Our embedding and a single layer of transformer encoder, recalled in Appendix A, share the same type of inductive bias, i.e, aggregating features relying on similarity weights. We now clarify their relationship. Our embedding is arguably simpler (see respectively size of attention and number of parameters in Table 1), and may compete in some settings with the transformer self-attention as illustrated in Section 5.
Shared reference versus self-attention. There is a correspondence between the values, attention matrix in the transformer and ϕ, P in Definition 3.1, yet also noticeable differences. On the one hand, Φz aligns a given sequence x with respect to a reference z, learned with or without supervision, and shared across the data set. Our weights are computed using optimal transport. On the other hand, a transformer encoder performs self-alignment: for a given xi, features are aggregated depending on a similarity score between xi and the elements of x only. The similarity score is a matrix product between queries Q and keys K matrices, learned with supervision and shared across the data set. In this regard, our work complements a recent line of research questioning the dotproduct, learned self-attention (Raganato et al., 2020; Weiqiu et al., 2020). Selfattention-like weights can also be obtained
with our embedding by computing P(x, zi)P(x, zi)> for each reference i. In that sense, our work is related to recent research on efficient self-attention (Wang et al., 2020; Choromanski et al., 2020), where a low-rank approximation of the self-attention matrix is computed.
Position smoothing and relative positional encoding. Transformers can add an absolute positional encoding to the input features (Vaswani et al., 2017). Yet, relative positional encoding (Dai et al., 2019) is a current standard for integrating positional information: the position offset between the query element and a given key can be injected in the attention score (Tsai et al., 2019), which is equivalent to our approach. The link between CKNs and our kernel, provided by this positional encoding, stands in line with recent works casting attention and convolution into a unified framework (Andreoli, 2019). In particular, Cordonnier et al. (2020) show that attention learns convolution in the setting of image classification: the kernel pattern is learned at the same time as the filters.
Multiple references and attention heads. In the transformer architecture, the succession of blocks composed of an attention layer followed by a fully-connected layer is called a head, with each head potentially focusing on different parts of the input. Successful architectures have a few heads in parallel. The outputs of the heads are then aggregated to output a final embedding. A layer of our embedding with non-linear kernel κ can be seen as such a block, with the references playing the role of the heads. As some recent works question the role of attention heads (Voita et al., 2019; Michel et al., 2019), exploring the content of our learned references z may provide another perspective on this question. More generally, visualization and interpretation of the learned references could be of interest for biological sequences.

heading: 5 EXPERIMENTS
text: We now show the effectiveness of our embedding OTKE in tasks where samples can be expressed as large sets with potentially few labels, such as in bioinformatics. We evaluate our embedding alone in unsupervised or supervised settings, or within a model in the supervised setting. We also consider NLP tasks involving shorter sequences and relatively more labels.

heading: 5.1 DATASETS, EXPERIMENTAL SETUP AND BASELINES
text: In unsupervised settings, we train a linear classifier with the cross entropy loss between true labels and predictions on top of the features provided by our embedding (where the references z and Nyström anchors w have been learned without supervision), or an unsupervised baseline. In supervised settings, the same model is initialized with our unsupervised method and then trained end-to-end (including z and w) by minimizing the same loss. We use an alternating optimization strategy to update the parameters for both SCOP and SST datasets, as used by Chen et al. (2019a;b). We train for 100 epochs with Adam on both data sets: the initial learning rate is 0.01, and get halved as long as there is no decrease in the validation loss for 5 epochs. The hyper-parameters we tuned include number of supports and references p, q, entropic regularization in OT ε, the bandwidth of Gaussian kernels and the regularization parameter of the linear classifier. The best values of ε and the bandwidth were found stable across tasks, while the regularization parameter needed to be more carefully cross-validated. Additional results and implementation details can be found in Appendix C.
Protein fold classification on SCOP 1.75. We follow the protocol described by Hou et al. (2019) for this important task in bioinformatics. The dataset contains 19, 245 sequences from 1, 195 different classes of fold (hence less than 20 labels in average per class). The sequence lengths vary from tens to thousands. Each element of a sequence is a 45-dimensional vector. The objective is to classify the sequences to fold classes, which corresponds to a multiclass classification problem. The features fed to the linear classifier are the output of our embedding with ϕ the Gaussian kernel mapping on k-mers (subsequences of length k) with k fixed to be 10, which is known to perform well in this task (Chen et al., 2019a). The number of anchor points for Nyström method is fixed to 1024 and 512 respectively for unsupervised and supervised setting. In the unsupervised setting, we compare our method to state-of-the-art unsupervised method for this task: CKN (Chen et al., 2019a), which performs a global mean pooling in contrast to the global adaptive pooling performed by our embedding. In the supervised setting, we compare the same model to the following supervised models: CKN, Recurrent Kernel Networks (RKN) (Chen et al., 2019b), a CNN with 10 convolutional layers named DeepSF (Hou et al., 2019), Rep the Set (Skianis et al., 2020) and Set Transformer (Lee et al., 2019), using the public implementations by their authors. Rep the Set and Set Transformer are used on the top of a convolutional layer of the same filter size as CKN to extract k-mer features. Their model hyper-parameters, weight decay and learning rate are tuned in the same way as for our models
(see Appendix for details). The default architecture of Set Transformer did not perform well due to overfitting. We thus used a shallower architecture with one Induced Set Attention Block (ISAB), one Pooling by Multihead Attention (PMA) and one linear layer, similar to the one-layer architectures of CKN and our model. The results are shown in Table 2.
Detection of chromatin profiles. Predicting the chromatin features such as transcription factor (TF) binding from raw genomic sequences has been studied extensively in recent years. CNNs with max pooling operations have been shown effective for this task. Here, we consider DeepSEA dataset (Zhou & Troyanskaya, 2015) consisting in simultaneously predicting 919 chromatin profiles, which can be formulated as a multi-label classification task. DeepSEA contains 4, 935, 024 DNA sequences of length 1000 and each of them is associated with 919 different labels (chromatin profiles). Each sequence is represented as a 1000× 4 binary matrix through one-hot encoding and the objective is to predict which profiles a given sequence possesses. As this problem is very imbalanced for each profile, learning an unsupervised model could require an extremely large number of parameters. We thus only consider our supervised embedding as an adaptive pooling layer and inject it into a deep neural network, between one convolutional layer and one fully connected layer, as detailed in Appendix C.4. In our embedding, ϕ is chosen to be identity and the positional encoding described in Section 3 is used. We compare our model to a state-of-the-art CNN with 3 convolutional layers and two fully-connected layers (Zhou & Troyanskaya, 2015). The results are shown in Table 3.
Sentiment analysis on Stanford Sentiment Treebank. SST-2 (Socher et al., 2013) belongs to the NLP GLUE benchmark (Wang et al., 2019) and consists in predicting whether a movie review is positive. The dataset contains 70, 042 reviews. The test predictions need to be submitted on the GLUE leaderboard, so that we remove a portion of the training set for validation purpose and report accuracies on the actual validation set used as a test set. Our model is one layer of our embedding with ϕ a Gaussian kernel mapping with 64 Nyström filters in the supervised setting, and a linear mapping in the unsupervised setting. The features used in our model and all baselines are word vectors with dimension 768 provided by the HuggingFace implementation (Wolf et al., 2019) of the transformer BERT (Devlin et al., 2019). State-of-the-art accuracies are usually obtained after supervised fine-tuning of pre-trained transformers. Training a linear model on pre-trained features after simple pooling (e.g, mean) also yields good results. [CLS], which denotes the BERT embedding used for classification, is also a common baseline. The results are shown in Table 4.

heading: 5.2 RESULTS AND DISCUSSION
text: In protein fold classification, our embedding outperforms all baselines in both unsupervised and supervised settings. Surprisingly, our unsupervised model also achieves better results than most supervised baselines. In contrast, Set Transformer does not perform well, possibly because its implementation was not designed for sets with varying sizes, and tasks with few annotations. In detection of chromatin profiles, our model (our embedding within a deep network) has fewer layers than state-of-the-art CNNs while outperforming them, which advocates for the use of attentionbased models for such applications. Our results also suggest that positional information is important
(Appendix C.4;C.2), and our Gaussian position encoding outperforms the sinusoidal one introduced in Vaswani et al. (2017). Note that in contrast to a typical transformer, which would have stored a 1000×1000 attention matrix, our attention score with a reference of size 64 is only 1000×64, which illustrates the discussion in Section 4. In NLP, an a priori less favorable setting since sequences are shorter and there are more data, our supervised embedding gets close to a strong state-of-the-art, i.e. a fully-trained transformer. We observed our method to be much faster than RepSet, as fast as Set Transformer, yet slower than ApproxRepSet (C.3). Using the OT plan as similarity score yields better accuracies than the dot-product between the input sets and the references (see Table 2; 4). Choice of parameters. This paragraph sums up the impact of hyper-parameter choices. Experiments justifying our claims can be found in Appendix C.
• Number of references q: for biological sequences, a single reference was found to be enough in the unsupervised case, see Table 11. In the supervised setting, Table 14 suggests that using q = 5 provides slightly better results but q = 1 remains a good baseline, and that the sensitivity to number of references is moderate. • Number of supports p in a reference: Table 11 and Table 14 suggest that the sensitivity of the model to the number of supports is also moderate. • Nyström anchors: an anchor can be seen as a neuron in a feed-forward neural network (see expression of ψ in 3.3). In unsupervised settings, the more anchors, the better the approximation of the kernel matrix. Then, the performance saturates, see Table 12. In supervised settings, the optimal number of anchors points is much smaller, as also observed by Chen et al. (2019a), Fig 6. • Bandwidth σ in gaussian kernel: σ was chosen as in Chen et al. (2019b) and we did not try to optimize it in this work, as it seemed to already provide good results. Nevertheless, slightly better results can be obtained when tuning this parameter, for instance in SST-2.
OTKE and self-supervised methods. Our approach should not be positioned against selfsupervision and instead brings complementary features: the OTKE may be plugged in state-ofthe-art models pre-trained on large unannotated corpus. For instance, on SCOP 1.75, we use ESM1 (Rives et al., 2019), pretrained on 250 millions protein sequences, with mean pooling followed by a linear classifier. As we do not have the computational ressources to fine-tune ESM1-t34, we only train a linear layer on top of the extracted features. Using the same model, we replace the mean pooling by our (unsupervised) OTKE layer, and also only train the linear layer. This results in accuracy improvements as showed in Table 5. While training huge self-supervised learning models on large datasets is very effective, ESM1-t34 admits more than 2500 times more parameters than our single-layer OTKE model (260k parameters versus 670M) and our single-layer OTKE outperforms smaller versions of ESM1 (43M parameters). Finally, self-supervised pre-training of a deep model including OTKE on large data sets would be interesting for fair comparison.
Multi-layer extension. Extending the OTKE to a multi-layer embedding is a natural yet not straightforward research direction: it is not clear how to find a right definition of intermediate feature aggregation in a multi-layer OTKE model. Note that for DeepSEA, our model with single-layer OTKE already outperforms a multi-layer CNN, which suggests that a multi-layer OTKE is not always needed.

heading: ACKNOWLEDGMENTS
text: JM, GM and DC were supported by the ERC grant number 714381 (SOLARIS project) and by ANR 3IA MIAI@Grenoble Alpes, (ANR-19-P3IA-0003). AA would like to acknowledge support from the ML and Optimisation joint research initiative with the fonds AXA pour la recherche and Kamet Ventures, a Google focused award, as well as funding by the French government under management of Agence Nationale de la Recherche as part of the “Investissements d’avenir” program, reference ANR-19-P3IA-0001 (PRAIRIE 3IA Institute). DC and GM thank Laurent Jacob, Louis Martin, François-Pierre Paty and Thomas Wolf for useful discussions.

heading: A ADDITIONAL BACKGROUND
text: This section provides some background on attention and transformers, Sinkhorn’s algorithm and the relationship between optimal transport based kernels and positive definite histogram kernels.
A.1 SINKHORN’S ALGORITHM: FAST COMPUTATION OF Pκ(x, z)
Without loss of generality, we consider here κ the linear kernel. We recall that Pκ(x, z) is the solution of an optimal transport problem, which can be efficiently solved by Sinkhorn’s algorithm (Peyré & Cuturi, 2019) involving matrix multiplications only. Specifically, Sinkhorn’s algorithm is an iterative matrix scaling method that takes the opposite of the pairwise similarity matrix K with entry Kij := 〈xi, zj〉 as input C and outputs the optimal transport plan Pκ(x, z) = Sinkhorn(K, ε). Each iteration step ` performs the following updates
u(`+1) = 1/n
Ev(`) and v(`+1) =
1/p
E>u(`) , (8)
where E = eK/ε. Then the matrix diag(u(`))Ediag(v(`)) converges to Pκ(x, z) when ` tends to ∞. However when ε becomes too small, some of the elements of a matrix product Ev or E>u become null and result in a division by 0. To overcome this numerical stability issue, computing the multipliers u and v is preferred (see e.g. (Peyré & Cuturi, 2019, Remark 4.23)). This algorithm can be easily adapted to a batch of data points x, and with possibly varying lengths via a mask vector masking on the padding positions of each data point x, leading to GPU-friendly computation. More importantly, all the operations above at each step are differentiable, which enables z to be optimized through back-propagation. Consequently, this module can be injected into any deep networks.
A.2 ATTENTION AND TRANSFORMERS
We clarify the concept of attention — a mechanism yielding a context-dependent embedding for each element of x — as a special case of non-local operations (Wang et al., 2017; Buades et al., 2011), so that it is easier to understand its relationship to the OTK. Let us assume we are given a set x ∈ X of length n. A non-local operation on x is a function Φ : X 7→ X such that
Φ(x)i = n∑ j=1 w(xi,xj)v(xj) = W(x) > i V(x),
where W(x)i denotes the i-th column of W(x), a weighting function, and V(x) = [v(x1), . . . , v(xn)]
>, an embedding. In contrast to operations on local neighborhood such as convolutions, non-local operations theoretically account for long range dependencies between elements in the set. In attention and the context of neural networks, w is a learned function reflecting the relevance of each other elements xj with respect to the element xi being embedded and given the task at hand. In the context of the paper, we compare to a type of attention coined as dot-product selfattention, which can typically be found in the encoder part of the transformer architecture (Vaswani et al., 2017). Transformers are neural network models relying mostly on a succession of an attention layer followed by a fully-connected layer. Transformers can be used in sequence-to-sequence tasks — in this setting, they have an encoder with self-attention and a decoder part with a variant of self-attention —, or in sequence to label tasks, with only the encoder part. The paper deals with the latter. The name self-attention means that the attention is computed using a dot-product of linear
transformations of xi and xj , and x attends to itself only. In its matrix formulation, dot-product self-attention is a non-local operation whose matching vector is
W(x)i = Softmax ( WQxix
>W>K√ dk
) ,
where WQ ∈ Rn×dk and WK ∈ Rn×dk are learned by the network. In order to know which xj are relevant to xi, the network computes scores between a query for xi (WQxi) and keys of all the elements of x (WKx). The softmax turns the scores into a weight vector in the simplex. Moreover, a linear mapping V(x) = WV x, the values, is also learned. WQ and WK are often shared (Kitaev et al., 2020). A drawback of such attention is that for a sequence of length n, the model has to store an attention matrix W with size O(n2). More details can be found in Vaswani et al. (2017).

heading: B PROOFS
text: B.1 PROOF OF LEMMA 3.1 Proof. First, since ∑n′ j=1 pP(x ′, z)jk = 1 for any k, we have
W2(x, z) 2 = n∑ i=1 p∑ k=1 P(x, z)ikd 2 κ(xi, zk)
= n∑ i=1 p∑ k=1 n′∑ j=1 pP(x′, z)jkP(x, z)ikd 2 κ(xi, zk)
= ‖u‖22,
with u a vector in Rnn′p whose entries are √ pP(x′, z)jkP(x, z)ikdκ(xi, zk) for i = 1, . . . , n, j = 1, . . . , n′ and k = 1, . . . , p. We can also rewrite W z2 (x,x ′) as an `2-norm of a vector v in
Rnn′p whose entries are √ pP(x′, z)jkP(x, z)ikdκ(xi,x ′ j). Then by Minkowski inequality for the `2-norm, we have
|W2(x, z)−W z2 (x,x′)| = |‖u‖2 − ‖v‖2| ≤ ‖u− v‖2
=  n∑ i=1 p∑ k=1 n′∑ j=1 pP(x′, z)jkP(x, z)ik(dκ(xi, zk)− dκ(xi,x′j))2 1/2
≤  n∑ i=1 p∑ k=1 n′∑ j=1 pP(x′, z)jkP(x, z)ikd 2 κ(x ′ j , zk) 1/2
=  p∑ k=1 n′∑ j=1 P(x′, z)jkd 2 κ(x ′ j , zk) 1/2 = W2(x ′, z),
where the second inequality is the triangle inequality for the distance dκ. Finally, we have
|W2(x,x′)−W z2 (x,x′)| ≤|W2(x,x′)−W2(x, z)|+ |W2(x, z)−W z2 (x,x′)| ≤W2(x′, z) +W2(x′, z) =2W2(x ′, z),
where the second inequality is the triangle inequality for the 2-Wasserstein distance. By symmetry, we also have |W2(x,x′)−W z2 (x,x′)| ≤ 2W2(x, z), which concludes the proof.
B.2 RELATIONSHIP BETWEEN W2 AND W z2 FOR MULTIPLE REFERENCES
Using the relation prooved in Appendix B.1, we can obtain a bound on the error term between W2 and W z2 for a data set of m samples (x 1, . . . ,xm) and q references (z1, . . . , zq)
E2 := 1 m2 m∑ i,j=1 |W2(xi,xj)−W z 1,...,zq 2 (x i,xj)|2 ≤ 4 mq m∑ i=1 q∑ j=1 W 22 (x i, zj). (9)
When q = 1, the right-hand term in the inequality is the objective to minimize in the Wasserstein barycenter problem (Cuturi & Doucet, 2013), which further explains why we considered it: Once W z2 is close to the Wasserstein distance W2, Kz will also be close to KOT. We extend here the bound in equation 6 in the case of one reference to the multiple-reference case. The approximate 2-Wasserstein distance W z2 (x,x ′) thus becomes
W z 1,...,zq 2 (x,x ′) :=
〈 1
q q∑ j=1 Pzj (x,x ′), d2κ(x,x ′)
〉1/2 = 1 q q∑ j=1 W z j 2 (x,x ′)2 1/2 . Then by Minkowski inequality for the `2-norm we have
|W2(x,x′)−W z 1,...,zq 2 (x,x ′)| = ∣∣∣∣∣∣∣ 1 q q∑ j=1 W2(x,x ′)2 1/2 − 1 q q∑ j=1 W z j 2 (x,x ′)2 1/2 ∣∣∣∣∣∣∣
≤ 1 q q∑ j=1 (W2(x,x ′)−W z j 2 (x,x ′))2 1/2 , and by equation 6 we have
|W2(x,x′)−W z 1,...,zq 2 (x,x ′)| ≤ 4 q q∑ j=1 min(W2(x, z j),W2(x ′, zj))2 1/2 . Finally the approximation error in terms of Frobenius is bounded by
E2 := 1 m2 m∑ i,j=1 |W2(xi,xj)−W z 1,...,zq 2 (x i,xj)|2 ≤ 4 mq m∑ i=1 q∑ j=1 W 22 (x i, zj).
In particular, when q = 1 that is the case of single reference, we have
E2 ≤ 4 m m∑ i=1 W 22 (x i, z),
where the right term equals to the objective of the Wasserstein barycenter problem, which justifies the choice of z when learning without supervision.

heading: C ADDITIONAL EXPERIMENTS AND SETUP DETAILS
text: This section contains additional experiments on CIFAR-10, whose purpose is to illustrate the kernel associated with our embedding with respect to other classical or optimal transport based kernels, and test our embedding on another data modality; additional results for the experiments of the main section; details on our setup, in particular hyper-parameter tuning for our methods and the baselines.
C.1 EXPERIMENTS ON KERNEL MATRICES (ONLY FOR SMALL DATA SETS).
Here, we compare the optimal transport kernel KOT (3) and its surrogate Kz (2) (with z learned without supervision) to common and other OT kernels. Although our embedding Φz is scalable, the exact kernel require the computation of Gram matrices. For this toy experiment, we therefore use
5000 samples only of CIFAR-10 (images with 32 × 32 pixels), encoded without supervision using a two-layer convolutional kernel network (Mairal, 2016). The resulting features are 3 × 3 patches living in Rd with d = 256 or 8192. KOT and Kz aggregate existing features depending on the ground cost defined by −κ (Gaussian kernel) given the computed weight matrix P. In that sense, we can say that these kernels work as an adaptive pooling. We therefore compare it to kernel matrices corresponding to mean pooling and no pooling at all (linear kernel). We also compare to a recent positive definite and fast optimal transport based kernel, the Sliced Wasserstein Kernel (Kolouri et al., 2016) with 10, 100 and 1000 projection directions. We add a positional encoding to it so as to have a fair comparison with our kernels. A linear classifier is trained from this matrices. Although we cannot prove that KOT is positive definite, the classifier trained on the kernel matrix converges when ε is not too small. The results can be seen on Table 6. Without positional information, our kernels do better than Mean pooling. When the positions are encoded, the Linear kernel is also outperformed. Note that including positions in Mean pooling and Linear kernel means interpolating between these two kernels: in the Linear kernel, only patches with same index are compared while in the Mean pooling, all patches are compared. All interpolations did worse than the Linear kernel. The runtimes illustrate the scalability of Kz.
C.2 CIFAR-10
Here, we test our embedding on the same data modality: we use CIFAR-10 features, i.e., 60, 000 images with 32 × 32 pixels and 10 classes encoded using a two-layer CKN (Mairal, 2016), one of the baseline architectures for unsupervised learning of CIFAR-10, and evaluate on the standard test set. The very best configuration of the CKN yields a small number (3 × 3) of high-dimensional (16, 384) patches and an accuracy of 85.8%. We will illustrate our embedding on a configuration which performs slightly less but provides more patches (16×16), a setting for which it was designed. The input of our embedding are unsupervised features extracted from a 2-layer CKN with kernel sizes equal to 3 and 3, and Gaussian pooling size equal to 2 and 1. We consider the following configurations of the number of filters at each layer, to simulate two different input dimensions for our embedding:
• 64 filters at first and 256 at second layer, which yields a 16 × 16 × 256 representation for each image.
• 256 filters at first and 1024 at second layer, which yields a 16 × 16 × 1024 representation for each image.
Since the features are the output of a Gaussian embedding, κ in our embedding will be the linear kernel. The embedding is learned with one reference and various supports using K-means method described in Section 3, and compared to several classical pooling baselines, including the original CKN’s Gaussian pooling with pooling size equal to 6. The hyper-parameters are the entropic regularization ε and bandwidth for position encoding σpos. Their search grids are shown in Table 7 and the results in Table 8. Without supervision, the adaptive pooling of the CKN features by our embedding notably improves their performance. We notice that the position encoding is very important to this task, which substantially improves the performance of its counterpart without it.
C.3 PROTEIN FOLD RECOGNITION
Dataset description. Our protein fold recognition experiments consider the Structural Classification Of Proteins (SCOP) version 1.75 and 2.06. We follow the data preprocessing protocols in Hou et al. (2019), which yields a training and validation set composed of 14699 and 2013 sequences from SCOP 1.75, and a test set of 2533 sequences from SCOP 2.06. The resulting protein sequences belong to 1195 different folds, thus the problem is formulated as a multi-classification task. The input sequence is represented as a 45-dimensional vector at each amino acid. The vector consists of a 20-dimensional one-hot encoding of the sequence, a 20-dimensional position-specific scoring matrix (PSSM) representing the profile of amino acids, a 3-class secondary structure represented by a one-hot vector and a 2-class solvent accessibility. The lengths of the sequences are varying from tens to thousands.
Models setting and hyperparameters. We consider here the one-layer models followed by a global mean pooling for the baseline methods CKN (Chen et al., 2019a) and RKN (Chen et al., 2019b). We build our model on top of the one-layer CKN model, where κ can be seen as a Gaussian kernel on the k-mers in sequences. The only difference between our model and CKN is thus the pooling operation, which is given by our embedding introduced in Section 3. The bandwidth parameter of the Gaussian kernel κ on k-mers is fixed to 0.6 for unsupervised models and 0.5 for supervised models, the same as used in CKN which were selected by the accuracy on the validation set. The filter size k is fixed to 10 and different numbers of anchor points in Nyström for κ are considered in the experiments. The other hyperparameters for our embedding are the entropic regularization parameter ε, the number of supports in a reference p, the number of references q, the number of iterations for Sinkhorn’s algorithm and the regularization parameter λ in the linear classifier. The search grid for ε and λ is shown in Table 9 and they are selected by the accuracy on validation set. ε plays an important role in the performance and is observed to be stable for the same dataset. For this dataset, it is selected to be 0.5 for all the unsupervised and supervised models. The effect of other hyperparameters will be discussed below.
For the baseline methods, the accuracies of PSI-BLAST and DeepSF are taken from Hou et al. (2019). The hyperparameters for CKN and RKN can be found in Chen et al. (2019b). For Rep the Set (Skianis et al., 2020) and Set Transformer (Lee et al., 2019), we use the public implementations by the authors. These two models are used on the top of a convolutional layer of the same filter size as CKN to extract k-mer features. As the exact version of Rep the Set does not provide any implementation for back-propagation to a bottom layer of it, we consider the approximate version of
Rep the Set only, which also scales better to our dataset. The default architecture of Set Transformer did not perform well due to overfitting. We therefore used a shallower architecture with one ISAB, one PMA and one linear layer, similar to the one-layer architectures of CKN and our model. We tuned their model hyperparameters, weight decay and learning rate. The search grids for these hyperparameters are shown in Table 10.
Unsupervised embedding. The kernel embedding ϕ, which is infinite dimensional for the Gaussian kernel, is approximated with the Nyström method using K-means on 300000 k-mers extracted from the same training set as in Chen et al. (2019b). The reference measures are learned by using either K-means or Wasserstein to update centroids in 2-Wasserstein K-means on 3000 subsampled sequences for RAM-saving reason. We evaluate our model on top of features extracted from CKNs of different dimensions, representing the number of anchor points used to approximate κ. The number of iterations for Sinkhorn is fixed to 100 to ensure the convergence. The results for different combinations of q and p are provided in Table 11. Increasing the number of supports p can improve the performance and then saturate it when p is too large. On the other hand, increasing the number of references while keeping the embedding dimension (i.e. p × q) constant is not significantly helpful in this unsupervised setting. We also notice that Wasserstein Barycenter for learning the references does not outperform K-means, while the latter is faster in terms of computation.
Supervised embedding. Our supervised embedding is initialized with the unsupervised method and then trained in an alternating fashion which was also used for CKN: we use an Adam algorithm to update anchor points in Nyström and reference measures z, and the L-BFGS algorithm to optimize

heading: 1024 85.8/95.3/96.8

heading: 2048 86.6/95.9/97.2

heading: 3072 87.8/96.1/97.4
text: the classifier. The learning rate for Adam is initialized with 0.01 and halved as long as there is no decrease of the validation loss for 5 successive epochs. In practice, we notice that using a small number of Sinkhorn iterations can achieve similar performance to a large number of iteration, while being much faster to compute. We thus fix it to 10 throughout the experiments. The accuracy results are obtained by averaging on 10 runs with different seeds following the setting in Chen et al. (2019b). The results are shown in Table 13 with error bars. The effect of the number of supports q is similar to the unsupervised case, while increasing the number of references can indeed improve performance.
C.4 DETECTION OF CHROMATIN PROFILES
Dataset description. Predicting the functional effects of noncoding variants from only genomic sequences is a central task in human genetics. A fundamental step for this task is to simultaneously predict large-scale chromatin features from DNA sequences (Zhou & Troyanskaya, 2015). We consider here the DeepSEA dataset, which consists in simultaneously predicting 919 chromatin profiles including 690 transcription factor (TF) binding profiles for 160 different TFs, 125 DNase I sensitivity profiles and 104 histone-mark profiles. In total, there are 4.4 million, 8000 and 455024 samples for training, validation and test. Each sample consists of a 1000-bp DNA sequence from the human GRCh37 reference. Each sequence is represented as a 1000× 4 binary matrix using onehot encoding on DNA characters. The dataset is available at http://deepsea.princeton. edu/media/code/deepsea_train_bundle.v0.9.tar.gz. Note that the labels for each profile are very imbalanced in this task with few positive samples. For this reason, learning unsu-
pervised models could be intractable as they may require an extremely large number of parameters if junk or redundant sequences cannot be filtered out.
Model architecture and hyperparameters. For the above reason and fair comparison, we use here our supervised embedding as a module in Deep NNs. The architecture of our model is shown in Table 15. We use an Adam optimizer with initial learning rate equal to 0.01 and halved at epoch 1, 4, 8 for 15 epochs in total. The number of iterations for Sinkhorn is fixed to 30. The whole training process takes about 30 hours on a single GTX2080TI GPU. The dropout rate is selected to be 0.4 from the grid [0.1; 0.2; 0.3; 0.4; 0.5] and the weight decay is 1e-06, the same as Zhou & Troyanskaya (2015). The σpos for position encoding is selected to be 0.1, by the validation accuracy on the grid [0.05; 0.1; 0.2; 0.3; 0.4; 0.5]. The checkpoint with the best validation accuracy is used to evaluate on the test set. Area under ROC (auROC) and area under precision curve (auPRC), averaged over 919 chromatin profiles, are used to measure the performance. The hidden size d is chosen to be either 1024 or 1536.
Results and importance of position encoding. We compare our model to the state-of-the-art CNN model DeepSEA (Zhou & Troyanskaya, 2015) with 3 convolutional layers, whose best hyperparameters can be found in the corresponding paper. Our model outperforms DeepSEA, while requiring fewer layers. The positional information is known to be important in this task. To show the efficacy of our position encoding, we compare it to the sinusoidal encoding used in the original transformer (Vaswani et al., 2017). We observe that our encoding with properly tuned σpos requires fewer layers, while being interpretable from a kernel point of view. We also find that larger hidden size d performs better, as shown in Table 16. ROC and PR curves for all the chromatin profiles and stratified by transcription factors, DNase I-hypersensitive sites and histone-marks can also be found in Figure 2.
C.5 SST-2
Dataset description. The data set contains 67,349 training samples and 872 validation samples and can be found at https://gluebenchmark.com/tasks. The test set contains 1,821 samples for which the predictions need to be submitted on the GLUE leaderboard, with limited number of submissions. As a consequence, our training and validation set are extracted from the original training set (80% of the original training set is used for our training set and the remaining 20% is used for our validation set), and we report accuracies on the standard validation set, used as a test set. The reviews are padded with zeros when their length is shorter than the chosen sequence length (we choose 30 and 66, the latter being the maximum review length in the data set) and the BERT
implementation requires to add special tokens [CLS] and [SEP] at the beginning and the end of each review.
Model architecture and hyperparameters. In most transformers such as BERT, the embedding associated to the token [CLS] is used for classification and can be seen in some sense as an embedding of the review adapted to the task. The features we used are the word features provided by the BERT base-uncased version, available at https://huggingface.co/transformers/ pretrained_models.html. For this version, the dimension of the word features is 768. Our model is one layer of our embedding, with ϕ the Gaussian kernel mapping with varying number of Nyström filters in the supervised setting, and the Linear kernel in the unsupervised setting. We do not add positonnal encoding as it is already integrated in BERT features. In the unsupervised setting, the output features are used to train a large-scale linear classifier, a Pytorch linear layer. We choose the best hyper-parameters based on the accuracy of a validation set. In the supervised case, the parameters of the previous model, w and z, are optimized end-to-end. In this case, we tune the learning rate. In both case, we tune the entropic regularization parameter of optimal transport and the bandwidth of the Gaussian kernel. The parameters in the search grid are summed up in Table 18. The best entropic regularization and Gaussian kernel bandwidth are typically and respectively 3.0 and 0.5 for this data set. The supervised training process takes between half an hour for smaller models (typically 128 filters in w and 3 supports in z) and a few hours for larger models (256 filters and 100 supports) on a single GTX2080TI GPU. The hyper-parameters of the baselines were similarly tuned, see 19. Mean Pooling and [CLS] embedding did not require any tuning except for the regularization λ of the classifier, which followed the same grid as in Table 18.
Results and discussion. As explained in Section 5, our unsupervised embedding improves the BERT pre-trained features while still using a simple linear model as shown in Table 17, and its supervised counterpart enables to get even closer to the state-of-the art (for the BERT base-uncased model) accuracy, which is usually obtained after fine-tuning of the BERT model on the whole data set. This can be seen in Tables 20; 21. We also add a baseline consisting of one layer of classical self-attention, which did not do well hence was not reported in the main text.

Review summaries:
The authors propose a new way to aggregate the embeddings of elements in a set (or sequence) by comparing it with respect to (trainable) reference set(s) via Optimal Transport (OT). The motivation to build such a pooling operation is derived from self-attention and the authors suggest an OT spin to that (e.g., the different reference sets/measures can be thought of as different heads in attention). This is, however, done in a principled way with the help of kernel embeddings and not just ad-hoc using the transport plan as the attention matrix.