Keywords: Neural networks, representation learning, clustering, Wasserstein generative adversarial networks, Kantorovich duality, Monge optimal transport, Gaussian mixture model, deep generative models, auto-encoder, uncertainty estimation

A I S La dernière décennie a vu les réseaux de neurones devenir un outil de référence dans l'apprentissage statistique. En effet, cette technologie s'est imposée comme un outil essentiel pour des types de données aussi variés que les images, le texte, le son, etc. dans des contextes à grande échelle. Les succès des réseaux neuronaux s'étendent également à l'apprentissage par renforcement (jeux et robotique) et à l'apprentissage non-supervisé (analyse et génération de données), avec notamment une qualité inégalée pour l'imitation d'images avec les réseaux génératifs adversaires. Néanmoins, les réseaux neuronaux restent difficiles à interpréter en tant qu'estimateurs statistiques. L'objectif de cette thèse est d'atténuer cet inconvénient et d'accroître encore la portée des réseaux de neurones grâce à trois applications dites d'intelligence artificielle : (i) le partitionnement des données en groupes (clustering) grâce à un algorithme qu'on propose, (ii) la recherche des coordonnées pertinentes d'un signal avec une notion qu'on a baptisé « la pire distance de Wasserstein » et (iii) la prédiction d'un résultat munie de l'estimation d'une incertitude associée en revisitant et en étendant les méthodes connues.

A travers ces trois contributions, nous nous proposons de répondre à trois questions sur les représentations : (i) Comment représentons-nous les objets qui appartiennent aux groupes que nous essayons de former ? (ii) Comment représenter ce qui fait qu'un objet ressemble au reste des objets de son type ? (iii) Comment représenter une incertitude associée à une prédiction automatique ? Ainsi, ce travail présente des moyens de s'éloigner de l'utilisation supervisée classique du deep learning (le domaine de l'apprentissage statistique des réseaux de neurones) avec la volonté d'interpréter ces mystérieuses « boîtes noires » grâce à de nouveaux outils algorithmiques et statistiques. Nous avons veillé à ce que notre utilisation des réseaux de neurones soit la plus interprétable possible pour une meilleure compréhension des données en jeu, au-delà des décisions automatiques.

Mots-clefs

Réseaux de neurones, apprentissage de représentations, partionnement de données (clustering), réseaux génératifs adversaires de Wasserstein, dualité de Kantorovich, transport optimal de Monge, modèle de mélange gaussien, modèles génératifs profonds, auto-encodeur, estimation de l'incertitude i

A B S T R A C T I N E N G L I S H

The last decade has seen neural networks become a reference tool in statistical learning. Indeed, this technology has established itself as an essential tool for data types as varied as images, text, sound, etc. in large scales settings. The success of neural networks also extends to reinforcement learning (games and robotics) and unsupervised learning (data analysis and generation) including unparalleled quality for image imitation with generative adversary networks. Nevertheless, neural networks remain difficult to interpret as statistical estimators. The aim of this thesis is to mitigate this drawback and further enhance the scope of neural networks through three so-called artificial intelligence applications: (i) partitioning data into groups (clustering) thanks to a proposed algorithm, (ii) finding the relevant coordinates of a signal with a notion that we dub "the worst Wasserstein distance" and (iii) predicting a result while estimating an associated uncertainty by revisiting and extending known methods.

Through these three contributions, we propose to answer three questions about representations: (i) How do we represent objects that belong to the groups we are trying to form? (ii) How do we represent what makes an object resemble its kind? (iii) How do we represent an uncertainty associated with an automatic prediction? Thus, this work presents ways to get away from the classic supervised use of deep learning (the statisical learning field related to neural networks) with the desire to interpret these mysterious so-called black boxes thanks to new algorithmic and statistical tools. We tried hard to make sure that our use of neural networks is as interpretable as possible for a better understanding of the data at stake, beyond automatic decisions.

Through perception and experience, each human being gathers data all the time and eventually process that data into representations that become our grasp onto the universe for interacting. We represent ideas and concepts for thinking and language and even get pre-representations of things we do not know. Although artificial intelligence has become a very popular expression, one can remark that representation is a key prerequisite for any of processing that we consider whether we talk about our intelligence as human beings our artificial intelligence for automatic machines. On one side of the spectrum of possibilities, excellent representations do not even require further processing and simple decision mechanisms on raw data are enough for good results. On the other side of the spectrum of possibilities, poor representations actually do require sophisticated expert-based or advanced statistical processing to output good results by taking into account domain-specific knowledge. During the 2010s, some major technological obstacles were crossed allowing the development of the so-called deep learning marked by the coming of a new era where the frontier between representation and processing becomes very blurry.

In many machine learning fields, defining a clear and sound objective is always key to produce good research but unfortunately, several so-called more intelligent tasks are impossible to define that well which leads us in this dissertation in statistics to insist on the statistical representation side rather than on the statistical processing side. Therefore, our main subject is representation: (i) of data among its groups in the first chapter, (ii) of data among its characteristics in the second chapter and (iii) of predictions with uncertainty in the third chapter. At the crossroads of three different fields namely statistics, deep learning and optimal transport that recently gained much scientific attention, our effort leverages that wealth of existing research to tackle representations in three contributions with different contexts: Wasserstein Clustering (DiWaC and GeWaC) to identify groups among data, Infinitesimal Wasserstein Maximal Distortion (InWaMaDi) to highlight relevant characteristics of data and Hypothesis of an Uncertainty Model (HUM) to estimate both supervised predictions and some uncertainty information.

M A C H I N E L E A R N I N G L A N D S C A P E

This dissertation focuses on neural networks (a. k. a. deep techniques) that enjoy recent tremendous sucess as a technology but much of the current described work is applicable to other tools like decision trees [START_REF] Breiman | Classification and regression trees[END_REF] or kernel-based methods [START_REF] Andrew | An Introduction to Support Vector Machines and Other Kernelbased Learning Methods[END_REF]. That being said, deep learning changed scientists traditions about data: usually researchers separated feature-extraction and automaticdecision making tasks into two jobs. Nowadays, this frontier becomes merely a blurred line as best systems are built by doing feature extraction and automatic decision simultaneously along layers (hence the "deep" adjective as more layers give more sophisticated systems). Mixing feature extraction and automatic decision in layers amplified by much more computation power than ever let neural networks become fancy again. The whole machine learning scientific community beyond those who uses neural networks refer to these low-financial-support periods as "artificial intelligence winters". In the 1990s and 2000s, even in a major conference named Neural Information Processing Systems, popular methods like kernel-based ones were more popular than neural networks at a time when the keyphrase "deep learning" was merely confidential. As announced in the seminal keynote at the International Joint Conference on neural networks in 2011, neural networks were back again in the performance leaderboards top methods thanks to recent hardware considerable computational improvements. Then, between 2011 and 2014, neural networks algorithms beat state-of-the-art records in several fields such as image processing, computer vision, speech recognition, machine translation with almost only artificial neural networks scientists and without researchers from specific domain expertise as explained by [START_REF] Ng | Deep learning[END_REF] which gave the surprising hope for an increasing ease in a growing list of applications.

During the last decade, a pleasing end-to-end paradigm emerged and says that systems should not be trained sequentially (or even independently) but rather simultaneously as a whole [START_REF] Ng | The state of artificial intelligence[END_REF]. This intuitive belief that consists in training several layers of a neural network at once is widespread for better empirical results. Unfortunately, doing an end-to-end training also means dealing with black boxes as intermediate neural networks layers that are difficult to interpret (and are even not identifiable whereas other less efficient methods still provide interpretation ease). End-to-end training seems to be encouraged for better than other known styles of training empirically and especially in Deep Learning [START_REF] Bojarski | End to end learning for self-driving cars[END_REF].

Nevertheless, this statement must be handled with care for pragmatic reasons beyond industrial scalability and loss of interpretable modularity as [START_REF] Glasmachers | Limits of end-to-end learning[END_REF] points out there are also some other effects: feeding a deep neural network with the concatenation of raw data and some non-deep-learning algorithms outputs is often hard to beat. For example, in video recognition [START_REF] Schmid | Active Large-scale Learning for Visual Recognition[END_REF][START_REF] Crasto | MARS: Motion-Augmented RGB Stream for Action Recognition[END_REF], it is recommended to augment the raw video voxels input with optical flow (which is a processed version of the same raw video pixels but for motion estimation). Indeed, because of the data starvation phenomenon (a. k. a. over-fitting), using off-the-shelves pretrained algorithms is a simplistic form of transfer learning combining knowledge (and sometimes data) from the current and from the previous tasks. Thus, pragmatically, it is sometimes useful not to follow the end-to-end-training approach just for the sake of it: sequentially trained and/or optimized modules can work very well and still provide the easiest interpretation for what each module does. One too naive end-to-endtraining approach would end up with a black box trained from scratch.

In this dissertation, there is a will to emphasize our scientific need to crack in deep learning black boxes (end-to-end or not) because that's how better data understanding gets in, beyond automatic decisions.

Dimensions

In this dissertation, we will consider large scale settings so we must be specific about what is considered large. Throughout the machine-learning-related fields, we can consider: (number of classes in classification, the output space dimension for regression and beyond, the approximate number of nodes in a grammar tree, or number of atoms of a chemical molecule graph...).

as emphasized by [START_REF] Harchaoui | Large-scale learning for image classification[END_REF] in the concept of machine learning cuboid. For each edge of this cuboid, we have direct optimization implications for feasible computations and best results so far:

• N ≫ 1 stochastic-gradient-based optimization algorithms are more suitable than in-memory alternatives because we only need a few data points (minibatches) at the same time per iteration;

• D ≫ 1 some further analysis should be conducted: dimensionality re- duction and domain-specific knowledge must be used at once for fighting against the curse of dimensionality phenomenon;

• K ≫ 1 one-class-vs-rest strategies are preferred rather than one-class-vs-one strategies for computational reasons. Indeed in a one-class-vs-rest strategy, we only need to combine K decisions separating each class with all the rest whereas in a one-class-vs-one strategy, we would have

K×(K-1) 2
decisions separating all combinations of classes pair.

Throughout this work, we are mainly interested in large N and large D data configurations without considering large K issues. For example, in clustering settings, K should be small because otherwise it defeats the data analysis purpose: having too many clusters does not help human beings to understand data. Although the large N and large D case has already been successfully investigated in the supervised classification context, at the beginning of this work (in 2016) little research had been conducted but we observe that this key preoccupation is finally entitled to scientific attention today.

Along these three dimensions of data analysis in machine learning, this Ph.D. dissertation proposes new (or revisited) representations: (i) simplifying the cardinality axis of N thanks to clustering in our first contribution, (ii) an attempt to better understand data at a coordinate level for a local dimensionality axis of D relevance assessment in our second contribution, and (iii) re-interpreting the output axis of K through uncertainty estimation in our third contribution.

Epistemology

Epistemology is the theory of knowledge [START_REF] Newman | Descartes' epistemology[END_REF][START_REF] Ahmad | The epistemology of Ibn Khald[END_REF]. In particular, in machine learning, one epistemology has consequences on our beliefs, opinions, justifications and finally our scientific methodology in our studies. Several epistemological ways to describe machine learning exist and for this dissertation we choose one with probabilistic perspectives [START_REF] Murphy | Machine Learning, a Probabilistic Perspective[END_REF] because of its ease of sophisticated interpretation for insights. We believe data is coming from a phenomenon that we call Nature that we represent by an idealized probabilistic distribution associated with a random and often multivariate variable (or a pair or a tuple of this random multivariate variables). In practice, we consider datasets as extracts of Nature. An (annotated) dataset is some collection of independent realizations that are identically (sampled) distributed from what we call Nature. Of course, this statement is falsifiable [START_REF] Bernard | Introduction à l'étude de la médecine expérimentale[END_REF] and maybe counter-intuitive but computer science, statistical learning, machine learning, data science and all these young sciences mixing mathematics and programming are just a few decades old, compared to more established several centuries old (or even millenia old) fields such as mathematics, biology, physics, chemistry, medecine etc.

Mathematically, it is convenient to choose that epistemology (rather than an other one) in order to introduce the notion of generalization capabilities of machine learning systems. Indeed, with other epistemology, datasets have a higher status and generalization becomes ill-defined (with the ad hoc notion of training error and testing generalization error). In this work, we accept that datasets empirical distributions (sum of Dirac distributions) are approximated and noisy versions of an idealized (probably smoother) distribution.

Different kinds of machine learning

Generally, for a given project, the job of a machine learning scientist (or data scientist depending on the name given by economic trends) is decomposed in two phases often in a loop: (i) training time (previous to or interleaved with a validation time) to match/imitate/reproduce the phenomenon in the presence of groundtruth information (labels, reward) for learning a model, (ii) test time (or execution time of the system we have just built) without access to groundtruth information because we are using our trained model. Following pioneers in machine learning [START_REF] Lecun | What's Wrong With Deep Learning[END_REF], we can roughly separate the machine learning landscape in three depending on what is accessible during training, validating and testing times: first supervised learning, second unsupervised learning and third reinforcement learning as this dissertation can find applications in all of these three machine learning fields.

For the purposes of notation, univariate functions are here generalized to the multivariate case by applying the associated univariate function to each entry and then concatenating everything such that the function output has the same shape as the input. For example, z ∈ R K is a vector whose K ∈ N * coordinates z (k) are indexed by k, then log(z) = log(z (1) ), . . . , log(z (k) ), . . . , log(z (K) )

⊤ .

Supervised Learning

Nature provides a pair of (input, output) random variables (x, y).

(x, y) ∼ Nature (1) collected in a training labelled dataset (or Nature extract as stated above): dataset = (x 1 , y 1 ), . . . , (x i , y i ), . . . , (x N , y N ) [START_REF] Chan | Pcanet: A simple deep learning baseline for image classification[END_REF] where each (x i , y i ) is a realization of (x, y) ∼ Nature On the one hand, input x often represents a question in various forms such as: a vector of numbers (categories, integers or floating decimal numbers), an image or a video (made of pixels, voxels in channels and beyond [START_REF] Ponce | Computer vision: a modern approach[END_REF]), a sound (its waveform or its time-frequency representation [START_REF] Li | Fundamentals of speech recognition[END_REF][START_REF] Mallat | A wavelet tour of signal processing: the sparse way[END_REF]), a gene (its ATGC or RNA-seq representation [START_REF] Barillot | Computational systems biology of cancer[END_REF]), a chemical molecule (its 3D graphical representation [START_REF] Zaslavskiy | Graph matching and its application in computer vision and computational biology R ´esum[END_REF]) etc. On the other hand, output y represents the answer of the question x in two main classes of problems: regression in which we deal real numbers and classification dealing with categories and integers. Of course, these kinds of separations are limited but somewhat useful to describe the main problems.

Many supervised learning problems share the same kind of optimization objective: min F E (x,y)∼Nature ℓ y, F (x) (3)

where ℓ y, F (x) measures the discrepancy of predicting F (x) from an input x instead of the ground truth label y. This equation Eq. (3) behaves like an aggregation of errors (when y = F (x)) summed up into one value (the lower, the better) over data. In layman's terms, ℓ(y, F (x)) is how much the system is punished for a mistake during training and is preferably zero for no error: a perfect y = F (x) scenario. More precisely, statistical learning becomes the task of finding parameters θ = θ F of function F that has the right structure (tree, random forest, linear or kernel-based support vector machines, neural networks etc.) that minimizes L:

min θ F L(θ F ) (4)
L(θ F ) = E (x,y)∼Nature ℓ y, F (x)

Indeed, the formulation of Eq. ( 4) has the merit of generalizing almost all supervised learning problems.

Beyond the scope of this dissertation, there is a considerable amount of scientific works about regularization notably inherited from Lagrangian optimization [START_REF] Boyd | Convex Optimization[END_REF] and similar techniques. In a nutshell, our Eq. ( 3) is still valid to fit in this kind of research by simply replacing the current function ℓ by function l: l(y, F (x)) = ℓ(y, F (x)) + λΩ(F )

(5)

where Ω gives the predictor function F some desirable properties (see the books of [START_REF] Bonnans | Numerical optimization: theoretical and practical aspects[END_REF] and of [START_REF] Boyd | Convex Optimization[END_REF] for further details) with a relative importance given by λ ∈ R + in order to provide generalization capabilities coping with the fact that we only have access to a limited training dataset instead of Nature itself in practice. In summary, supervised learning is finding a function F that maps x to y based on a training dataset assuming that such an idealized function F * exists (sometimes we only need F * to only be a relation and not necessarily a function):

y ≃ F * (x) (6)
In order to build an estimator F of the desired decision function F * , we usually solve an optimization problem:

F = min F L(F ) (7)
where the loss function L is a proxy of all errors that we want to ideally minimize in one value L(F ) measuring all aggregated discrepancies between the ground truth y and prediction F (x).

In practice, we only have access to a limited amount of annotated (x, y) data that we call training data to fit our prediction function F , so the L(F )s are estimated by approximations L(F ) as if the available data was all the data in the universe: almost as if the Nature smooth distribution was replaced by the dataset empirical distribution. The hope consists in saying that at test time, for a new and unseen input x coming from the same (x, y) distribution as in training time but where the true output y is unknown, we can predict an estimated output ŷ = F (x) that is close to the true output y.

Learning is possible because at training time, we have access to both the input x and output y in a dataset. Back to our philosophical considerations about epistemology, we may consider that the supervised learning task consists in compressing the relationship between the training pairs (x, y) in a fitted predictor F that once trained, one only needs x to recover a lossy version of y ≃ F (x).

During the "compression process of learning" we accept some loss in information that we sacrifice to get better compression rate in general non-lab conditions. This way, we can study many machine learning tools such linear dot products, tree, random forest of trees, non-linear kernel evaluations, vanilla neural networks, convolutional or recurrent neural networks etc. and their respective algorithms, computations and structures with the same compression-flavored point of view.

In supervised learning we have:

x and y at training time x without y at testing time that we estimate Looking at supervised learning as a compression problem is interesting for understanding recurring trade-off through out this scientific literature: models should be sophisticated enough to recover outputs information from inputs without too much loss of information (complexity and number shoud go up) while still being sufficiently sober (i. e. not too sophisticated,(complexity and number goes down) otherwise generalization capabilities dramatically drop down while compression is given up.

Unsupervised Learning

In unsupervised learning, Nature does not provide any more information than the data x itself.

x ∼ Nature (8)

This machine learning field is useful for data analysis that has more broader scientific purposes than the industrial applications of supervised learning.

x without any labeled y information that we still estimate to get structure from data for the sake of interpretable knowledge discovery Imitating data (Generative Adversarial Networks GANs for example) [START_REF] Goodfellow | Generative adversarial nets[END_REF] and clustering [START_REF] Jain | Data clustering: 50 years beyond K-means[END_REF] are two unsupervised learning tools to exhibit data analysis as a fundamentally intelligent tool for scientists, intelligent etymologically meaning from latin understanding the underlying structure of data. Unsupervised learning is often ill-posed which remains a mystery because the same automatic tasks can be evaluated subjectively many times by different human beings with still some consistency (e. g. for clustering) which gives hopes for improvements. More mathematically, indeed, we use the Hadamard definition [START_REF] Maz | Jacques Hadamard: a universal mathematician[END_REF] for a well-posed problem and we understand that all these three points cannot apply in clustering objectives optimization:

1. a solution exists, 2. the solution is unique, 3. the solution's behaviour changes continuously with the initial conditions.

In this dissertation, clustering is tackled in a chapter as the computation of the maximal optimal transports between groups (e. g. clusters) while leveraging the GANs scientific literature in terms of numerical tools. Furthermore, another chapter is trying to make use of enhanced Wasserstein distance among distributions for relevant features weighting of data without explicit supervision from any task but rather the likeness of each point with the remaining dataset points that we express with the notion of worst optimal transport once again thanks to its powerful mathematics and algorithmics machineries.

Reinforcement Learning

Reinforcement Learning [START_REF] Sutton | Reinforcement learning: An introduction[END_REF] is an area of machine learning about how agents (say robots) sequentially observing environment inputs x from Nature could take the best sequence of actions y in order to maximize some untimely reward r also given by Nature and not necessarily after each action while maintaining a state representation s (or z) of both history and environment.

(x, r) ∼ Nature, but r is not always given (i. e. we often have r = 0 ) (

For example, playing automatically chess, Go, cards are famous applications associated with artificial intelligence victories in controlled settings face to expert human beings.

In unsupervised learning, we have:

x but no supervised action y is given, only some reward r once in a while during training time

x with some reward r once in a while at testing time and we estimate the best sequence of actions y

In reinforcement learning, taking into account uncertainty estimation is certainly helpful for the distentangling the recurring exploration / exploitation dilemma [START_REF] Sutton | Reinforcement learning: An introduction[END_REF]. This opens up a potentially large range of possible applications for our contribution dealing with uncertainty.

N E U R A L N E T W O R K S

Since their introduction in computer science [START_REF] Rosenbaltt | The perceptron-a perciving and recognizing automation[END_REF], artificial neural networks loosely inspired by the biological neurons did not stop fascinating researchers until today. In a nutshell, a neural network implements a function F from R D to R K by the composition of L ∈ N * layers (or sub-functions) (F ℓ ) ℓ=1,...,L ) of the form:

(∀ℓ ∈ 1, L -1 ) F ℓ = a • Linear ℓ (10) F L = Linear L
where:

• The non-linear element-wise function a called activation is often chosen among the hyperbolic tangent function, the positive part function (or rectified linear unit ReLU), the sigmoid function (or logit function).

• The Linear ℓ s functions are matrix-vector product linear operators in the form of:

(∀x ∈ R i ℓ ) Linear ℓ (x) = M ℓ x + b ℓ (11)
for

M ℓ ∈ R o ℓ ×i ℓ and b ℓ ∈ R o ℓ with (i ℓ , o ℓ ) ∈ N * × N * and (∀ℓ ∈ 1, L - 1 ) o ℓ = i ℓ+1 (i 1 = D and o L = K)
Thanks to the universal approximation theorem of [START_REF] Hornik | Approximation capabilities of multilayer feedforward networks[END_REF], we only need mild conditions on a (unbounded and non-constant) for the set of functions expressed in this form to be dense over the set of Lebesgue-integrable functions. This means that each time, we have an objective function to minimize over a set a function, there exists a neural network F that is able to approximate the optimal solution arbitrarily well. In terms of computer science, this is good news because instead of minimizing over a set of functions, we can reasonably minimize over a set of parameters approximating that function we are looking for. Thus, we transformed a functional optimization problem into a numerical optimization problem over the matrices M ℓ s, the biases b ℓ s but also the number of layers L and the input/output parameters (i ℓ , o ℓ ). The initial enthusiasm provoked by this statement was dampened through decades because this is only an existence theorem that does not provide a way in itself to find these parameters. Moreover it turns out that the statistical estimation of these parameters is difficult due to the over-fitting phenomenon [START_REF] Scholkopf | Learning With Kernels: Support Vector Machines, Regularization, Optimization, and Beyond[END_REF] (a. k. a. data starvation which is a complementary metaphor). Optimization is also slow (especially in the 1960s, 1970s, 1980s and even 1990s compared to nowadays).

Recently, neural networks became suddenly more plausibly useful since storage and computations were getting much faster and cheaper. Indeed, on top of a dramatic computational speed improvement, fast random data access is also crucial for realistic real-world applications. During the 1990s, even with industrial-level quality results, scientists in neural networks for machine learning did not catch the whole Research community world wide attention at first. As data storage and collection problems have been nicely solved thanks to the decrease of hardware's costs (and thus the increase of available computational power especially with the rise of CPUs and GPUs parallelism and distribution), increase of data access speed and software solutions (like HDFS [START_REF] Shvachko | The Hadoop distributed file system[END_REF] and Spark [START_REF] Zaharia | Spark: Cluster computing with working sets[END_REF]) that all appeared in an era around the 2000s called Big data. In the 2010s, these tools for manipulating data were key for large-scale training and execution of machine learning engines [START_REF] Castelluccio | AI rising. Strategic Finance[END_REF]. Another factor of scientific success is the extensive use of world-wide open source repository which pioneered in terms of reproducibility and best practice sharing. On the theoretical part, the mild conditions of the universal approximation theorems [START_REF] Cybenko | Approximation by superpositions of a sigmoidal function[END_REF][START_REF] Hornik | Approximation capabilities of multilayer feedforward networks[END_REF][START_REF] Gao | Degrees of freedom in deep neural networks[END_REF] only gives approximation ability up to our statistical estimation ability (hence the need of big cardinality datasets compared to dimensionality). Otherwise, neural networks notoriously suffer from data starvation (a. k. a. overfitting) and the neural networks need appropriate structures, computations and optimization procedures to inject enough knowledge into the systems in order to get the tremendous sucess we benefit today. Still, in spite of these recent and great improvements, there is a lack of solid theoretical grounds (especially compared to its Reproducing Kernel Hilbert Space RKHS counterpart) and many open questions remain unanswered as of today although we can mention the works of [START_REF] Vidal | Mathematics of deep learning[END_REF] and [START_REF] Arora | Generalization and equilibrium in generative adversarial nets (GANs)[END_REF].

Input Data

Thanks to the universal approximation theorem presented previouly, it is reasonable to parse several kinds of input data to see how neural networks can be fed. Indeed, neural networks universal approximation theorems do only provide an existence result of a desired function but no explicit way to get or estimate it. Even with this important theorem, scientists sill have to work to build an adequate deep learning structure to cope with the estimation problem of an almost unreachably ideal neural networks parameters. Even if we had a procedure that could reach that ideal neural networks, this would be fitted on training data only which is not a guarantee for good results on unseen yet data.

Historically [START_REF] Rosenbaltt | The perceptron-a perciving and recognizing automation[END_REF], input data x ∈ R D is a vector fed to a function implemented by a one-layer perceptron: a composition of a matrix-vector product and an element-wise sigmoid function. At the same time, stochastic optimization tools were revisited, implemented for large scale settings with the associated theoretical background provided by the Robbins-Monro theorem [START_REF] Robbins | A stochastic approximation method[END_REF] (we will come back on it later). The idea of composition for more sophisticated neural networks came fast with the introduction of multilayer perceptron [START_REF] Rosenblatt | Principles of Neurodynamics: Perceptrons and the Theory of Brain Mechanisms[END_REF] along with the back-propagation algorithm which gracefully adapts the differenciation chain rule for efficient algorithms (the forward step is the evaluation of the neural network function and the backward step is the computation of its gradient with respect its parameter). Nowadays, thanks to much engineering progress, more and more artificial intelligence promises are kept.

Since the dawn of the 1990s, new kinds of information media became available on computers storage systems with increasing sophistication and some dramatic research-to-product type of improvements emerged:

T E X T Word processor softwares rapidly replaced typewriters and made natural language processing possible and one can cite the old Reuters text dataset [START_REF] Hayes | Construe-TIS: A System for Content-based Indexing of a Database of News Stories[END_REF] for example. In terms of artificial intelligence, it seems that this medium along with genetics data are the most difficult one;

I M A G E Likewise, digital recording of photographs made the singleton dataset of Lena [START_REF] Roberts | Picture Coding Using Pseudo-Random Noise[END_REF] in image processing grow from only one to several thousands images size datasets with MNIST [LeCun et al., 1989b] and to millions images size datasets such as ImageNet [START_REF] Fei-Fei | Imagenet: crowdsourcing, benchmarking and other cool things[END_REF] and beyond.

The associated tremendous research progress made possible real-world applications in everyday applications;

S O U N D Early speech recognition systems also benefit from datasets collection since one of the oldest: TIMIT [START_REF] Zue | Speech database development at MIT: Timit and beyond[END_REF], even the music processing got its own MusicNet dataset [START_REF] Thickstun | Invariances and Data Augmentation for Supervised Music Transcription[END_REF];

V I D E O the recent 2018 YouTube 8 millions videos dataset [START_REF] Abu-El-Haija | YouTube-8M: A Large-Scale Video Classification Benchmark[END_REF] seem promising for same kind of quantitative-qualitative upward gap and improvements as image in a near future.

There is an interesting hypothesis to maybe understand why the text medium is so hard to manipulate within statistical frameworks compared to the other ones although it has been the first to be digitized with enormous ever-growing quantity data: close features in the other media are more clearly dependent (two neighboring pixels in images and even voxels in videos are highly dependent like consecutive sound samples are) and this dependence fades away with longer horizons but unfortunately text features (words or even letters) have stronger and longer range interactions that suggests to see them through grammar in spite of the scientifically obsolete and harmful but widespread saying from Frederick Jelinek in 1985:

"Every time we fire a phonetician/linguist, the performance of our system goes up."

Combining strong statistical tools and linguistics is probably the best alternative for future natural language representations.

For all these different media (or information supports), a major research pattern can be analyzed: most data contain redundant information so discarding stuttered information is useful in order to manipulate the relevant infortmation only. For example, spatial data like images, have highly dependent close pixels and thus intuitively small-sized convolution kernels seem appropriate to decorrelate the redundant information. Another example is sequence data or smooth time series data where neighboring data (with respect to the sequence or temporal axes) should also be decorrelated thus 1D-convolutions [START_REF] Zhang | Character-level convolutional networks for text classification[END_REF] or recurrence seem appropriate as much as recurrent neural networks can [START_REF] Murakami | Gesture recognition using recurrent neural networks[END_REF]Taguchi, 1991, Hochreiter and[START_REF] Hochreiter | Long short-term memory[END_REF] as described in Fig. 1. This huge research effort allowed for large scale applications (both in terms of training size or cardinality N and data size or dimensionality D and in the end the output size K). Meanwhile, image and sound processing were increasingly developping methods using convolutions when LeCun et al. [1989b,a] applied those convolutions for image classification with both iterated compositions and backpropagation which turned out to be also useful for sound and even text processing [START_REF] Lecun | Convolutional networks for images, speech, and time series[END_REF]. For time series, there is still a scientific debate among the deep learning community about whether recurrent or convolutional neural networks (which abreviate to RNN and CNN respectively). Meanwhile the sequential nature of words enumeration in text is tackled as time series but there is probably more hidden structure yet to get from grammar (as pioneering work from [START_REF] Socher | Parsing natural scenes and natural language with recursive neural networks[END_REF] pointed out since 2011).

These combined breakthroughs allowed ambitious real-world applications throughout the 1990s and the 2000s. Today, new constibutions with improved engineering, structure, optimization and regularization tools are still referring to early works with even homage (e. g. GoogLeNet by [START_REF] Szegedy | Going Deeper with Convolutions[END_REF] referring to LeNet by LeCun et al. [1989b]). Image and Speech recognition is embedded in many everyday products, even video recognition begins to have reliable industrial applications. Natural language processing also gets impressive translation results1 but there is still room for improvements. For all these different media, it seems that the same phenomenon occurs: unleashing clean dataset with thorough engineering effort astonishingly helps the scientific community to bring back high quality prototypes and ultimately products. Indeed, if we measure the time delay between a dataset release and available products and the best example is the AdaBoost implementation for face detection with Haar features (MIT-CMU frontal faces dataset [START_REF] Sung | MIT+ CMU frontal face dataset a, b and c[END_REF]] and prize-wining paper [START_REF] Viola | Robust real-time face detection[END_REF] less than 3 years after). In computer vision, the same phenomenon occured with more and more available large cardinality datasets and thus ready-to-use products as Fig. 2 shows. Now we present a table of tips to adapt unconstrained vanilla neural networks functions to constrained custom ones without the need of tediously maintaining the constraints (which is not recommended for stochastic gradient descent optimization, the weapon of choice to face large scale datasets according to the neural networks literature [START_REF] Lecun | Efficient backprop[END_REF]).

Mathematically, the universal approximation theorem allows to parse a very large class of function (more precisely we only need the Lebesgue-integrability [START_REF] Hornik | Approximation capabilities of multilayer feedforward networks[END_REF]):

x ∈ R D and F (x) ∈ R K (12)
and Table 1 enumerates many kinds of functions and their corresponding implementations thanks to the unconstrained original neural network function F .

Property Implementation

Positivity exp(F (x)) or F (x) 2 or even max(0, F (x))

Boundness between m and

M (m < M) m + (M -m)σ(F (x)) where σ(z) = 1 1+exp(-z) which is related to the SoftMax function when K = 2 Probability Vector (i. e. in K-dimensional simplex) SoftMax(F (x)) where SoftMax(z) k = exp(z k ) ∑ K ℓ=1 exp(z ℓ )
which is related to the multilogit model [START_REF] Hastie | The elements of statistical learning: data mining, inference and prediction[END_REF] and the logsumexp trick

Positivity and (semi)-definiteness matrix C(x) × C(x) ⊤ where C(x) is a lower triangular free matrix with positive diagonal entries (even bounded for more stability in practice) which is related to the Cholesky decomposition [START_REF] Golub | Matrix computations[END_REF] 1-Lipschitz function Online power iteration on each matrix-vector product inside the neural network that implements the function F which is described in the spectral normalization work by [START_REF] Miyato | Spectral normalization for generative adversarial networks[END_REF] Bijection (one-to-one function)

Composition of layers such as

x ⊤ [:d] , s(x [:d] ) × x [d:] + t(x [:d] ) ⊤ ⊤
(with pythonic indexation of coordinates) where s and t are regular neural network functions. This technique also gives to the log-determinant of the Jacobi matrix without too much computation burden thanks to the original paper of [START_REF] Dinh | Density estimation using real nvp[END_REF] Recursive Function

If t is time and y 0 an initialization:

y t+1 = G(y t , F (x t ))
which is pedagogically well presented by [START_REF] Karpathy | The Unreasonable Effectiveness of Recurrent Neural Networks[END_REF] and [START_REF] Chakraborty | Trapping effect analysis of AlGaN/InGaN/GaN Heterostructure by conductance frequency measurement[END_REF] with improved LSTMs variants by [START_REF] Hochreiter | Long short-term memory[END_REF] and GRUs [START_REF] Cho | Learning phrase representations using RNN encoder-decoder for statistical machine translation[END_REF] Table 1: Implementations for several types of Functions with respect to their inputs nature and functional properties

O P T I M A L T R A N S P O R T

The optimal transport research field has proven to be crucial in redefining our modern world as we know it, across a stunningly wide range of applications, since its French birth in the XVIIIth century [Monge, 1781], with scientists from very different backgrounds and application areas such as [START_REF] Rabin | Wasserstein barycenter and its application to texture mixing[END_REF] in image processing, [START_REF] Courty | Optimal Transport for Domain Adaptation[END_REF] in near-general data domain adaptation, [START_REF] Abouchar | Air Transport Demand, Congestion Costs, and the Theory of Optimal Airport Use[END_REF] for airports management, decisive World War II military battle victories [START_REF] Smolinski | The Scale of Soviet Industrial Establishments[END_REF], breakthrough innovation in modern economy [START_REF] Galbraith | The pragmatism of John Kenneth Galbraith[END_REF] and flabbergasting futuristic industrial revolutions such as semi-automatic objects generative design [START_REF] Shu | 3D Design Using Generative Adversarial Networks and Physics-based Validation[END_REF]. After this long one-sentence celebration, we now briefly review the basics of such a prolific mathematical offspring.

Formulations

In a data space X equipped with a metric c, we want to measure a distance between two piles of data that is related to that metric. We are willingly using the vague word pile because thanks to the notion of Dirac distibutions we have access to both smooth densities, empirical distributions and a wide variety of distributions in general. Mathematically, we thus manipulate two distributions µ and ν with associated two variables x ∼ µ and y ∼ ν both living in X . On those distributions, we want to compute the quantity W c (µ, ν) measuring how much different the piles µ and ν are. To that end, three equivalent formulations exist for that same quantity:

M O N G E F O R M U L AT I O N W c (µ, ν) = inf T * (µ)=ν E x∼µ c (x, T(x)) (13)
where T * (µ) denotes the push forward of µ by transport map T. For real understanding, we refer the reader looking for details to academic textbooks on Probability such as Random Measures, Theory and Applications [START_REF] Kallenberg | Random measures, theory and applications[END_REF].

WA S S E R S T E I N F O R M U L AT I O N W c (µ, ν) = inf γ∈Γ(µ,ν) E (x,y)∼γ c (x, y) (14)
where Γ(µ, ν) denotes the set of couplings γ based on the two distributions µ, ν: the collection of all probability measures on X × X with marginals µ and ν in order to maintain two coupling properties:

∀y ∈ X , E x∼µ [γ(x, y)] = ν(y) (15) and ∀x ∈ X , E y∼ν [γ(x, y)] = µ(x) (16) K A N T O R O V I C H -R U B I N S T E I N F O R M U L AT I O N (for the L 2 euclidean distance cost: c(x, y) = x -y 2 ) W(µ, ν) = sup C∈Lip-1 E x∼µ [C(x)] -E y∼ν [C(y)] (17)
where Lip-1 is the 1-Lipschitz functions set (from X ⊂ R D to R). More general formulas exist in Real analysis and Probability by [START_REF] Dudley | Real analysis and probability[END_REF], but changing the euclidean distance cost is possible up to how difficult redefining the Lipschitz property of C is.

As a scientist apprentice, one can notice that for producing research work, it seems that the Monge formulation is more amenable to intuition, the Wasserstein formulation is more suitable for probabilistic and geometric perspectives and the Kantorovich-Rubinstein formulation is more convenient for devising algorithms thanks the decoupling of µ and ν in the formulas which makes computations easier (simple expectations that a usual Monte Carlo estimation can handle).

Algorithms

Since the middle of the XXth century, three major algorithmical tools emerged:

1. Hungarian Discrete Method [START_REF] Kuhn | The hungarian method for the assignment problem[END_REF] 2. Entropy-Regularized Sinkhorn Fixed-Point Method [START_REF] Cuturi | Sinkhorn distances: Lightspeed computation of optimal transport[END_REF][START_REF] Genevay | Entropy-regularized optimal transport for machine learning[END_REF] 3. Wasserstein Generative Adversarial Networks [START_REF] Goodfellow | NIPS 2016 tutorial: Generative adversarial networks[END_REF] but we are well aware that with a certain amount of pragmatism it is beyond the scope of this humble state of the art to present the great mathematical achievements in optimal transport. Daring to write a summary of such a huge mathematical and on-going research field is unsettling and we refer the reader to three great references for best and rather exhaustive optimal transport overview: (i) Optimal Transport for Applied Mathematicians by [START_REF] Santambrogio | Optimal transport for applied mathematicians[END_REF] on the general scientific culture side, (ii) Computational Optimal Transport by [START_REF] Peyré | Computational optimal transport[END_REF] on the statistical and programming side, (iii) Optimal Transport: Old and New by [START_REF] Villani | Optimal transport: old and new[END_REF] on the probabilistic and theoretical side. Historically, [START_REF] Kuhn | The hungarian method for the assignment problem[END_REF] found a cubic complexity algorithm to solve an optimal transport in the discrete case which was internationally widespread as the so-called Hungarian Method. Linear cost optimization of one-to-one assignments between two sets of elements is reccuring in many real-world applications as briefly enumerated above. Indeed, once a pairwise assignment cost matrix is given, minimizing the associated cost sum with respect to the best possible assignment map boils down to the ticking of one matrix entry per row and per column (with some additional dummy entries for coping with the rectangular matrix case i. e. different cardinalities for the two sets at hand). Beyond this seminal successful attempt to cast optimal transport as a linear problem, some other works pushed the analysis further with network flows [START_REF] Ahuja | Network flows[END_REF], with graph theory angle [START_REF] Goldberg | Finding Minimum-Cost Circulations by Canceling Negative Cycles[END_REF], then Dynamic Programming mixed with fluid mechanics reasoning came in with [START_REF] Benamou | A computational fluid mechanics solution to the Monge-Kantorovich mass transfer problem[END_REF] for improved computational speed. Special discrete-continuous distributions cases were also efficiently tackled by [START_REF] Mérigot | A multiscale approach to optimal transport[END_REF] with some exceedingly fast convergence thanks to [START_REF] Lévy | A Numerical Algorithm for L2 Semi-Discrete Optimal Transport in 3D[END_REF]. First in industrial logistics, these approaches were used in a surprisingly wide range of applications from e. g. worker to work assignments, airplane to airport assignments, to communication protocol load balancing systems etc.

Recently, [START_REF] Cuturi | Sinkhorn distances: Lightspeed computation of optimal transport[END_REF] revisited entropy regularized transport to efficiently solve almost the same optimal transport problem with the Sinkhorn iterative fixed-point-type algorithm, which is especially handy when polynomial complexity is not realistic in large scale settings. The idea is that they are willing to trade some approximated optimal transport due to transport entropic regularization for realistic speed and doable computations. Surprisingly, even exploiting the entropy-regularized properties of the Sinkhorn (and thus non-optimal) transport itself has value in many applications such as robust finance [START_REF] March | Multidimensional martingale optimal transport[END_REF], ranking [Vert], photo album summarization [START_REF] Liu | Lsmi-sinkhorn: Semi-supervised squared-loss mutual information estimation with optimal transport[END_REF]... This is explainable because traditionally, regularizations schemes are meant to make numerical and stability problems vanish. That fixed-point Sinkhorn theorem gives extremely fast convergence rate of transport entropy regularized over Wasserstein distances computations: less than a dozen of iterations are enough in practice. This approach on top of dramatic computational accelerations makes it indispensable both for theoretical analysis and for many real-world applications. In spite of diligent progress for Wasserstein Generative Adversarial Networks (as we will describe later), the work accomplished by [START_REF] Genevay | Entropy-regularized optimal transport for machine learning[END_REF] still presents Sinkhornbased techniques as a great mathematical and numerical alternative for efficient optimal-transport-related solutions in machine learning.

In order to imitate high dimensional data, the principle of the milestone work of [START_REF] Goodfellow | Generative adversarial nets[END_REF] about Generative Adversarial Networks (GAN) is prototypically new for a fascinating worldwide series of research papers. As illustrated in Fig. 3, from a pseudo-random generator, we can sample some low dimensional noise that is transformed thanks to a variable generator function to get generated data within the original data space. The role of the critic function is to estimate a divergence between the real data and the generated data distributions and the generator's role is to minimize it. The most common metaphor for this mainstream press acclaimed technique is considering the generator as a forger trying to fool the detective embodied by the critic within an adversarial objective. The detective wants to distinguish generated and real data and the forger wants to produce generated data that are indistinguishable when compared to real data. More mathematically, it turns out that there is indeed a link between that min-max optimization and a Nash equilibrium as emphasized by [START_REF] Fedus | Many paths to equilibrium: Gans do not need to decrease a divergence at every step[END_REF]. Following the seminal Generative Adversarial Networks (GANs) work of [START_REF] Goodfellow | Generative adversarial nets[END_REF], [START_REF] Arjovsky | Wasserstein generative adversarial networks[END_REF] used the Kantorovich-Rubinstein formulation in order to imitate data extending GAN from Jensen-Shannon divergence minimization to Wassertein distance minimization between generated and real data distributions which paves the way for revisiting optimal transport in the context of unsupervised learning. Answering the high research expectations for GANs, the recent contribution of [START_REF] Miyato | Spectral normalization for generative adversarial networks[END_REF] called spectral normalization had a tremendous impact. Indeed, revisiting the power iteration numerical recipe [START_REF] Press | Numerical recipes 3rd edition: The art of scientific computing[END_REF] at each linear or convolutional steps to elegantly enforces the Lipschitz property required by the Kantorovich-Rubinstein duality allows variation constraints withtout unstable stochastic gradient projection techniques beyond Wasserstein distances and optimal transport. Until the spectral normalization technique, neural networks had a tendency to implement function that are not regular and enforcing the Lipschitz property (i. e. constraining the variations of the functions implemented this way) gives a beneficial regularization effect beyond Wasserstein distance estimation. This spectral normalization technique gave so much optimization stability that facing unheard-of large scale settings is made possible and provides the extraordinary images imitation results of the BigGAN approach [START_REF] Brock | Large scale GAN training for high fidelity natural image synthesis[END_REF]. This engineering achievement is also convincing thanks the already-mastered residual convolutional neural networks [START_REF] He | Deep residual learning for image recognition[END_REF] ResNet tool coming from the supervised image classification research. As a matter of fact, it is fair to say that ResNet [START_REF] He | Deep residual learning for image recognition[END_REF] gave the fascinating ability to neural networks of handling recursively raw, moderately pre-processed and highly pre-processed data at each layer by short-cutting the traditional successive-layers structure thanks to additional cross-layer connections. All combined, it made BigGAN [START_REF] Brock | Large scale GAN training for high fidelity natural image synthesis[END_REF] impressive rendering results possible, neatly fighting against the curse of dimensionality effect on the learned parameters side confronted with exceedlingly large scale unsupervised conditions for both data cardinality N (number of images in the dataset) and dimensionality D (number of pixels for the image high resolution) under unsupervised conditions.

R D R d

Critic

R E P R E S E N TAT I O N S

Knowing how to represent data is understanding data and the underlying structure beneath it (and vice versa). As is, data has in general too much dimensionality to be plotted (D > 3) and finding a useful projections seems to interleave regular dimensionality reduction tehniques and clustering (i. e. the task of building groups as we will see later).

Since 2012, several researchers identified what large scale settings for computations and memory meant not only for applied mathematics but also for the worldwide economy describing it at the Big Data era: Among them Zikoupoulos and Eaton [2016], [START_REF] Peters | The Age of Big Data[END_REF], [START_REF] Jordan | On statistics, computation and scalability[END_REF]. Consider N, the number of elements of a database (the cardinality), D the size of each element (the dimensionality), then we observe two training conditions or regimes for machine learning algorithms:

B I G D ATA R E G I M E N D ≫ 1 Statistical
theorems behave nicely but software programming was difficult before data storage and computations speed dramatically improved (even on the software side for parallelization)

S M A L L D ATA R E G I M E N D ≪ 1 or N ≃ D
As we will see, the curse of dimensionality make things very difficult in terms of statistics but computations are easy.

In practice, software and hardware issues of the big data regime have been solved in the 2000s and the beginning 2010s at an industry-quality level. Domain specific and statistical expertise are still required for small data regime. Indeed with images for example, hierarchical convolutions organized layers in neural networks has a decorrelating effect which reduces the impact of the huge dimensionality of data: intermediate results (called feature maps) are still very big on the first layers but the number of parameters has been substantially decreased thanks to the small-sized convolution kernels. For other domains, applying Convolutional Neural Networks (CNN) did work but mastery of each domain must not be ignored. In other words, injecting enough knowledge into an automatic system is basically diminishing the dimensionality (i. e. removing redundancy with respect to the task at hand) which is good news towards coping with the curse of dimensionality: this transforms a difficult small data regime into an easier big data regime.

Big data and neural networks

We previouly established that a convincing decrease of dimensionality D thanks to the appropriate representation tools is key for good empirical results in machine learning. Now it is time to explain how to optimize an objective function under large cardinality conditions. Since the 1950s, some early work mixing statistics and optimization, [START_REF] Robbins | A stochastic approximation method[END_REF] allowed large cardinality training dataset thanks to stochastic optimization. One interesting aspect of the Robbins-Monro theorem is that the objective (nor the full-gradient) does not need to be evaluated directly anymore but only a biased-free estimate of its gradient with respect to the learned parameters. In the end, the large scale constraint is relieved because only randomly picked mini-batches of data are needed instead of the whole dataset during the optimization. This early mathematical finding paved the way for contemporary large scale learning fulfilled ambitions [START_REF] Bonnans | Numerical optimization: theoretical and practical aspects[END_REF] beyond deep learning. Several scientific avenues have been taken with sometimes sophisticated algorithmical tools [START_REF] Bertsekas | Nonlinear programming[END_REF] with new programming context (e. g. distributed systems [START_REF] Hendrikx | Accelerated decentralized optimization with local updates for smooth and strongly convex objectives[END_REF] or limited memory systems [START_REF] Defazio | SAGA: A fast incremental gradient method with support for non-strongly convex composite objectives[END_REF]).

In a nutshell, the Robbins-Monro theorem allows to manipulate an idealized loss function (over all the data of the universe) without the need to compute its values nor its gradients with respect to its parameters as long as one can provide a biased-free estimator of the loss function gradient needed for the stochastic gradient descent optimization scheme. Surprisingly, it turns out that artificial neural networks actually look like biological neural networks for learning (but much less for its structure and running behavior inspite of 1950s predictions). Indeed, stochastic gradient descent follows a Hebbian rule 2 :

θ t+1 = θ t -α t f t ( 18
)
where t is an iteration index, θ the parameters that are learned, α is a learning rate (often constant) and finally f t is an incremental progress computed thanks to a measured error. Numerically, f t is the gradient of an objective L to minimize:

f t = ∇ θ L(θ t ) (19)
2 Donald O. Hebb was an influential neuro-psychologist from the 1950s which can be approximated by a biased-free estimate of ∇ θ L(θ t ) according to the Robbins-Monro theorem. For a general objective to minimize such as L(θ F ) = E (x,y)∼Nature ℓ y, F (x) the hebbian rule has a replaced stochastic gradient ft instead of the "true" one f t :

θ t+1 = θ t -α t ĝt (20)
where

ft = 1 B B ∑ b=1 ∇ θ ℓ y i b , F (x i b ) (21)
and

i b ∼ U N (1, N) is an uniform index parsing the N-cardinality training dataset in mini-batches of size B.
It is noteworthy to recall that there is a simple case where the best learning rate is given in closed form: the online estimation of a mean random multivariable which is related to the least squares problem:

µ t+1 = µ t + 1 t + 1 (x t+1 -µ t ) (22)
From a practitioner point of view, large cardinality problems have been solved both on the software/hardware side and on the optimization side. The dimensionality issues still remain and is in essence more specific to each application we tackle.

The Curse of Dimensionality

Enhancing data structure in the data representation is an efficient way to cope with the curse of dimensionality issue that we focus on in this section. As said above, in images and speech, convolutions provided enormous improvements in terms of accuracy for the given task. Word ordering structure is today considered as standard since the rise of embeddings techniques with word2vec [START_REF] Mikolov | Distributed representations ofwords and phrases and their compositionality[END_REF] and variants even with the recent BERT method. The fundamental idea is to transform observed occurrences of words ordering on large corpora into a regression problem providing relevant words representations into real-valued vectors. The situation is much more difficult in genomics despite worldwide sincere attention [START_REF] Barillot | Computational systems biology of cancer[END_REF]. Indeed 1 human is DNA-represented by D ≃ 10 12 nucleotides and the worldwide population is N ≃ 10 9 which corresponds to a small data regime once again. Even in diagnostics problem healthy v. s. ill detection or classification problem, in theory a balanced and worldwide-sized annotated DNA-dataset is not enough for classification nor interpretation which is quiet an embarassing (theoretical) scenario. We did not yet leverage enough desirable and yet unknown DNA structural information among those nucleotides to get an easier big data regime. Now we know that high cardinality is not a problem anymore if dimensionality is not too high keeping the problem under the big data regime. In small data regime when dimensionlity is too high, devising algorithms is difficult because of the well-known Curse of Dimensionality 3 that we mentioned earlier. To explain that phenomenon dubbed the Curse of Dimensionality, a classic example considers the volume of a sphere of unit radius in dimension D which follows In Fig 4, first we can see that volume V(D) reaches its maximum for D = 5 and right after, that same volume V(D) dramatically decreases towards 0 as D grows (lim D→∞ V(D) = 0) which is not intuitive and rather deceiving because of how we experience spatial neighborhood notions in our 2D and 3D living environment as human beings (we would expect that volume to grow indefinitely as it does for D = 1, 2, 3, 4). For example, in low dimensionality regimes (D = 1, 2 or 3), if a point lies close to the origin inside an unit ball neighborhood, then the volume to exhaustively parse for finding it is reasonably big (V(D) ≃ 2, 3.14 or 4.18 respectively). But for high dimensionality regimes (e. g. D = 30), the corresponding volume to parse gets extremely small (V(D) =≃ 2 × 10 -5 ) which is unsettling: we would think that looking for a point inside a higher dimensional sphere would be larger but this is not true. In the end, this means that rudimentary notions like distances or even similarities behave unexpectedly in such high dimensionality data regimes. Another pedagical example considers the ratio R(D) between a 0.9-radius and 1-radius balls' volumes (R(D) = 0.9 D ). In Fig. 5, that R(D) ratio also exhibits an embarassing phenomenon with respect to our intuition: as the dimensionality D increases, the volume ratio R(D) goes to zero 3 The keyphrase Curse of Dimensionality was first mentioned by [START_REF] Bellman | Dynamic Programming[END_REF] to describe the need of efficient algorithmic tools like Dynamic Programming [START_REF] Dasgupta | Algorithms[END_REF] to efficiently explore huge discrete solutions spaces, historically later it also concerned many kinds of large data spaces (lim D→∞ V(D) = 0) which means that in high dimensionality regimes the orange 0.1-peel occupies almost all the entire orange volume in layman's terms. For pedagogical reasons, we may present a metaphor that we dubbed the "crinkled paper in a room". Indeed, in general, real data distribution seem to behave like a crinkled sheet of paper with low intrinsic manifold dimensionality (say d = 2) inside of a large space of dimensionality D (say a d = 2 crinkled sheet of paper living in a room of D = 3 dimensions to help our 3D-living creatures intuition). In this metaphor, revealing the data structure and focusing on the independent variables ruling the data boils down to how we iron or un-crinkle that sheet of paper to explicitly have access to the real degrees of freedom parsing the data.

V(D) = π k k ! if D = 2k is even, or V(D) = 2(k!)(4π) k (2k + 1)! if D = 2k + 1 is odd (23)
Interestingly, following this metaphor helps us understand than the GAN approach is doing the opposite: starting from a low dimensional uniform law (the ironed sheet of paper) transformed by the generator into a higher dimensional and much more sophisticated law. GAN optimization is basically learning how to un-iron or crinkle a clean sheet of paper into something close to the real data manifold.

Unfortunately, in terms of measures theories, neither the Riemann nor Lebesgue measures allocate weight to that sheet of paper (crinkled or not, it is respectively not defined and zero) which may lead to explore new theoretical research avenues and maybe studying other probability foundations like maybe not yet sufficiently explored in machine learning based on the Hausdorff measure [START_REF] Abbott | Hausdorff Measures[END_REF] or using probability freed from measure theory like the acclaimed attempt of [START_REF] Breiman | Probability, volume 7 of classics in applied mathematics[END_REF]. Another way to cope with this difficulties is to abandon the notion of probability distributions to embrace energy-based models as initiated by [START_REF] Lecun | A tutorial on energy-based learning[END_REF] which require less assumptions to model data and allows to generalize a distribution into an energy (through a negexponentiation from an energy analogy from Physics but without the constraint of having a unit measure over space).

On the kernel-based techniques [START_REF] Shawe-Taylor | Kernel methods for pattern analysis[END_REF] side during the 1990s and 2000s, the dimensionality problems disappear because once a kernel similarity matrix is built, one does almost not need the data anymore to operate analysis (with supervised kernel-based support vector machines or unsupervised spectral clustering). With kernels, dimensionality problems are avoided thanks to the kernel similarity matrix. All dimensionality-wise considerations are relegated to the crucial kernel definition. Unfortunately, such techniques are limited due to the inherent quadratic memory complexity of such pairwise structures implied by the similarity matrix.

Meanwhile, on the Model-based side, naive approaches do not succeed to achieve good results because of that high dimensionality: they ultimately loose their specific model selection capabilities because of over-parametrization to match that high data dimensionality. Indeed, in reasonable dimensionality conditions, the main advantage of such techniques is their ease for probabilistic interpretation and model selection. Over-parametrization (or over-fitting) can be seen as a data starvation phenomenon: a large number of parameters to fit would require a huge amount of data to get reliable estimations which in practice leads to poor performance. We see that once again, the interesting factor is the ratio N D of cardinality N over D and not one of them without considering the other one.

In the context clustering of large scale dimensionality, parcimonious and clusterwise representations [Bouveyron and Brunet-Saumard, 2014b] circumvent these high dimensionality problems and still keep the appealing probabilistic properties of model-based clustering without sacrificing accuracy. One can remark that in the supervised classification literature, sparsity (and even structured sparsity [START_REF] Jenatton | Structured sparsity-inducing norms: Statistical and algorithmic properties with applications to neuroimaging[END_REF]) also did cope with dimensionality problems that are similar in essence.

Model-based clustering algorithms are popular because they are renowned for their probabilistic foundations and their flexibility [START_REF] Duda | Pattern classification[END_REF]. Indeed, even for non-statisticians, the possibility to output meaningful probabilities is intuitive and principled. The main drawback of mixture-based and model-based methods for clustering is the lack of richness (in Kleinberg's sense see [START_REF] Kleinberg | An impossibility theorem for clustering[END_REF] but we will come back on it later) due to necessary distribution assumptions that may not be necessarily true for real data which justifies our attempt to alleviate this limitation thanks to the functional expressivity of neural networks.

One fundamental machine learning hypothesis is recurring in the literature [START_REF] Murphy | Machine Learning, a Probabilistic Perspective[END_REF][START_REF] Duda | Pattern classification[END_REF][START_REF] Bishop | Pattern recognition[END_REF]: real data live in a low-dimensional (of dimension D) manifold in a much higher dimensional space (of dimension D and d ≪ D). A classic pedagogical example consists of the independent and uniform sampling of each pixel of an image: there is no realistic chance to produce a convincing photograph! This means that even sophisticated mathematical object such as photographs lie on a manifold of lower intrinsic dimensionality than the number of pixels multiplied by the number of channels. In a reverse fasion, this has been confirmed by the DCGAN work of [START_REF] Radford | Unsupervised representation learning with deep convolutional generative adversarial networks[END_REF] that is able to generate D = 3 × 256 × 256 ≃ 2 × 10 5 convincing DCGAN images from a random uniform variable made of d = 100 independent coordinates). Thanks to this low-manifold-dimensionality hypothesis for data in mind, it is reasonable to investigate some dimensionality reduction techniques.

Dimensionality Reduction

At the beginning of the XXth century, Principal Component Analysis (PCA) was invented [Pearson, 1901]. This technique finds an optimal linear (or affine) projection with respect to compression/decompression quadratic reconstruction error. This algorithm gave birth to two more recent ones: (i) its kernelized extension [START_REF] Schölkopf | Nonlinear component analysis as a kernel eigenvalue problem[END_REF]] (euclidean distances can be expressed with dot products that are in turn replaced by kernel evaluations in a Reproducing Kernel Hilbert Space RKHS following the well-known kernel trick) and (ii) autoencoders [START_REF] Kramer | Nonlinear principal component analysis using autoassociative neural networks[END_REF][START_REF] Bourlard | Auto-association by multilayer perceptrons and singular value decomposition[END_REF][START_REF] Vincent | Stacked denoising autoencoders: Learning useful representations in a deep network with a local denoising criterion[END_REF] which replaces compression and decompression by one neural network each. Fig. 7, 8 and 9 show PCA, (Gaussian) kernel-PCA and AE dealing with nonlinearly-structured yet simple 2D distributions shown in Fig. 6. It turns out that (up to the Gaussian parameter of our kernel-PCA), the non-linearity improvements of PCA in two different variants namely kernel-PCA and AE does help ironing the crinkled distribution of interest (to follow our metaphor in section 1.5.2). Dimension reduction approaches such as principal components analysis (PCA) or even autoe-encoders (AE) [START_REF] Vincent | Stacked denoising autoencoders: Learning useful representations in a deep network with a local denoising criterion[END_REF] may help for clustering but, as is, are not designed with a clustering mindset which causes poor results in practice. These global dimension reduction techniques are pragmatic but ignore information which is discriminant for separating clusters. Indeed, clusters are usually living in different sub-spaces between clusters if there exist. Back in the original data space, there is no reason to find an easy-to-find common linear sub-space that is discriminant enough to separate all the classes at the same time. For example, if clustering is taken into account while reducing the dimensionality, then one solution could be to divide the reduced space into as many zones as clusters such that they do not overlap while still reducing the dimensionality. Thus, we avoid generic approaches because they cannot afford by themselves to capture these subtleties in the data. One must combine dimensionality reduction and clustering. Clustering could be looked at as an extremely simplified version of the data by just keeping the index of the cluster the data belong to.

The high dimensionality clustering literature [START_REF] Bouveyron | High-Dimensional Data Clustering[END_REF] tends to show that clustering and dimensionality should be done a the same time (i.e. in an end-to-end fashion) and not sequentially. Indeed, on the one hand, doing clustering first for huge dimensionality data is computationally difficult for obvious reasons and also statistically difficult because of the curse of dimensionality.

On the other hand, doing dimensionality reduction first looses hidden clusterwise information about data. One major issue of this work is precisely trying to tackle this "chicken or egg" problem.

Vanilla Auto-encoders alone do not allow to specify a precise probabilistic structure for a low-dimensional representation. This limits their combination with model-based clustering techniques. Furthermore, optimizing an auto-encoder and a Gaussian mixture generally implies the use of a trade-off hyper-parameter to combine these two objectives. This hyper-parameter is possibly hard to tune as cross-validation is not an option in our unsupervised settings as no validation score can be used by definition.

The problem of learning representations from data in an unsupervised manner is a long-standing problem in machine learning [START_REF] Bengio | Representation learning: A review and new perspectives[END_REF], LeCun et al., 2015]. Principal Components analysis (PCA) and auto-encoders (AE) which can be seen as non-linear extension of PCA [Baldi and [START_REF] Curran Associates | Neural networks and principal component analysis: Learning from examples without local minima[END_REF] have been used for representing faces [START_REF] Turk | Face recognition using eigenfaces[END_REF] or to produce a hierarchy of features [START_REF] Chan | Pcanet: A simple deep learning baseline for image classification[END_REF]. Other techniques have been used such as sparse coding [START_REF] Mairal | Sparse representation for color image restoration[END_REF] where the representation of one image is a linear combination of a few elements in a dictionary of features. More recently [START_REF] Bojanowski | Unsupervised learning by predicting noise[END_REF] learned features unsupervisedly by a procedure that consists in mapping a large collection of images to noise vectors through a deep convolutional neural networks.

Clustering and dimensionality reduction are interleaved. The importance of finding a suitable representation for unsupervised tasks was first highlighted by [START_REF] Chang | On using principal components before separating a mixture of two multivariate normal distributions[END_REF], who showed that embeddings based on principal component analysis were often unfit for clustering purposes so we suggest the idea of learning both clustering and dimensionality at the same time in an end-to-end deep learning fahsion. In a more model-based literature [START_REF] Bouveyron | Model-Based Clustering and Classification for Data Science: With Applications in R[END_REF], combining clustering and dimensionality reduction simultaneously also proved more successful than separating dimensionality reduction and clustering sequentially, which in turn, was already more successful than doing only one of them for both results. This means that both deep learning and Bayesian literatures tend to show a certain symbiosis between clustering and dimensionality reduction towards data analysis and understanding.

In the context of linear embeddings (that offers dimensionality reduction), the main approach was to combine linear discriminant analysis with the k-Means (k-Means) algorithm (DisKMeans) [De la [START_REF] De La Torre | Discriminative cluster analysis[END_REF] or more generally a Gaussian Mixture Model (Fisher-EM) [START_REF] Bouveyron | Simultaneous model-based clustering and visualization in the fisher discriminative subspace[END_REF]. Much less research is available in relation to non-linear embeddings. [START_REF] Archambeau | Manifold constrained variational mixtures[END_REF] however proposed to use manifold learning in combination with GMM. Combining clustering with representation learning has been done with deep learning techniques in the past. An early attempt was explored by [START_REF] Trigeorgis | A deep semi-nmf model for learning hidden representations[END_REF] who used a deep semi-non-negative-matrix-factorization (NMF) model to specifically factorize the input into multiple stacking factors which are initialized and updated layer by layer with k-Means on the last layer.

Neural networks have proven successful in the context of supervised classification and even regression [START_REF] Goodfellow | Deep Learning[END_REF]. Indeed, their ability to transform data such that the frontiers between classes are hyperplanes in the classification setting have made them very popular. In spite of the non-convexity of their optimization scheme, today, they are superior to convex machineries such as Support Vector Machines even for kernelized ones in almost every domain in Computer Vision and sound processing for example. The idea of having learned features has already been tackled by [START_REF] Chen | Deep learning with nonparametric clustering[END_REF] used Deep Belief Networks together with maximum-margin clustering. [START_REF] Wang | Learning a task-specific deep architecture for clustering[END_REF] jointly optimized a sparse coding objective and a clustering loss. Eventually, all these recent approaches have been empirically outperformed by auto-encoders' style machineries.

When it comes to compressing data while limiting loss of reconstruction information, auto-encoders have proved efficient [START_REF] Vincent | Stacked denoising autoencoders: Learning useful representations in a deep network with a local denoising criterion[END_REF]. Briefly, an auto-encoder is a neural network made of two parts: (i) the encoder maps the data in a low-dimension space, (ii) the decoder maps them back to the original space. An auto-encoder is trained to reconstruct the data in the original space (usually in a least squares fashion but it could be any differentiable metric). At the end, if the reconstruction error is low, then codes resulting from the encoder (also called "bottleneck") have compressed data without loosing too much information (because by construction, it is possible to rebuild data from codes thanks to the decoder). The main assumption behind this technique is that the input data space of high dimensionality contains structure that could be successfully embedded in a lower-dimensionality manifold [START_REF] Alain | What regularized auto-encoders learn from the datagenerating distribution[END_REF]Bengio, 2014, Sonoda and[START_REF] Sonoda | Decoding stacked denoising autoencoders[END_REF] and the code space plays that embedding role.

Generative Adversarial Networks (GAN) [START_REF] Goodfellow | Generative adversarial nets[END_REF]] establish a min-max game between a generator neural network on one side and a discriminator or critic neural network on the other side in order to generate data (from random noise) that the critic cannot distinguish from the real data. From that influential work emerged Adversarial Auto-Encoders (AAE) by [START_REF] Makhzani | Adversarial autoencoders[END_REF], Wasserstein Auto-Encoders (WAE) by [START_REF] Tolstikhin | Wasserstein auto-encoders[END_REF] and Adversarially Learned Inference (ALI) by [START_REF] Dumoulin | Adversarially learned inference[END_REF]. In a few words, thrice are turning an auto-encoder into a generative model. They are trained in different ways that put an arbitrary fixed prior distribution in the code space. For the clustering chapter, we were initially inspired by these approaches with learned mixture distribution instead of a fixed prior one.

D I S S E R TAT I O N O U T L I N E

This thesis is illustrating the process of learning representations using neural networks and optimal transport through three applications:

C L U S T E R I N G joint work with Pierre-Alexandre Mattei, Andrés Almansa and Charles

Bouveyron It is about unsupervisedly representing large-scale datasets in groups. This offers data tools to get a better intimate knowledge of the data with whereas the usual deep learning supervised classification algorithms do not so easily unless tedious manual annotation is already paid for at least on training data;

U N S U P E R V I S E D F E AT U R E I M P O R TA N C E
It consists in analyzing data at a coordinates level with a wide of applications from pure data understanding to background/foreground image segmentation in an unsupervised manner (the only remaining supervision being a pile of images containing the same semantic class of content). In this work, we only propose an attempt to accomplish that desirable goal and we provide a sound theoretical framework to do it;

P R E D I C T I O N W I T H U N C E R TA I N T Y
We insist on a better interpretation of supervisedly trained neural networks output in terms of uncertainty (especially for classification probabilities) towards a simple yet efficient way to improve uncertainty estimation in such supervised learning scenarios. In this contribution, sensitive applications can be a little more reliably envisioned when simple industrial constraints or more complex health, security and even justice issues are involved.

At the end of this mauscript, we provide an appendix Generative Adversarial Networks Initialization with Auto-Encoders. This is an heuristic for practical initialization of GAN training with tips revisiting pre-training traditions from the 1990s but for contemporary machine learning tools.

In this state of the art chapter, we introduced the challenges related to thrice learning representations, optimal transport and neural networks and the next chapters will be devoted to applications in clustering, unsupervised feature importance extraction and supervised uncertainty estimation thanks to the aforementioned tools. Some ongoing and future work are envisioned as a conclusion and an appendix presents practical yet effective techniques for training GANs that we gathered from experience.

Clustering

Clustering is one of the oldest unsupervised learning task [START_REF] Jain | Data clustering: 50 years beyond K-means[END_REF]. Clustering [START_REF] Duda | Pattern classification[END_REF] is the task of making groups without the need of any manual annotations. Along with dimensionality reduction, clustering is a desirable goal in data analysis, visualization and is often a preliminary step in many algorithms for example in computer vision [START_REF] Ponce | Computer vision: a modern approach[END_REF] and natural language processing [START_REF] Goldberg | Neural network methods for natural language processing[END_REF]. Clustering and more generally data analysis does not only consist in pre-processing steps, it is about helping us (as human beings) understanding the underlying structure of data at hand. Meanwhile, the computer vision field has recently witnessed major progress thanks to end-to-end deep-learning systems since the seminal work of LeCun et al. [1990] and more recently of [START_REF] Krizhevsky | Imagenet classification with deep convolutional neural networks[END_REF]. Most of the work however has been carried out in a supervised context. Our effort leverages that wealth of existing research but in an unsupervised framework.

While optimal transport [START_REF] Villani | Optimal transport: old and new[END_REF] have gained recent attention especially for generating data (i.e. imitating data) in large scale settings (large both in terms of dataset cardinality N and dimensionality D) with Generative Adversarial Networks (GAN) originated by [START_REF] Arjovsky | Wasserstein generative adversarial networks[END_REF] and Sinkhorn divergences by [START_REF] Genevay | Entropy-regularized optimal transport for machine learning[END_REF], we chose to ignore imitating capabilities and just use this literature to algorithmically manipulate Wasserstein distances. For example, this considerable amount of anterior work gives us a significant ease for optimization with helpful tools such as stochastic gradient descent.

The purpose of this research is to build a linear-complexity algorithms that use non-linear embeddings into code spaces. Indeed, in the clustering literature, one can distinguish two kinds of clustering algorithms with respect to their computation and memory complexity as function of the cardinality N. On one side, we have linear algorithms such as k-Means (k-Means) and Gaussian Mixture Models (GMM), which usually work directly on the data (i.e. without any medium such as embeddings and transformed version of the raw data). On the other side, we also have quadratic and cubic algorithms such as hierarchical clustering [START_REF] Duda | Pattern classification[END_REF] and spectral clustering [START_REF] Ng | On spectral clustering: Analysis and an algorithm[END_REF][START_REF] Zelnik-Manor | Self-tuning spectral clustering[END_REF][START_REF] Luxburg | A tutorial on spectral clustering[END_REF] that use pairwise similarities to emphasize the latent clustering structure lying on the data. Now, we describe some statistical problems related to the clustering task, and we will enumerate some famous clustering algorithms.

Clustering is an ill-posed problem

In general, there exists no clear, objective means of defining a "good clustering". For a fixed number of groups, [START_REF] Kleinberg | An impossibility theorem for clustering[END_REF] presents three desirable properties for a given clustering algorithm, namely: and he prooved the impossible existence of such an algorithm featuring all these three properties simultaneously. In other words, his clustering impossibility theorem tells us that clustering is an ill-posed problem. To add insult to injury, when data representation (or embedding) is involved, that clustering task becomes all the more unclear because the underlying metric is allowed to change arbitrarily. Indeed, the algorithm could expand the distance between points in the embedding space that initially were located near to each another, which would break initial pairwise "distance" constraints between the initial points and would inevitably violate the internal structure relating data. Well aware about these difficulties, we decided to try anyway following our scientific predecessors as clustering is useful in practice "as is".

Taking advantage of the abundant optimal transport literature with probabilities (relaxing hard clusters memberships definition to prefer probabilities) and also the neural networks literature (which successfully handles arbitrary classes shapes in supervised contexts) make our efforts reasonable towards a useful clustering algorithm for practioners. Let's review the clustering axioms of [START_REF] Kleinberg | An impossibility theorem for clustering[END_REF]. First, once a metric is chosen, scale-invariance can be given for free thanks to the geometric optimal transport interpretations i. e. all 1-Wasserstein distance would be multiplied like the data accordingly without changing the optimization results. Second, Wasserstein distances operate on all pairs of distributions to the contrary of the Kullback-Leibler divergence (which requires common support which explains the use of infinite support distributions for the models like the Gaussian), and thus no cluster shapes assumption is required, thus richness would be fulfilled by the functional expressivity power provided by neural networks. Third, unfortunately, we would not be able to achieve consistency because our models and algorithms strongly depend on the euclidean distance. In fact, this third consistency property could be partially reached thanks to a generalized notion of Wasserstein distance defined by the maximum Wasserstein distance when the metric parses a family of distances (which makes the maximum of them still a distance) but this would require further scientific work that we just skimmed in our unsupervised feature importance contribution.

Richness (i. e. free clusters shapes robustness) is obtained thanks to an intermediate space that we call embedding space or code space which has lower dimensionality that the data space. If theoretical or practical tools are given to navigate between these spaces, back and forth without loosing too much information between data and codes thanks to encoder and decoder functions, then, model-based distribution assumptions can be made on the code side which gives us the richness property in return on the data side like in the work of [START_REF] Jiang | Variational deep embedding: A generative approach to clustering[END_REF]. Indeed, at the beginning of this thesis (mid-2016), our first intuition was to put a mixture distribution (which is a typical model-based idea) at the bottleneck of an auto-encoder (which is a typical deep learning unsupervised tool) with the hope of gathering the best of these two universes: probabilistic ease for model selection coming from model-based legacy on one side and rich representations with neural networks coming from deep learning legacy on the other side. successful in several papers that we briefly present here. The first work we saw doing clustering in the code space of an auto-encoder is the one of [START_REF] Song | Deep auto-encoder based clustering[END_REF] and Yang et al. [2016a]. The idea is to cope with large data space dimensionality D ≫ 1 thanks to an intermediate code space of lower dimensionality d ≪ D. More precisely, [START_REF] Song | Deep auto-encoder based clustering[END_REF] considered a k-Means-regularized auto-encoder loss to get a code space that is more easily clustered with k-Means namely their loss is the sum of the reconstruction and the k-Means residual with a chosen hyper-parameter combining both. This philosophy is the one adopted for our own approach but with a mixture of distributions ∑ K k=1 π k × p k of distributions p k weighted by proportions π. We use a code space where clusters memberships (not data themselves but encoded versions of them) are easily computed. In our preliminary experiments, we found that optimizing the k-Means objective (online) when doing joint clustering and feature learning did not work well. We believe this is because it creates high magnitude gradients for points that are far away from cluster centers. Moreover there are sharp discontinuities at cluster boundaries whereas GMM diminishes that effect thanks to low density/probability values for far points. This empirical conclusion seems to confirm what [START_REF] Xie | Unsupervised deep embedding for clustering analysis[END_REF] also observed.

In a similar spirit, [START_REF] Huang | Deep embedding network for clustering[END_REF] have developped a locality-preserving and group-sparsity constraints method to handle the clustering. Yang et al. [2016b] alternate between supervised classification and feature learning through Convolutional neural networks (CNN) for images clustering. They significantly improved the state-of-the-art but their method is limited by its intrinsically quadratic complexity. In a similar spirit, [START_REF] Xie | Unsupervised deep embedding for clustering analysis[END_REF] embrace the t-SNE framework [START_REF] Maaten | Visualizing data using t-sne[END_REF] in a clustering context through an autoencoder in a non-model-based fashion. But doing t-SNE first and then clustering is not a good idea because of the same loss of useful cluster-wise information that occurred with PCA or auto-encoders.

All these works tend to show that simultaneous representation learning (by dimensionality reduction for example) and clustering actually do help each other. The reader will find an excellent review in the work of [START_REF] Aljalbout | Clustering with deep learning: Taxonomy and new methods[END_REF].

Clustering in Large-Scale Cardinality Regimes k-Means and Mixture Models have been studied in large scale cardinality settings [Bottou andBengio, 1995, Cappé and[START_REF] Cappé | On-line expectation-maximization algorithm for latent data models[END_REF] but these algorithms work thanks to strong distribution assumptions (like cluster-wise Gaussian clusters shapes for GMM) directly in the original data space (i.e. without embeddings). Agglomerative clustering methods greedily use a square similarity matrix to fusion data into clusters but the building of that N × N matrix is undoable for large cardinality N. Spectral clustering [START_REF] Zelnik-Manor | Self-tuning spectral clustering[END_REF] works with very mild (or no) cluster shapes assumptions thanks to the kernel trick which gives access to high (and even infinite) dimensional representations space at the cost of a square similarity matrix once again which inevitably blocks the way leading to large cardinality datasets. Nevertheless, [START_REF] Choromanska | Fast spectral clustering via the nyström method[END_REF] found a way to gracefuly alleviate this problem through to the the Nyström method that only demands the computations of only few entries of that nonstorable square similarity matrix thanks to a low-rank approximation (which is justified by the low intrinsic manifold dimensionality hypothesis related to our 2D crinkled sheet of paper in a 3D room metaphor). The present work is an attempt to provide a scalable method with the mildest possible cluster shapes assumptions thanks to neural networks that already have these two desirable properties in the supervised context: scalability and mild data assumptions.

The universal approximation theorems [START_REF] Hanin | Approximating continuous functions by relu nets of minimal width[END_REF] allow neural networks to achieve richness (in the sense of [START_REF] Kleinberg | An impossibility theorem for clustering[END_REF]) in theory. But in practice, that richness requires the estimation of a lot of parameters which is not reliable unless we have a large cardinality dataset during training compared to the dimensionality as explained above in section 1.5. In this case, large cardinality datasets are handled thanks to stochastic gradient optimization [START_REF] Bach | Beyond stochastic gradient descent for large-scale machine learning[END_REF]. Indeed, considerable research in supervised classification has been conducted based on these foundations during the last decades but this work's challenge is about extending this success to unsupervised classification (a. k. a. clustering).

Generative approaches produce a model in the form of a synthetic data distribution that is supposed to be close to the original data distribution with respect to a criterion such as the Kullback-Leibler divergence (which is equivalent to maximizing the likelihood as explained by the Pattern Recognition and Machine Learning textbook of [START_REF] Bishop | Pattern recognition[END_REF]) typically optimized with Expectation-Maximization [START_REF] Dempster | Maximum likelihood from incomplete data via the EM algorithm[END_REF]. Parameters and hyper-parameters are two different things: parameters are optimized whereas hyper-parameters are imposed before optimization and can be selected after optimization among a set of optimized models (i. e. model selection). One considerable advantage of generative techniques over others is that hyper-parameter selection is made easy through model selection thanks to principled mathematical (often Bayesian) foundations. Inddeed, building such a model for clustering gives strong tools to evaluate generalization capabilities (with famous criteria such as Akaike Information Criterion AIC, Bayesian Information Criterion BIC summarized by [START_REF] Duda | Pattern classification[END_REF] or even Integrated Completed Likelihood ICL [START_REF] Biernacki | Assessing a mixture model for clustering with the integrated completed likelihood[END_REF] etc.).

Discriminative methods for clustering were initially inherited from supervised classification these last two or three decades. They are also extended to unsupervised classification (a. k. a. clustering). In clustering these discriminative approaches would not build a model that would fit the data but would rather separate the ouput classes or groups from each other (e.g. in a one-vs-one or one-vs-rest manner) focusing on the boundaries of the groups rather than on the groups themselves. Spectral Clustering [Von [START_REF] Luxburg | A tutorial on spectral clustering[END_REF] or DIFFRAC [Bach and Harchaoui, 2008] are two examples of such techniques.

k-Means solves an (Optimal) Transport Problem

We take a close look at the k-Means loss for data (x i ) i=1...N into K groups:

min σ,µ 1 N N ∑ i=1 x i -µ σ(i) 2 2
(24) that we optimize over the centroids µ = (µ k ) k=1...K and the assignments function

σ from 1, N ⊂ N to 1, K ⊂ N.
We studied k-Means which is probably both the oldest and most famous clustering algorithm [START_REF] Jain | Data clustering: 50 years beyond K-means[END_REF] and we realized that it tries to efficiently solve an optimal transport problem. Put differently, the global minimum of the k-Means loss satisfies the optimal transport problem of choosing a limited K number centroids (µ k ) k=1,...,K such that the data empirical distribution p = 1 N ∑ N i=1 δ P x i on the one hand and the centroids weighted distribution q = ∑ K k=1 π k × δ µ k on the other hand would be the closest possible in the 2-Wasserstein sense associated to the squared euclidean distance, although usually, we use the 1-Wasserstein distance associated with the plain and simple euclidean distance:

W c 2 (p, q) = min γ∈Γ(p,q) E (x,m)∼γ x -m 2 2 (25)
and for discrete distributions:

W c 2 (p, q) = 1 N K ∑ k=1 π k × E x∼p x -µ k 2 2
(26) which is the k-Means loss and of course we have a link between the proportions π and the assignments σ: π k = #{σ(i)=k| i∈ 1,N } N as for discrete distributions, optimal transport plans are degenerated [START_REF] Peyré | Computational optimal transport[END_REF]. This interesting link between clustering and optimal transport encouraged us to investigate further between these two scientific literatures: clustering and optimal transport. Generalizing this observation to more sophisticated distributions thanks to Generative Adversarial Networks could lead to having a fatter support distribution than just Diracs:

min σ,(θ p k ) k=1...K W( K ∑ k=1 π k × p k , 1 N N ∑ i=1 ×δ x i ) (27)
where θ p k parametrized the kth cluster generator distribution p k (that was previously reduced to a Dirac distribution located on the kth centroid µ k ). In practice, we can have latent variables z ∼ N (0 d , I d ) and we define in y k ∼ p k by y = G k (z) as we will develop later. That initial idea we had is encouraging in a sense that it suggests an interesting mixing between optimal transport, neural networks and model-based clustering. We think that there is an opportunity here to mention that k-Means is not a regular special case of Gaussian Mixture Models through Expectation Maximization (EM-GMM) although this false idea is widespread. It is true that when neg-exponentiated, the k-Means looks like the negative likelihood of a Gaussian mixture fitted on the data with identity (or equally proportional to) covariance matrices and equal proportions but there is major difference: memberships probabilities in EM-GMM are not degenerated but they are in k-Means. k-Means can still be seen as an example of the Expectation-Maximization technique but these distinctions have been clearly made by [START_REF] Celeux | A classification EM algorithm for clustering and two stochastic versions[END_REF]. The link between k-Means and optimal transport is stronger than the one between k-Means and EM-GMM and also more fruitful in terms of open research.

Unsupervised Feature Importance

As suggested earlier when we analyzed the consequences of the impossibility theorem by [START_REF] Kleinberg | An impossibility theorem for clustering[END_REF], metric invariance is an interesting subject. To handle this difficulty related to the choice of the metric, one can work with a set of metrics because we know that the metric defined by the upper bound evaluation over a set of ditances is also a distance. To the best of our knowledge, this angle has not been tackled by the research community to study unsupervised feature importance extraction.

In supervised learning, [START_REF] Breiman | Random forests[END_REF] proposed routines based on permutation and mean decrease in impurity but much less work has been done in the unsupervised context. This is probably due to the fact two interleaved problems remain: metric learning, feature selection which makes our task ill-posed. We actually suggest that it is worth trying to improve regular and generic euclidean approaches.

Uncertain Predictions

In this chapter devoted to uncertainty in supervised learning, we borrow scientific items from Bayesian and frequentist scientific communities. Indeed, we believe that mixing both scientific cultures is beneficial in general and for uncertainty estimation in particular: probabilistic interpretations provided by Bayesian formulations and function expressivity fitted to large scale data provided by the large frequentist fauna of algorithms which is of course not limited to deep neural networks.

Bayesian and Frequentist scientists in Statistical Learning

There is a so-called rivalry in automatic statistical learning between Bayesians on one side and frequentists on the other side that is quiet disturbing when seen from a young scientific point of view. Indeed, we can laughably notice that Bayesian scientists are allowed to use histograms of frequencies and frequentist scientists are allowed to use the Bayes rule. This kind of debate is often sterile but it is probably a characteristic of still young non-unified sciences with varying names accross trends. In this section, we briefly describe the specific issues taken into account by these two commonly separated communities.

The Bayesian framework [START_REF] Barber | Bayesian Reasoning and Machine Learning[END_REF] is characterized by its use of parametrized probability laws which allows to benefit from interpretation ease (including uncertainty) when predicting information of producing description about data. The preferred statistical tool is usually the so-called bayes rule, hence the name of this scientific community. The choice of the modelled distributions carries interpretation and knowledge that is elegantly injected into the trained systems.

In contrast, within the frequentist approach [START_REF] Bishop | Pattern recognition[END_REF], we assume that uncertainty is inherenty present due to the randomness coming from repeatable experiments producing empirical observations. Hence, many machine learning problems tackled with a frequentist point of view is a statistical estimation problem based on observed data. If one could accept caricatures, then we would say that Bayesian statistical learning is a principled probabilistic and thus interpretable framework (hence their appealing reputation but at the cost of often wrong model assumptions) which offers extraordinary research avenues such as model selection without extra data and meaningful probability interpretations whereas frequentist statistical learning produce good-results-oriented blackboxes [START_REF] Neal | Bayesian Learning for Neural Networks[END_REF] with impressive recognition rates resculpting our modern world.

Recently, in a Ph.D. thesis [START_REF]Uncertainty in Deep Learning[END_REF], dared the idea of taking advantage of both worlds in a Bayesian deep learning approach (although the merit was certainly to revisit such a counter-intuitive approach that in fact dates back to at least in the 1990s by [START_REF] Neal | Bayesian Learning for Neural Networks[END_REF] or even by [START_REF] Bishop | Mixture density networks[END_REF]). Today in 2020, there is a still a controversial debate on this subject: Bayesian machine learning injects some knowledge through a distribution prior (not always a Gaussian prior even if we can recognize this is the best studied distribution and the most frequently used) for inputs and outputs of statistical predictors. This guiding of the machine learning at both optimization and prediction steps assumes some knowledge about input and output data but most of the time that knowledge is not existing hence the debate. Frequentist neural networks with all their parameters and Bayesian statistics with all their parametrized distributions might seem uneasy to coexist in the same unique machine learning method.

Sources of Uncertainty

Beyond frequentist and Bayesian considerations, according to recent research work [START_REF]Uncertainty in Deep Learning[END_REF] (we recommend the reader to read this Ph.D. thesis for details and comprehensive bibliography about prediction with uncertainty), uncertainty can be broken down into several facets: [START_REF] Srivastava | Improving neural networks with dropout[END_REF]. Overfitting is bad extrapolation in essence. Thus, neural networks are particularly prone to extrapolation issues like extrapolation uncertainty.

E X T R A P O L AT I O N U N C E R TA I N T Y
A L E AT O R I C U N C E R TA I N T Y Some noise could have been introduced in training data (some wrong coordinates, some wrong labels etc.). Well-studied statistical machineries like linear and kernel-based support vector machines [START_REF] Andrew | An Introduction to Support Vector Machines and Other Kernelbased Learning Methods[END_REF] proposed the notion of margin to cope with some part of that uncertainty but it seems that linking that margin to probability estimation is still an open research problem even with a logistic regression loss instead of a hyperplane in practice;

E P I S T E M I C U N C E R TA I N T Y
The initial machine learning problem might be illdefined: many solutions can solve the problem, so we do not know which one to choose objectively. For example, the butterfly effect (in layman's terms) is a source of epistemic uncertainty against meteorological forecasts meaning that uncertainty is inherently related to the studied laws of physics as a scientific standing point (see for example [START_REF] Epstein | Stochastic dynamic prediction[END_REF]). Thus, the problem modeling could be unsufficient which introduces randomness.

There exists more precise analysis for describing these interleaved sources of uncertainty and readers may find more exhaustive research in the work of [START_REF] Kennedy | Bayesian calibration of computer models[END_REF]. But all these phenomena boil down to how systems should and could handle the possible lack of confidence tainting automatic predictions.

WA S S E R S T E I N C L U S T E R I N G

Man Gave Names to All the Animals Bob Dylan in Slow Train Coming, 1979 Joint work with Charles Bouveyron, Pierre-Alexandre Mattei and Andrés Almansa

A B S T R A C T
Clustering is partitionning the data into groups. Deep approaches for clustering are promising for extending the success of neural networks beyond the limits of supervised classification. In this chapter, clustering is tackled at the crossroad of several literatures: auto-encoders, generative adversarial networks (GANs) for optimal transport, statistical mixture models. Two criteria are studied: (i) the first generatively minimizes the Wasserstein distance between data and clusterseparated generated data inspired by the GANs success and (ii) the second discriminatively maximizes over all partitions the Wasserstein distances between the associated groups. These two mechanisms are compatible with model selection according to a Wasserstein criterion measured on held-out validation data. Competitive results are achieved on benchmark datasets such as images, sparse and dense data, with the benefits of selecting the number of groups which promises interesting further research.

K E Y W O R D S
Clustering, Neural Networks, Wasserstein Generative Adversarial Networks, Mixture Model, Deep Generative Models, Optimal Transport, Discriminative Clustering, Generative Clustering

I N T R O D U C T I O N

In a new environment, the first thing a human being (even as a child) does is mentally grouping the elements of the surroundings and putting names on those groups. The ambition of this work is to tackle the problem of clustering within the deep learning framework. While the vast majority of the model-based clustering approaches focused on the Kullback-Leibler divergence (closely related to maximum likelihood up to an additional constant term [START_REF] Bishop | Pattern recognition[END_REF]), we investigate in this chapter an attempt to operate with the Wasserstein distance in lieu of Kullback-Leibler divergence for clustering with mixture distributions.

In the first part of this chapter, we take some generative ingredients coming from the machine learning literature to build our clustering technique as chronologically, we got inspired by Generative Adversarial Networks (GAN) remarkable results [START_REF] Goodfellow | NIPS 2016 tutorial: Generative adversarial networks[END_REF]. In the second part of this chapter, on the contrary, we choose more discriminative ingredients while still using some of the statistical and algorithmical tools empowered by the experience gained previously. Recently and independently from us, we have seen that [START_REF] Mukherjee | ClusterGAN: Latent Space Clustering in Generative Adversarial Networks[END_REF] successfully built a better generative clustering work than us but our approach is slightly different: Jensen-Shanon divergence from the original GAN [START_REF] Goodfellow | Generative adversarial nets[END_REF] whereas we used the Wasserstein distance in our case thanks to the work of [START_REF] Arjovsky | Wasserstein generative adversarial networks[END_REF].

In the discriminative clustering context, [START_REF] Bouveyron | On the estimation of the latent discriminative subspace in the Fisher-EM algorithm[END_REF] revisited an ancient approach from [START_REF] Fisher | The use of multiple measurements in taxonomic problems[END_REF]: separating clusters with a maximum Kullback-Leibler divergence inter-cluster hence the discriminative point of view. Likewise and with the same spirit, we propose to maximally separate data into clusters which is exactly adopting a discriminative angle to clustering instead of generating closest possible clusters to data which corresponds to a generative angle. Here, we explore new ways to do clustering with inter-cluster Wasserstein distances maximization moving away from the conventional literature by borrowing ideas from generative model-based approaches and discriminative methods.

Related Work

Authors such as [START_REF] Flamary | Wasserstein discriminant analysis[END_REF] actually did go down that appealing research road for Wasserstein-distance-based discriminative clustering with simultaneous linear dimensionality reduction. Indeed, [START_REF] Flamary | Wasserstein discriminant analysis[END_REF] adapted the old latent discriminant analysis of [START_REF] Fisher | The use of multiple measurements in taxonomic problems[END_REF] but with the Wasserstein distance instead of a sum of distances approach that is closer to a Kullback-Leibler mindset while still being both discriminative and model-based. The choice of divergence (even beyond Kullback-Leibler and Wasserstein) and the choice between generative and discriminative approaches are scientifically intruiguing and produce different algorithms in practice as far as clustering is concerned.

The ambition of this work is to tackle clustering within the deep learning framework because of its large success in supervised learning that we want to inherit in unsupervised classification (a. k. a. clustering). More precisely, our motivation does not come from the popularity of the so-called deep learning approaches but is rather nurtured by the functional expressivity that neural networks and the fact that the associated scientific community is profuse in free high-quality toolboxes, which is, of course, impossible to distangle from the deep learning tremendous popularity.

To achieve these goals, we take advantage of the recent Wasserstein Generative Adversarial Networks research to estimate these aforementioned generative and discriminative criteria. The smoothness implied by the deep learning optimization gradual procedures made us think that soft memberships probabilities should be prefered over discrete categorization output. Handling the Wasserstein distance meant dealing with transport plans which are uneasy objects intuitively in terms of software programming. Thanks to the Kantorovich-Rubinstein formulation, we end up with uncoupled lossed for the distributions which is easier but at the price of an adversarial min-max optimization for the generative case which is notoriously difficult to monitor in practice (although much research attention simplified it). For the discriminative case, maximizing the Wasserstein distance estimated by the Kantorovich-Runbinstein maximization duality is especially suitable for optimization stability. Clearly, we do benefit from algorithmic tools from the adversarial neural networks for the generative part of this work but we still do benefit from this tools for the discriminative easier algorithm construction, especially for critics (or potential) lipchitzian functions. Indeed, for the discriminative algorithm, we surprisingly avoid undergoing the difficulty and instability of a min-max (adversarial) optimization with just an overall maximization instead.

At the very start of this thesis in 2016, in our preliminary experiments, a Gaussian mixture model trained with Expectation-Maximization [START_REF] Dempster | Maximum likelihood from incomplete data via the EM algorithm[END_REF] on codes coming from a vanilla MLP1 auto-encoder without convolutions is able to reach above 80% of unsupervised clustering accuracy on the famous digits MNIST images dataset2 directly on raw pixels. Encouraged by this surprisingly good result, we pursue our efforts towards an attempt to take advantage of that empirical fact in a sound framework.

While supervised classification has been a long-standing problem for many decades until recently thanks to computational hardware dramatic improvements and a considerable research effort towards statistical tools and algorithms, unsupervised classification (a. k. a. clustering) is still a difficult area but taking advantage of the supervised findings proved efficient in terms of research. In fact, although we can only confirm that clustering is an ill-posed problem as explained earlier due to the Kleinberg's impossibility theorem [START_REF] Kleinberg | An impossibility theorem for clustering[END_REF], we still find it desirable with the same issues that the supervised classification research had to solve in easier settings because supervised classification has a clearer objective. Beyond ill-posed problems, the way we tackle model-based clustering in this work has at least two issues, namely model selection and mixture parametrization that we brielfy describe here before presenting two techniques, one generative which we improved thanks into a discriminative one.

Model Selection

In supervised techniques, after training different models with some different hyper-parameters, we can find th best set of hyper-parameters thanks to (cross-) validation: measuring the different associated accuracy scores on a labeled heldout dataset. Unfortunately, this is not a solution in our unsupervised clustering task where no labels exist (never: neither at training, nor validation and of course not at testing stages). Nevertheless, we have two ways to circumvent this problem [START_REF] Bouveyron | Model-Based Clustering and Classification for Data Science: With Applications in R[END_REF]:

W I T H A H E L D -O U T A N D U N L A B E L E D D ATA S E T our model-based techniques
compare distributions within a statistically meaningful quantity (divergence or distance minimization or maximization): one representing the data and another one representing a fitted model. A natural way to select some hyper-parameters sets among several trained models is to measure the same training quantity and same model but with different data corresponding to the validation dataset. In this work, we develop two techniques: DiWaC (which corresponds to a discriminative approach) and GeWaC (which corresponds to a generative one) and these objectives measure a model-data fitting Wasserstein score and allow us to measure a model-data fitting score on unseen and held-out data. In supervised classification, the validation accuracy in measured from the comparison between pre-annotated labels and predictions from a model optimized with training data. In our unsupervised clustering case, we allow ourselves to loosely adapt the held-out validation score expression because we simply measure some Wasserstein distances between model distributions and validation empirical data distributions (without labels otherwise this would not be realistic) to check under-and over-fitting phenomena;

W I T H T O U T A N Y O T H E R D ATA S E T T H A N T H E I N I T I A L T R A I N I N G S E T
a vast literature exists with Kullback-Leibler divergences and Likelihood Maximization for selecting a good model among many [START_REF] Schwarz | Estimating the dimension of a model[END_REF][START_REF] Biernacki | Assessing a mixture model for clustering with the integrated completed likelihood[END_REF] using the number of model parameters following an Occam's razor for Bayesian Machine Learning to compensate for extra data validation but to the best of our knowledge, we do not know nonlikelihood-based techniques (here we use Wasserstein distances instead).

We clearly chose the first approach with some held-out unlabeled data with our Wasserstein distances. Now we explore how we manipulate mixtures distributions thanks to the reparametrization trick.

Reparametrization Trick for a Mixture

In order to handle distributions and more precisely being able to differentiate our objective functions with respect to the parameters of these distributions for learning purposes, we adapt the Reparametrization Trick from [START_REF] Kingma | Variational dropout and the local reparameterization trick[END_REF] first used in a variational context3 . Indeed, we have chosen a probabilistic framework from which probability distribution parameters are estimated. Usually, in deep learning, first, the optimization process finds functional parameters for output prediction thanks to the minimization (or maximization) of an objective funcion and the notion of gradient gets easy under mild (sub)-differentiability conditions towards training optimization. In our case, this is different: we do not optimize functional parameters for an output prediction but we do find density parameters for a distribution which plays the role of the output prediction in a probabilistic fashion. In less than a decade, the variational (see the surveys accomplished by [START_REF] Kingma | An introduction to variational autoencoders[END_REF]) and optimal transport (see the book of [START_REF] Peyré | Computational optimal transport[END_REF]) literatures in machine learning got much attention for that estimation scenario with probability distribution parameters. Instead of using some sophisticated mathematical tools for differentiating some objective over some distribution parameters, many researchers used the Reparametrization Trick thanks to [START_REF] Kingma | Variational dropout and the local reparameterization trick[END_REF]: transforming a known pseudo-random noise as an estimator of the theoreretical random variable for optimization reasons which makes differentiation easier. [START_REF] Graves | Stochastic backpropagation through mixture density distributions[END_REF] gives a more comprehensive overview of the problem.

Indeed, [START_REF] Kingma | Variational dropout and the local reparameterization trick[END_REF] allow to directly specify a prior distribution over the code space of a variational auto-encoder (VAE). Inference is done using a stochastic gradient variational bayesian (SGVB) method, based on a reparametrization of the variational lower bound. In this work, we will revisit and adapt this technique called the Reparametrization Trick for our distributions settings. Deep generative models for clustering may be built using a mixture model as prior distribution. This approach was recently explored by [START_REF] Dilokthanakul | Deep unsupervised clustering with gaussian mixture variational autoencoders[END_REF] and [START_REF] Jiang | Variational deep embedding: A generative approach to clustering[END_REF] who used a Gaussian mixture prior.

In practice, the Reparametrization Trick consists in manipulating a random variable z coming from a parametrized distribution (say Gaussian of mean a and covariance matrix B = CC ⊤ such that z ∼ N (a, B)) thanks to a default random generator (here ǫ ∼ N (0, I)) that is transformed through (z = a + Cǫ which simulates z ∼ N (a, B)) in order to get a differentiable version of the random variable we needed. This technique can used in many different ways according to the literature [START_REF] Graves | Stochastic backpropagation through mixture density distributions[END_REF][START_REF] Doersch | Tutorial on variational autoencoders[END_REF][START_REF] Blei | Variational inference: A review for statisticians[END_REF].

G E N E R AT I V E WA S S E R S T E I N C L U S T E R I N G

Contemporary to a frenetic rythm of papers about Generative Adversarial Networks (GAN) since their spectacular birth [START_REF] Goodfellow | Generative adversarial nets[END_REF], we admit we have been deeply influenced especially by the Wasserstein-based approaches [START_REF] Arjovsky | Wasserstein generative adversarial networks[END_REF] on the one hand and the French aura associated to optimal transport since the acclaimed work of [START_REF] Villani | Optimal transport: old and new[END_REF] on the other hand. Indeed, [START_REF] Arjovsky | Wasserstein generative adversarial networks[END_REF] proposed a praiseworthy attempt to outline this particularly fast-growing scientific landscape of GANs initiated by [START_REF] Goodfellow | NIPS 2016 tutorial: Generative adversarial networks[END_REF] by insisting on the mathematical quantities minimized between real and generated beyond the detective-forger metaphor describing the revisited min-max adversarial optimization method (the forger being the generator function and the detective being the discriminator or critic function).

At the same time, we wanted to accomplish some model-based clustering contributions thanks to the statistical legacy summarized by Bouveyron and Brunet-Saumard [2014b] and later in the commendable book of [START_REF] Bouveyron | Model-Based Clustering and Classification for Data Science: With Applications in R[END_REF]. As clustering and data imitation are both unsupervised tasks, we wondered since 2016, how was it possible to accomplish clustering with data imitation tools such as GANs. Indeed, we investigate in this part a generative approach by mimicking data that matches best real data with respect to the Wasserstein distance.

Deep Generative Models for Clustering

In the recent statistical learning literature, there is a significant trend towards better deep generative models (DGM) based on different inference procedures: (i) likelihood, (ii) GANs.

First, the goal of the likelihood-based approach is to minimize the Kullback-Leibler divergence between the original data distribution and the parametrized model data distribution (which is equivalent to maximizing the likelihood [START_REF] Duda | Pattern classification[END_REF]). A recent example of that kind is the work of [START_REF] Dinh | Density estimation using real nvp[END_REF] that had to introduce neural networks bijections in order to apply a revisited change of variable formula on the likelihood mainly because the Kullback-Leibler divergence requires same distribution supports and thus a bijection is almost mandatory to parse the entire data space (including the vast empty data zones exactly like Gaussian mixture where each Gaussian component support is the entire space).

Second, inverting the problem of accessing to the parametrized model data distribution is to sample generated data from it which what [START_REF] Goodfellow | Generative adversarial nets[END_REF] astutely proposed with GANs: they minimize a divergence between real and generated data over the generator function parameters implemented by a neural network thanks to a critic function that evaluates that divergence. Since then, the Wasserstein distance [START_REF] Gulrajani | Improved training of wasserstein gans[END_REF] through the Wasserstein GAN (WGAN) was also proposed among other divergences summarized by a survey done by [START_REF] Goodfellow | NIPS 2016 tutorial: Generative adversarial networks[END_REF] gives a glimpse of what is offered by that thriving on-going scientific literature.

We aim at clustering a dataset of N points x 1 , ..., x i , ..., x N samples of the random variable x living in a space X (say R D ) into K homogeneous groups. We suppose there exists a latent code space Z of a low dimension d (say d = 10 which is low compared to the original data dimensionality D: d ≪ D) such that there is a mapping D between Z to X connecting the random variable x in X and its latent counterpart z in Z. We also assume that z follows a mixture M of rather common distributions as components (say Gaussian).

At the heart of our Generative Wasserstein Clustering (GeWaC) model, there is an auto-encoder made of: (i) an encoder network E (parametrized by θ E ) and (ii) a decoder network D (parametrized by θ D ). That auto-encoder plays the role of a two-ways bridge between the data space and a code (or latent) space which more akin to clustering alleviating the curse of dimensionality thanks to its lower dimensionality (d instead of D).

Our model consists in saying that the data have been generated as follows. The clustering variable c c ∼ Cat(π)

corresponds to a categorical random variable among K ≥ 2 clusters with prior proportions defined in vector π (of K positive scalars which sums to 1). In this generative process, once the cluster k is chosen, one can generate a code in Z which implies a mixture marginal for z:

z|c = k ∼ g k (.; ξ k ) z ∼ K ∑ k=1 π k g k (.; ξ k ) (29)
for the probability distributions (g k ) K k=1 of each of the K components parametrized by ξ k . Ultimately a point in X is generated:

x|z ∼ δ D(z) (30)
To translate into statistics the assumption that the code data lie extremely close to a low-dimensional manifold, we should have written N (D(z), σ 2 I p ) instead of a Dirac δ D(z) located on a decoded data point D(z) and them further assume that σ → 0. In this context, the posterior probability needed to cluster x is given by

P(c = k|x) = E z∼p(z|x) [P(c = k|z)] (31)
Although any parametric density functions can be used for each mixture component g k , we restrict ourselves in this work to densities allowing the use of the reparameterization trick [Kingma andWelling, 2013, Kingma et al., 2015] which has been described above and has a central role in our adversarial optimization as we will see later.

In the following, we will illustrate our methodology using a mixture of Gaussians in the code space, i.e. we choose:

g k (.; ξ k ) = N (.; µ k , Σ k ) (32)
(a Gaussian distribution with mean µ k and covariance matrix Σ k ) where (π k ) k=1,...,K , (µ k ) k=1,...,K , and (Σ k ) k=1,...,K are the mixture parameters stored in θ M . Note that, beyond the Gaussian mixture prior that we consider here, our approach could be extended to any mixture of reparametrizable distributions: one might for example consider a mixture of von Mises as the prior distribution, in order to obtain interesting visualizations on a hyper-sphere, such as the ones of [START_REF] Davidson | Hyperspherical variational auto-encoders[END_REF]. Our graphical model displayed in Fig. 11 summarizes the chosen data generation modelling. In order to fit our generative model, we generate a random variable D(z) to match x in terms of Wasserstein distance in an adversarial fashion. We only have access to samples of x (through the dataset) and D(z) (through a pseudo-random generator for z and the decoder or generator D), which is a scenario where the WGANs proved successful as Fig. 12 summarizes the modeling proposed here.

The posterior probability p(z|x) in Eq. ( 31) is hard to compute because of the nonlinearity of the decoder. But, as in the work of [START_REF] Kingma | Auto-encoding variational bayes[END_REF], we can approximate it using an inference network q(z|x) built according to the encoder E . As emphasized by [START_REF] Kingma | Variational dropout and the local reparameterization trick[END_REF], minimizing the Kullback-Leibler divergence between the true posterior and q(z|x) leads to minimizing a penalized quadratic auto-encoder loss. Since σ → 0, the dominating term in this loss will precisely be the loss of a vanilla auto-encoder which is what we do in practice for the sake of simplicity. Eventually, we can compute an approximation of P(c = k|x) by simply replacing the true posterior by the approximation, which leads to the maximum-a-posteriori (MAP) rule in code space:

P(c = k|x) ≃ π k g k (E (x); ξ k ) ∑ K k ′ =1 π k ′ g k ′ (E (x); ξ k ) (33)
following the Bayes formula.

After several attempts, we came up with two fruitful strategies to make our system work: first, the Concatenation Trick from [START_REF] Dumoulin | Adversarially learned inference[END_REF] to ensure some consistency between code and data spaces and second, the Reparametrization Trick from [START_REF] Kingma | Variational dropout and the local reparameterization trick[END_REF] to handle parametrized mixture distributions described earlier.

Concatenation Trick

Our concatenation approach consists in operating a clustering that is consistent both in the embedding code and input data spaces as described in Fig. 12: • on the one hand, there is one space (called "code space") with real encoded data (random variable E (x)) next to a generated mixture distribution z ∼ M;

• on the other hand, there is an other space (called "input space") of real data at hand (random variable x) coexisting with generated decoded signal D(z).

This mechanism allows our GAN technique to bring real and generated data distributions together. We use a small dimensionality code space to benefit from the natural probabilistic interpretation of mixture (e. g. Gaussian) for clustering in which each code mixture component corresponds to a cluster. To guarantee a minimum level of consistency for our encoder/decoder system bridging between our two aforementioned spaces, we use a technique that we dub the "concatenation trick" proposed by [START_REF] Dumoulin | Adversarially learned inference[END_REF]. This is necessary to make sure the encoder E and decoder D functions are reciprocal mathematically almost everywhere on the data and code manifolds at hand (see [START_REF] Donahue | Adversarial feature learning[END_REF] for proof).

On the data side, if we only bring together the distributions of x and D(z), then this would not have been enough because clustering would take place for E (x) and z. On the code side, if we only bring together the distributions of E (x) and z, then this would not have been enough either because true data x would not be even concerned anymore which would be equivalent to disconnecting the codes from data. Concatenating data and codes solves these two problems. Fortunately, in other settings (data imitation and compression), the concatenation trick formulation of [START_REF] Dumoulin | Adversarially learned inference[END_REF] nicely fits ours, thanks to the idea of bringing together the distributions of

• x = a(x) = x ⊤ , E (x) ⊤ ⊤ ∼ p considered as real with x representing data; • ỹ = b(z) = D(z) ⊤ , z ⊤ ⊤ ∼ q considered as generated with z sampled over a parametrized mixture q = M(π, µ, Σ) = ∑ K k=1 π k N (µ k , Σ k )
. by minimizing the Wasserstein distance between p and q over all parameters (mixture, encoder and decoder). This way, we elegantly get the almost everywhere reciprocity between the encoder E and decoder D on the data manifold adapting the work of [START_REF] Dumoulin | Adversarially learned inference[END_REF] and also [START_REF] Donahue | Adversarial feature learning[END_REF].

The nature of the chosen GAN (Wasserstein GAN or other) is not crucial here as [START_REF] Arjovsky | Wasserstein generative adversarial networks[END_REF] explained that the original GAN [START_REF] Goodfellow | Generative adversarial nets[END_REF] uses the Jensen-Shannon divergence and WGAN uses the Wasserstein distance (which is a fortiori a divergence) but [START_REF] Dumoulin | Adversarially learned inference[END_REF] empirically show the improvement attributed to the concatenation trick vs. without it in terms of image rendering. Interestingly, concatenation becomes a key element in our case to make the whole system work because of our clustering goal (that [START_REF] Dumoulin | Adversarially learned inference[END_REF] do not have) to maintain data-code consistency through the encoder/decoder reciprocity. Intuitively, if the concatenated variables in the previous bulletted list have similar distributions with respect to the chosen GAN-specific divergence at hand, then the marginals are close too. Thus, our use of the concatenation trick gets handy because it ensures that our pair of encoder/decoder functions (E , D) still behaves in a reciprocal fashion in spite of the mixed GAN and mixture framework in a surprisingly stable manner in terms of optimization.

In fact, we first proposed a model without concatenation which meant the adversarial WGAN optimization did not involved the encoder E as only the distribution of x and D(z) were being put together: this was not principled especially because the clustering decision rule given in Eq. ( 33) involves E which forced us to retrain a final post-processing with an adhoc autoencoder loss with all parameters fixed except the ones of E for what seems to be an update with respect to the decoder D. The contribution of [START_REF] Dumoulin | Adversarially learned inference[END_REF] made all these considerations disappear in an easy yet principled fashion. Now, we minimize the Wasserstein distance between these two augmented data associated distributions thanks to the Kantorovich-Rubinstein duality with a WGAN [START_REF] Arjovsky | Wasserstein generative adversarial networks[END_REF][START_REF] Salimans | Improved techniques for training gans[END_REF][START_REF] Gulrajani | Improved training of wasserstein gans[END_REF][START_REF] Miyato | Spectral normalization for generative adversarial networks[END_REF] :

max ∇C ≤1 E x∼p C a(x) -E z∼q C b(z) (34)
with C called the critic and implemented in practice by a neural network of parameters θ C and constrained as 1-Lipschitzian by a chosen method among weight-clipping [START_REF] Arjovsky | Wasserstein generative adversarial networks[END_REF], augmented Lagrangian [START_REF] Gulrajani | Improved training of wasserstein gans[END_REF] or with more stability from online power iteration [START_REF] Miyato | Spectral normalization for generative adversarial networks[END_REF].

Assuming that the codes z come from a mixture of Gaussians with full covariance matrices Σ k ∈ R d×d (d ≪ D), for each k among the K components the corresponding variable z k follows:

y k = b(z k ) = D(S k × e + µ k ) ⊤ , (S k × e + µ k ) ⊤ ⊤ ( 35
)
where e ∼ N (0 d , I d ) and S k is a Cholesky decomposition-inspired representation (with non-zeros only in the lower-triangular part and strictly positive diagonal entries guaranteed by exponentials first and then affine-transformed sigmoids for better eigenvalues amplitude control) of the full covariance [START_REF] Kingma | Auto-encoding variational bayes[END_REF][START_REF] Kingma | Variational dropout and the local reparameterization trick[END_REF] for unconstrained optimization which is much easier and well studied in stochastic gradient settings.

Σ k = S k × S ⊤ k such that the transformed random variable S k × e + µ k behaves as if it comes from N (µ k , S k × S ⊤ k ) = N (µ k , Σ ⊤ k ) in the spirit of the Reparametrization Trick [Kingma
All the equations above meet in:

L(θ E , θ D , θ M , θ C ) = E x∼p C x ⊤ , E (x) ⊤ ⊤ (36) - K ∑ k=1 π k × E e∼N (0,I d ) [C ([D(S k × e +µ k ) ⊤ , (S k × e + µ k ) ⊤ ⊤
and thus we optimize:

min θ E ,θ D ,θ M max θ C L(θ E , θ D , θ M , θ C ) (37) Indeed, max θ C L(θ E , θ D , θ M , θ C ) corresponds to the Wasserstein distance between real data x ⊤ , E (x) ⊤ ⊤ and generated data D(S k × e + µ k ) ⊤ , (S k × e + µ k ) ⊤ ⊤
by the mixture M and decoder D.

Originally Wasserstein GAN uses a simple fixed distribution (Gaussian or uniform) for the random noise generator that is transformed by a neural network called the generator to fit the data distribution in the Wasserstein sense. Here, we use a tunable mixture distribution instead for clustering purposes. In a regular Wasserstein GAN, the fixed distribution has no parameter taking part in the data generation mechanism. For that reason among other vocabulary reasons, in our GeWaC algorithm, we call generator the union of the parametrized mixture distribution M associated with the decoder neural network D (that brings the noise generated from the mixture into fake data).

The mechanism that we dub concatenation trick proposed by [START_REF] Dumoulin | Adversarially learned inference[END_REF] naturally enforces the encoder and the decoder being reciprocal which can be proven in a way that is very close to what [START_REF] Donahue | Adversarial feature learning[END_REF] did for interested readers.

Algorithm

Our GeWaC algorithm can be decomposed in three successive steps:

Auto-Encoder Initialization

We train a classic auto-encoder (E , D) (for an encoder E and decoder D):

(θ 0 E , θ 0 D ) = arg min θ E ,θ D E x∼p D • E (x) -x 2 2 (38)

EM-based Gaussian mixture initialization

We fit a Gaussian mixture model

M SoftMax(α), µ, (S k × S ⊤ k ) k=1,...,K (39)
parametrized by θ M = (α k , µ k , S k ) k=1,...,K with the Expectation-Maximization algorithm [START_REF] Dempster | Maximum likelihood from incomplete data via the EM algorithm[END_REF] on the encoded data:

θ 0 M = arg max θ M E x∼p log K ∑ k=1 π k × N (µ k , S k × S ⊤ k ) E 0 (x) (40)
where the covariance matrices are parametrized by Σ k = S k × S ⊤ k to save tedious symmetry and eigenvalues signs constraints and the proportions are parametrized by π = SoftMax(α) for easily imposing 1-sum positive constraints on the proportions.

Critic Initialization

The critic function C role is to estimate the Wasserstein distance in order to get good gradient estimation for the rest of the parameters. The previous steps intialized a generator process that we should evaluate first before taking the gradient from it for the other parameters:

θC = arg max θ C L(θ 0 E , θ 0 D , θ 0 M , θ C ) (41)
from Eq. ( 37) which is optimized thanks to the algorithm 1 calling the algorithm 2 but without 3.

Clustered Data Generation

The previous steps made this final step well initialized to optimize all the parameters thanks to the algorithms 1 calling both 2 and 3:

( θE , θD , θM ) = arg min θ E ,θ D ,θ M max θ C L(θ E , θ D , θ M , θ C ) (42)
Our GeWaC technique is "end-to-end trainable". Steps 1, 2 and 3 are just reasonable initializations for step 4. Finally, we use the MAP rule in the code space (Bayes formula) in order to finally cluster the data points which makes our simple approach after training time particularly fit for clustering. For training this generative clustering, the optimization is orchestrated by algorithm 1 that calls Wasserstein distance estimation updates in iterations involving only critics neural networks parameters which is described by algorithm 2 and minimizes the Wasserstein distance estimates with respect to all parameters in algorithm 3 excluding the ones of the critics.

Algorithm 1 Optimization algorithm

1: while θ E , θ D and θ M have not converged do 2: Sample a mini-batch of size B × K from the dataset x i,k i = 1, . . . , B k = 1, . . . , K 3:
Compute the critic evaluation on the mini-batch of points and codes concatenation

a i,k ← C([x ⊤ i,k , E (x i,k ) ⊤ ] ⊤ ) i = 1, . . . , B k = 1, . . . , K 4
b i,k ← C([D(S k × e + µ k ) ⊤ , (S k × e + µ k ) ⊤ ] ⊤ ) i = 1, . . . , B k = 1, . . . , K 4: Compute w ← 1 B × K B ∑ i=1 K ∑ k=1 a i,k - K ∑ k=1 π k 1 B B ∑ i=1 b i,k
5: Perform a gradient ascent step with w over θ C

Model Selection for GeWaC

Once we get our trained data generator, we can measure the Wasserstein distance between generated data and some held-out validation data distributions (that we can sample from) to check under/over-fitting and ultimately choose the number of classes or the architecture of neurons and layers. There is one subtlety though: we must measure non-augmented data Wassertein distance between non-augmented generated data and decoded mixture noise with recomputed proportions from validation memberships probabilities means. Otherwise, the encoder-decoder part of our systems will not be fairly compared: by removing the dimensionality augmentation, we make possible the selection of the coding space dimensionality for example. Thus, algorithmically, transforming our training procedure into a validation one consists in keeping algorithms 1 and 2 and leaving 3 and the augmentation parts. We could consider our generative attempt as a natural extension of what the Expectation-Maximization for the Gaussian mixture model algorithm does with the Kullback-Leibler divergence but with stochastic gradient descent for the Wasserstein Distance. One difficulty appears though from mixture distributions being problematic because of the partial discreteness of the parameter Algorithm 3 Wasserstein Minimization Step 1: Free encoder, decoder and mixture gradients accumulators 2: Sample some B × K Gaussian noise codes from a pseudorandom generator

e i,k ∼ N (0, I d ) i = 1, . . . , B k = 1, . . . , K 3:
Compute the critic evaluation on the mini-batch of decoded noise and its original version concatenation

b i,k ← C([D(S k × e + µ k ) ⊤ , (S k × e + µ k ) ⊤ ] ⊤ ) i = 1, . . . , B k = 1, . . . , K 4: Compute w ← 1 B × K B ∑ i=1 K ∑ k=1 a i,k - K ∑ k=1 π k 1 B B ∑ i=1 b i,k 5: Perform a gradient descent step with w over (θ E , θ D , θ M )
space [START_REF] Graves | Stochastic backpropagation through mixture density distributions[END_REF]. Specifically, mixture weights are particularly difficult to optimize. Consequently, [START_REF] Dilokthanakul | Deep unsupervised clustering with gaussian mixture variational autoencoders[END_REF] assume that these weights are known beforehand, which is not always the case in real-world problems. [START_REF] Jiang | Variational deep embedding: A generative approach to clustering[END_REF] propose a method for optimizing these weights but does not provide empirical evidence on imbalanced data to support this scheme. On our side, we did manage to get reasonable results but without completely being fair about proportions: if we trust our careful initialization procedures and let proportions have very low learning rate, then in practice, it is as if we freeze these proportions to a constant vector value although we did not realize it at first while being deceived by apparently good results except when initial proportions are not valid.

Collapsing Effects

Since its introduction of GANs, [START_REF] Goodfellow | Generative adversarial nets[END_REF] warned the reader about what they called the "Helvetica scenario" in which their generator is trained too often compared to the not-enough-updated discriminator (a word that is translated by critic since Wassertein GANs [START_REF] Arjovsky | Wasserstein generative adversarial networks[END_REF]). Indeed, a generator can intuitively be good at generating data from a specific region of space without being able to generalize to other space zones of data as the discriminator is fooled when comparing good localized generated data compared to real data. This is a kind of a spatial unsupervised over-fitting that is commonly defined in supervised learning. We end up with a GAN that gets unable to generate the same variety of data as real data hence the collapsing effect phenomenon name.

In our GeWaC technique, we observed a similar problem on the mixture side. Indeed, because of the richness of potential functions expressed by neural networks, only one obviously non-clustering mode of a mixture is enough to go through the decoder and parse the whole data manifold space. We decided to give up this technique for that reason to prefer a discriminative approach that does not need to produce data. In fact, it would be interesting to adapt what [START_REF] Mukherjee | ClusterGAN: Latent Space Clustering in Generative Adversarial Networks[END_REF] recently did but for the Wasserstein distance in order to see if their one-hot canonical latent augmentation clustering encoding circumevents our generative problems.

D I S C R I M I N AT I V E WA S S E R S T E I N C L U S T E R I N G

In this section, we do clustering with neural networks thanks to a discriminative (instead of generative) objective optimization: clusters distributions form mixture components like [START_REF] Fraley | Model-based clustering, discriminant analysis, and density estimation[END_REF], [START_REF] Bouveyron | Model-Based Clustering and Classification for Data Science: With Applications in R[END_REF] did in the past. The discriminative aspect comes from the fact that we are trying to optimally separate the data clusters in terms of Wasserstein distance (benefiting from the recent rising of deep learning scientific techniques since [START_REF] Arjovsky | Wasserstein generative adversarial networks[END_REF], [START_REF] Gulrajani | Improved training of wasserstein gans[END_REF], [START_REF] Miyato | Spectral normalization for generative adversarial networks[END_REF] but also [START_REF] Genevay | Entropy-regularized optimal transport for machine learning[END_REF]) in a one-versusrest fashion.

In a discriminative clustering, there is a fundamental limitation of Kullback-Leibler divergence techniques: distributions must share same (infinite) support. Indeed, if density supports are not the same, the Kullback-Leibler divergence is not defined (the logarithm of a zero probability being -∞) so in these circumstances we can (must) use infinite support densities such as Gaussians or related to artificially separate clusters that have mathematically same support which is not natural: How come separated clusters share same model density supports? The unsatisfactory answer consists in having low density separation zones between them. From that perspective, Wasserstein distances are better because they are well-defined for non-equal density supports which constitutes its main advantage thanks to its geometric properties. Clusters can now mathematically have model densities without any overlap.

We propose an algorithm to perform unsupervised classification (a.k.a. clustering) within this framework that we call "DiWaC" for Discriminative Wasserstein Clustering. Using the idea of separating clusters as a discriminative objective function is not new (see for example Spectral Clustering [Von [START_REF] Luxburg | A tutorial on spectral clustering[END_REF] or DIFFRAC [Bach and Harchaoui, 2008]) but the fact that we handle distributions allows to enable both out-of-sample clustering and model selection at the same time which makes our technique appealing even in large scale settings.

This work, to the best of our knowledge, is the first that maximizes Wasserstein distances between clusters in a discriminative manner. The advantage of maximizing in our case is that the Kantorovich-Rubinstein formulation makes it an overall maximization which is good news in terms of programming and convergence ease compared to usual (and painful to monitor and debug) min-max-type of optimization in Generative Adversarial Networks.

DiWaC is built at the crossroads of the auto-encoders, generative adversarial networks, optimal transport and statistical mixture models literatures. We recall that an auto-encoder is a neural network made up of two parts: on the one hand, an encoder that transforms the initial data into a smaller code space followed by a decoder that sends the codes back to the initial data space by trying to reconstruct them approximately in the sense of a quadratic loss (for a conventional auto-encoder), a Kullback-Leibler divergence (for a variational auto-encoder [START_REF] Kingma | Auto-encoding variational bayes[END_REF]) or the Wasserstein distance (for an adversarial auto-encoder [START_REF] Tolstikhin | Wasserstein auto-encoders[END_REF]).

By presenting clustering in a discriminative fashion, we end up by defining a good data partitionning in clusters as one whose components are as far apart as possible from each other. Hence, mathematically, we maximize the weighted sum of Wasserstein's distances between each cluster components and all others. Thanks to Kantorovich's formulation of Wasserstein's distances, the optimization of this criterion is a maximization (without minimization) on probabilities and critics functions (also called potential in the optimal transport literature). Thus, we benefit from the algorithmic tools coming from adversarial neural networks, especially for the critics Lipschitzian functions [START_REF] Miyato | Spectral normalization for generative adversarial networks[END_REF], without suffering from the hardship of an adversarial optimization (as there is only maximization and no minimization any more).

Wasserstein Distances between Clusters

Concretely, clustering is the task of gathering the data x in K ≥ 2, K ∈ N well separated groups but here we relax the hard notion of group into soft memberships through:

τ(x) = [P(c = 1|x), . . . P(c = k|x), . . . , P(c = K|x)] ⊤ (43)
The number K of groups can be found through model selection which is explained later (section 2.3.5) Through clustering we infer a function to describe hidden structure from unlabeled data. We aim at clustering a dataset of N samples x 1 , ..., x i , ..., x N of the random variable x (of distribution p) living in a space X (say X = R D ) in K groups (or clusters). Usually, in the literature in Machine Learning (a fortiori including Deep Learning), scientists minimize probability divergences, namely: Kullback-Leibler divergence which is equivalent to maximizing the likelihood in many cases (with a wide range of tools from logistic regression i. e. cross-entropy loss for classification to Expectation-Maximization for clustering), the Jensen-Shannon divergence in Generative Adversarial Networks [START_REF] Goodfellow | Generative adversarial nets[END_REF] as remarked and extended to Wasserstein distance (which is thus also a divergence with nicer and geometric properties) by [START_REF] Arjovsky | Wasserstein generative adversarial networks[END_REF]. The originality of the current work in clustering settings where the grouping are unknown lies on the maximization of divergences between groups instead of the usual minimization of divergences. This makes our contribution close to Discriminative Latent Models [Bouveyron and Brunet-Saumard, 2014a], Fisher Expectation Maximization [START_REF] Bouveyron | Simultaneous model-based clustering and visualization in the fisher discriminative subspace[END_REF] and Spectral Clustering [Von [START_REF] Luxburg | A tutorial on spectral clustering[END_REF] and we took some ideas from this vein of research (especially about how we handle proportions to normalize our objectives as we will see). The Wasserstein distance is also a divergence and thus has better properties than the Kullback-Leibler which constitute the main motivation of this work: mainly symmetry and geometric interpretations such as the triangle inequality and being defined even when the compared distributions do not have the same support.

With this research background in mind, we choose to still benefit from the Wasserstein generative adversarial networks (WGAN) without being adversarial. Indeed, we just see WGAN as an algorithmic tool to manipulate Wasserstein distances for large scale datasets thanks to its neural networks stochastic optimization. Thus a probabilistic configuration is needed and we choose the one of mixture models which actually stood the test of time in several research milestones done by Machine Learning pioneers like [START_REF] Dempster | Maximum likelihood from incomplete data via the EM algorithm[END_REF], [START_REF] Lloyd | Least squares quantization in PCM[END_REF], [START_REF] Blei | Latent dirichlet allocation[END_REF], [START_REF] Jain | Data clustering: 50 years beyond K-means[END_REF]. Let us model our data as coming from a mixture distribution of K unknown bu separated components p 1 , . . . , p k , . . . , p K weighted by proportions π

1 , • • • , π k , • • • , π K : x ∼ p = K ∑ k=1 π k × p k (44)
which uses the classical data generation model of mixture models [START_REF] Bouveyron | Model-Based Clustering and Classification for Data Science: With Applications in R[END_REF]]:

• Pick a cluster index k from the multinomial distribution of parameters π

• Pick a data point x from the data component distribution p k

Identifying the different components distributions p k s (and the proportions π k s) is recasting clustering in a probabilistic manner. Based on the presented mixture model, we try to maximize the Wasserstein distances between the different p k s themselves.

Meanwhile, thanks to the Bayes formula, we can get:

p k (x) = p(x) × τ k (x) π k (45)
and similarly, we define:

pk (x) = p(x) × 1 -τ k (x) 1 -π k (46)
the distribution of the remaining points without the kth group. A reasonable objective appears to consist in simultaneously maximizing all inter-cluster Wasserstein distances W(p k , pk ). Indeed, distributions p k s with highly overlapping support correspond to low inter-cluster Wasserstein distances because each component corresponds to a cluster. In contrast, well separated distributions p k s correspond to bigger inter-cluster Wasserstein distances W(p k , p k ′ ) and W(p k , pk ) (with k = k ′ ).

Unormalized and Normalized sum of inter-Cluster Wasserstein Distances

We previously established that inter-cluster Wasserstein distances maximization for building a training objective could be interesting for clustering. In this part of the current work, we present two attempts: one simple sum and one weighted sum to combine the inter-cluster Wasserstein distances.

Definition 1. Unormalized sum of inter-cluster Wasserstein distances

L u (p 1 , . . . , p k , . . . , p K ) = 1 K K ∑ k=1 W(p k , pk ) (47)
Definition 2. Normalized sum of Wasserstein inter-cluster distances L n (p 1 , . . . , p k , . . . , p K , π 1 , . . . , π k , . . . ,

π K ) = K ∑ k=1 π k × (1 -π k ) × W(p k , pk ) (48)
In these two objectives definitions, we maximize over p 1 , . . . , p k , . . . , p K and π 1 , . . . , π k , . . . , π K verifying:

• π 1 + π 2 = 1 with π ∈ R 2 + • p = π 1 × p 1 + π 2 × p 2

Sanity Checks

Thanks to a theoretical example with an obvious clustering outcome one could wish, we now try to understand the interests of choosing one of the two defined objectives. In R 2 , let us take a data distribution made of Gaussian and Dirac components: This almost Gaussian mixture case (corrupted by a down-weighted and farlocated Dirac distribution), it is still interesting to analyze what would happen when clustering with K = 2 groups making Eq. ( 47) ending up with the maximization of W(p 1 , p 2 ) (p 1 and p 2 being unknown).

p = 1 α + β + 1 × (αN (0 2 , σ × I 2 ) + βN (b, σ × I 2 ) + δ c ) (49
We consider two candidate clusterings to understand how compatible is the objective what we would expect as a good solution.

Example 1. Ideally, a satisfactory solution would be:

p good 1 = 1 α + 1 × (αN (0 2 , I 2 ) + δ c ) and p good 2 = N (b, I 2 ) ( 50
)
with proportions

π good 1 = α + 1 α + β + 1 ≃ 0.33 and π good 2 = β α + β + 1 ≃ 0.66 (51)
(switching the indices 1 and 2 does not break any generality). Indeed we believe that the oulier Dirac distribution δ c should be neglectable (thanks to the coefficients α and β being much bigger than one).

In other words, this would correspond to each cluster being associated with a single Gaussian beacuse the contribution of the Dirac in the mixture behaves like an outlier. Thus, we investigate here a rudimentary sanity check towards oultliers robustness.

Example 2. There is a bad clustering candidate that unfortunately gets a better score with respect to Eq. ( 47):

p bad 1 = 1 α + β × (αN (0 2 , I 2 ) + βN (b, I 2 )) and p bad 2 = δ c (52)
with proportions

π bad 1 = α + β α + β + 1 ≃ 0.99 and π bad 2 = 1 α + β + 1 ≃ 3.3 × 10 -5 (53)
In the Normalized Spectral Clustering literature (well explained by Von [START_REF] Luxburg | A tutorial on spectral clustering[END_REF]), a similar normalization is used for more robust clustering with respect to outliers.

Indeed for Eq. ( 47), the good and bad solutions give approximately (thanks to the chosen caricatural coefficients α, β, σ, m and M):

L u (p good 1 , p good 2 ) ≃ 9 (54) L u (p bad 1 , p bad 2 ) ≃ 100 (55)
which selects the bad candidate. This prooves the sensitivity of the unormalized objective in L u with respect to the outlier Dirac distribution in Eq. ( 49).

In objectives represented in L u and L n , the clusters are separated but L n avoids degenerate clustering candidates thanks to the π k × (1π k ) term. In our sanity check example, the previous bad clustering is smashed out by the normalization in L n , the first cluster occupies more than π 1 ≃ 99% of the data leaving the second singleton set with less than π 2 ≃ 0.003%) whereas the unormalized L u objective rewards a bad degenerate candidate solution. Indeed for the normalized objective L n gives approximately:

L n (p good 1 , p good 2 , π good 1 , π good 2 ) ≃ 3.92 (56)
L n (p bad 1 , p bad 2 , π bad 1 , π bad 2 ) ≃ 5.9 × 10 -3 (57)

selecting the expected good candidate clustering.

One-vs-One and One-vs-Rest Strategies

There is also a more principled way to look at our normalized objective Eq. ( 48) and the π k × (1π k ) term. In supervised classification (say logistic regression nicely explained by [START_REF] Hastie | The elements of statistical learning: data mining, inference and prediction[END_REF]), we have one-vs-one and one-vs-rest strategies that we use here. More concretely, in a one-vs-one strategy, it would be reasonable to consider the probability π k × (1π k ′ ) of choosing one point from cluster k and the other point from cluster k ′ in an independent fashion. Now we imagine two data generation models to add some theoretical justification to our objective Eq. ( 48). First, we can imagine a generation model:

(if k = k ′ the case is obvious since W(p k , p k ) = 0)
The observed Wasserstein distance is a random variable whose mean is:

L OvO n (p 1 , . . . , p k , . . . , p K , π 1 , . . . , π k , . . . , π K ) = K ∑ k=1 K ∑ k ′ =1 π k × π k ′ × W(p k , p k ′ ) (58)
which is the mean one-vs-one inter-Wasserstein distance between clusters selected by a multinomial distribution of parameter π.

The same reasoning can be done with a second and slightly different data generation model:

1. Sample two clusters indices (k, k ′ ) independently from the multinomial distribution of parameters π 2. Observe the Wasserstein distance W(p k , pk ) between the kth component and the rest

The observed Wasserstein distance is a random variable whose mean is now:

L OvR n (p 1 , . . . , p k , . . . , p K , π 1 , . . . , π k , . . . , π K ) = K ∑ k=1 π k × (1 -π k ) × W(p k , pk ) (59)
which is the mean one-vs-rest inter-Wasserstein distance between clusters selected by a multinomial distribution of parameter π. This one-vs-rest approach is preferred rather than the one-vs-one for combinatorial reasons as K grows. One can notice that the previous normalized Wasserstin sum definition is equal to that one-vs-one point of view definition L OvR n = L n but presented in a different fashion.

Decoupling the compared distributions in the Wasserstein distance computation is easier with the euclidean ℓ 2 distance. Indeed, the Kantorovich-Rubinstein [START_REF] Dudley | Real analysis and probability[END_REF] as it is already used in WGAN [START_REF] Arjovsky | Wasserstein generative adversarial networks[END_REF] tells us a decoupled re-definition:

W(µ, ν) KR = sup C∈Lip-1 E x∼µ C(x) -E y∼ν C(y) (60)
where Lip 1 is the 1-Lipschitz functions set. Thus, we can write that for three distributions q 1 , q 2 , q 3 and proportions κ 1 , κ 2 (with κ 1 + κ 2 = 1):

W(κ 1 × q 1 + κ 2 × q 2 , q 3 ) = sup C∈Lip-1 E x∼κ 1 ×q 1 +κ 2 ×q 2 C(x) -E y∼q 3 C(y) (61)
but

E x∼κ 1 ×q 1 +κ 2 ×q 2 C(x) = R D C(x) (κ 1 × q 1 (x) + κ 2 × q 2 (x)) dx (62) = κ 1 × R D C(x) × q 1 (x)dx + κ 2 × R D C(x) × q 2 (x)dx = κ 1 × E x∼q 1 C(x) + κ 2 × E x∼q 2 C(x) (63)
and

E y∼q 3 C(y) = (κ 1 + κ 2 ) × E y∼q 3 C(y) (64) thus W(κ 1 × q 1 + κ 2 × q 2 , q 3 ) = sup C∈Lip-1 E x∼κ 1 ×q 1 +κ 2 ×q 2 C(x) -E y∼q 3 C(y) (65) = sup C∈Lip-1 κ 1 × E x∼q 1 C(x) -E y∼q 3 C(y) +κ 2 × E x∼q 2 C(x) -E y∼q 3 C(y) ≤ κ 1 × sup C∈Lip-1 E x∼q 1 C(x) -E y∼q 3 C(y) +κ 2 × sup C∈Lip-1 E x∼q 2 C(x) -E y∼q 3 C(y)
because the maximum of a sum is lower or equal than the sum of each term maxima (which is also true with supremum instead of maximum) which gives

W(κ 1 × q 1 + κ 2 × q 2 , q 3 ) ≤ κ 1 × W(q 1 , q 3 ) + κ 2 × W(q 2 , q 3 ) (66)
and finally implies that

L OvR n (p 1 , . . . , p K , π 1 , . . . , π K ) ≤ L OvO n (p 1 , . . . , p K , π 1 , . . . , π K ) (67)
and thus maximizing the lower complexity one-vs-rest objective is maximizing a lower bound of the one-vs-one objective which is usual in Machine Learning.

In this section, we tried to provide some theoretical justification for our training objective and it is time to tackle the deep learning angle of it.

Estimating Wasserstein Distances with Deep Learning

In terms of neural network optimizatin, if we sum up everything we described, membership probabilities τ k (x) = P(c = k|x) are the unknown function outputs of our clustering formulation which can be represented thanks to a positive function f chained with a sum-1 normalization layer:

∀k ∈ 1, K τ k (x) = f k (x) ∑ K ℓ=1 f ℓ (x) (68)
On the Kantorovich side, the critics C k are also neural networks (of parameters θ C k ) but with a 1-Lipschitz property gracefuly provided by spectral normalization for each linear (or convolutional) non-activation layers thanks to a simple yet efficient power iteration technique developped recently by [START_REF] Miyato | Spectral normalization for generative adversarial networks[END_REF] which gives enough stability to face large scale datasets in the work of [START_REF] Brock | Large scale GAN training for high fidelity natural image synthesis[END_REF]. Injecting all these formulas in L n = L OvR n , the maximization over all parameters θ (concatenating the θ f k s and the θ C k s) becomes (in several lines):

L OvR n (p 1 , . . . , p K , π 1 , . . . , π K ) = E x∼p K ∑ k=1 τ k (x) -π k × C k (x) (69)
while proportions π k are equal to E x∼p P(c = k|x) which is kept and main- tained à la online k-Means [START_REF] Bottou | Convergence Properties of the K-Means Algorithms[END_REF] in their means updates (which can be seen like an Optimal Control [START_REF] Bertsekas | Dynamic programming and optimal control[END_REF] self-regulated closed loop mechanism). Suppose that an oracle gave us the optimal (π * , θ * E , θ * M ) and we keep them fixed, does adding a constant of the C k s critics functions could diverge to ±∞?

The answer is no because, by definition of the proportions in the mixture, we know that E x∼p [τ k (x)] = π k which means that the sign of τ k (x)π k cannot be constant (it can be always zero but that is an easy case for what we need to prove). Thus, adding a constant b k to C k will not change the objective. (and adding a constant multiplier a k to a valid 1-Lipschitz C k would violate that 1-Lipschitz constraint).

More generally, the fact that we maintain the null equalities

∀k ∈ 1, K E x∼p [τ k (x) -π k ] = 0 (70)
(thanks to [START_REF] Bottou | Convergence Properties of the K-Means Algorithms[END_REF]) and the Lipschitz property together prevent the critics C k s f rom diverging;

The objective of Eq. ( 69) favors τ k (x) to be far from π k thanks to the maximization: (i) if cluster memberships τ k (x) are too close to their means π k , then the objective would be close to zero (its lower bound because a positive linear combination of Wasserstein distances Eq. ( 69) is non-negative); (ii) when τ k (x) is high above its mean π k (bounded by 1), C k (x) will be high and when τ k (x) is low under its mean (bounded by 0), C k (x) will be low too. Accordingly, it is interesting to see that C k (x) can be seen as a relaxed decision function (i. e. high for points in k th cluster and low for the other clusters, but bounded in terms of variation due to its Lipschitzian property);

Overlapping Gaussian components of mixture M are avoided. Intuitively, if we take a region of the data space where τ k (x) and τ ℓ (x) (with k = ℓ) are high (meaning E (x) is on an overlapping zone between two Gaussian components k and ℓ), then the Wasserstein distances W(p k , pk ) and W(p ℓ , pℓ ) (parts of the objective sum) could have been even more maximized on this region because they are related to W(p k , p ℓ ). Thus our algorithm favors partitions over covers (in spite of our soft relaxation of memberships).

Algorithm

Algorithm 4 optimizes Eq. ( 69) which has only a maximization steps which provides a certain engineering ease compared to the min-max optimization scheme required for the GANs. Nevertheless, initialization routines turned out to be crucial in practice. For this highly non-convex optimization, suggest three successive careful initialization steps before the real opitmization.

The optimization of Eq. ( 69) that we want to accomplish has no guarantee to converge in a global extrema with our neural networks approach, which is why careful initialization is crucial. That is why our DiWaC algorithm can be decomposed in four successive steps (three for intializations and one for the real optimization). In theory, if we had a universal (i. e. regardless of the convexity) ideal optimizer at our disposal, then the three first steps would be useless but to avoid useless and spurious local minima, we do a three steps initialization for a final one.

1. Auto-encoder initialization We train a classic auto-encoder (E , D) (for an encoder E and decoder D parametrized respectively by θ E and θ D ):

(θ 0 E , θ 0 D ) = arg min θ E ,θ D E x∼p D • E (x) -x 2 2 (71)
using an auto-encoder initialization is interestingly reminiscent to 1990s ans 2000s pretraining for neural networks at a time when training time was a much more a burden than today (although [START_REF] Goodfellow | Deep Learning[END_REF] also mention better activation functions, weight initialization, variants of gradient descent and regularization methods, we empirically observed that pretraining is still benefitial for many scenarios both in terms of performance and sometimes even overall training time).

EM-based Gaussian mixture codes initialization

We fit a Gaussian mixture model M of proportions ω, means µ and covariance matrices S with the Expectation-Maximization algorithm [START_REF] Dempster | Maximum likelihood from incomplete data via the EM algorithm[END_REF] on the encoded data E (x):

θ 0 M = arg max θ M E x∼p log K ∑ k=1 ω k × N (E 0 (x); µ k , S k × S ⊤ k ) (72)
where θ M = (ω k , µ k , S k ) k=1,...,K . We also precise that this EM-GMM op- timization is itself initialized with k-Means++ [START_REF] Arthur | k-means++: The advantages of careful seeding[END_REF] 4 on the codes (or encoded inputs that lie at the bottleneck of the previous step's autoencoder).

Critics initialization

We define K critics functions (C k ) k=1,...,K implemented thanks to [START_REF] Miyato | Spectral normalization for generative adversarial networks[END_REF] and parametrized by θ C that would estimate the Wasserstein distances weighted sum with memberships probabilities provided by the previous step's EM. Indeed, we define clustering probabil-

ity functions τ k (x) = ω k ×N (E (x);µ k ,S k ) ∑ L ℓ=1 ω ℓ ×N (E (x);µ ℓ ,S ℓ )
to maximize over θ C :

θ 0 C = arg max θ C E x∼p K ∑ k=1 τ k (x) -π k × C k (x) (73)
Indeed, we recall that the critics functions C k s are just a convenient tool to estimate the Wasserstein distances. With bad critics, the gradient of the objective over the memberships probabilities parameters would be also wrong which is bad news especially when memberships functions are previously and nicely initialized. This suggests this warmup step 3 before the real clustering optimization.

Without that warm-up step 3 (right before the real step 4), we would deteriorate the previous clustering initializations quality provided by steps 1 and 2 (that is dubbed the "AE + GMM" initialization by [START_REF] Xie | Unsupervised deep embedding for clustering analysis[END_REF]) because the critics would not be trained enough to estimate that good clustering while still inflicting ignorant gradient steps on the clustering functions. To avoid such a bad scenario after steps 1 and 2, the last initialization step 3 corresponds to a warm-up before the real optimization step 4. This way, thrice the encoder, the mixture and the critics are reasonably well initialized and the real core step 4 can begin where all three are not free but just nicely initialized, relaxed and further optimized.

Core clustering

With the same objective, we do the core, final and wellinitialized optimization:

θ = arg max θ E x∼p K ∑ k=1 τ k (x) -π k × C k (x) (74)
as described in Algorithm 4.

Model Selection for DiWaC

The goal of model section is to check under/over-fitting and ultimately choose hyper-parameters such as the number of classes or the architecture of neurons and layers among several trained models with different hyper-parameters. In this work, our inter-cluster Wasserstein distance measured in Eq. ( 69) can be measured on some held-out validation data distributions. As stated earlier, critics are just a convenient tool to estimate the Wasserstein distances, so there are at least two ways to estimate the validation objective:

1. Directly (D) "as is": we simply apply Eq. ( 69) except for adapted proportions π measured on these held-out data:

π validation ← 1 N validation N validation ∑ i ′ =1 τ(x i ′ ) (75)
2. Re-Fitted (R) "with re-optimized critics": we use the same previously adapted proportions Eq. ( 75) and we optimize again the objective Eq. ( 69) but only with respect to the critics which corresponds to step 4, Eq. ( 73).

The goal is to refine the objective which corresponds to a separation power of the clustering in terms of interpretation.

(x i b ) b=1,...,B where i b ∼ U N (1, N) 9: Compute g i b k ← ω k × N (E (x i b ); µ k , S k ) ∑ K ℓ=1 ω ℓ × N (E (x i b ); µ ℓ , S ℓ ) 10: so that τ i b k ← g i b k ∑ K ℓ=1 g i b ℓ 11: Perform a gradient ascent step of 1 B ∑ B b=1 ∑ K k=1 (τ i b k - π k ) × C k (x i b ) to update θ 12:
Update proportions à la online k-means [START_REF] Bottou | Convergence Properties of the K-Means Algorithms[END_REF]

π ← π + 1 T + 1 × 1 B B ∑ b=1 τ i b -π 13: end while
Traditionally, model selection is accomplished according to the likelihood (or completed likelihood) related to the Kullback-Leibler divergence (see the Bayesian Information Criterion BIC [START_REF] Schwarz | Estimating the dimension of a model[END_REF] or Information Completed Likelihood ICL [START_REF] Biernacki | Assessing a mixture model for clustering with the integrated completed likelihood[END_REF] for examples of extrapolated generalization power measurements). One of the originality of our approach is that our model selection technique is done according to an other divergence which is the Wasserstein distance. This parallel allows us to use the historical likelihoodbased literature by replacing the well-known Kullback-Leibler divergence by the Wasserstein distance to maybe open new avenues of research for future investigations (e. g. adapting BIC and/or ICL beyond likelihood and Kullback-Leibler divergence).

In practice here, we train M times our model with different set of hyperparameters (different number of clusters, different neural networks structure, different mixture components structure etc.) Models ( θm ) m=1,...,M are now evaluated on held-out validation data this time instead of training data (in which they are trained as usual). Finally, we select the best set of parameters: the model indexed by m for which the objective (defined in Eq. ( 69) with adapted proportions) is maximum.

Changing the Metric beyond the Euclidean Distance for DiWaC

In terms of ill-posed problem, the axioms of [START_REF] Kleinberg | An impossibility theorem for clustering[END_REF] are almost satisfied: scale-invariance is valid, cluster-shapes-invariance is almost valid up to the expression power of the encoder and critic functions (a. k. a. the capacity of these neural networks in practice) that can be improved thanks to model selection and only metric-invariance remains but we have the intuition that it can be improved in future works thanks to the notion of worst metric among a large class of metrics as briefly evoked in the next chapter. Indeed, at the beginning of that clustering chapter, we mentioned that our clustering algorithm could not achieve consistency (i. e. metric invariance) because optimal transport required an initial and definitive commitment for the unique distance choice in the data space to build a Wasserstein distance for distributions over the data space. In fact, we could apply some ideas of the next contribution of this thesis to handle a very large class of distances at once which is possibly better than only one euclidean distance. Of course, according to the previously cited clustering nofree lunch theorem [START_REF] Kleinberg | An impossibility theorem for clustering[END_REF] 5 , we will never be able to cover all possible distances. Nevertheless, we could contemplate a solution where our hereby euclidean distance clustering is an initialization for a newer and better clustering technique for a large class of distances simulataneously in the future towards more consistency with the work of [START_REF] Kleinberg | An impossibility theorem for clustering[END_REF] in mind.

E X P E R I M E N T S

Implementation details and experimental setup

We did our experiments in Python by using the pyTorch [START_REF] Paszke | K. Pearson. Liii. on lines and planes of closest fit to systems of points in space[END_REF] and scikit-learn [START_REF] Pedregosa | Scikit-learn: Machine learning in Python[END_REF] libraries. with the same learning rate of 10 -5 with the Adam default optimization strategy [START_REF] Kingma | A method for stochastic optimization[END_REF] everywhere. k-Means and GMM are not easily compatible with large scale datasets which is why we took only a reasonable subset of large datasets for these initializations (which is why this is not crucial). In fact, the online learning of k-Means [START_REF] Bottou | Convergence Properties of the K-Means Algorithms[END_REF] is also possible for EM thanks to [START_REF] Cappé | On-line expectation-maximization algorithm for latent data models[END_REF] but we found our results already satisfying. In our preliminary experiments, the "AE + GMM" baseline (just the first 2 steps of our algorithm) performed poorly (and with bad reliability not reproducibility) without k-Means++ [START_REF] Arthur | k-means++: The advantages of careful seeding[END_REF] and Xavier neural networks weights initialization [START_REF] Glorot | Understanding the difficulty of training deep feedforward neural networks[END_REF]. Thus, all the experiments we report in this work use them.

For optimization reasons (unconstrained or implicitly constrained optimization is more stable than explicitly constrained optimization especially for stochastic gradient descent), we use two tricks:

• the SoftMax trick proportions are paremetrized by free logits that are converted into proportions through a SoftMax function 6 ;

• the Cholesky trick each covariance matrix is parametrized by its square root which always exists in the Cholesky decomposition sense for any covariance matrix (because it must be symmetric definite and positive, see [START_REF] Press | Numerical recipes 3rd edition: The art of scientific computing[END_REF]) and is a lower-triangular matrix whose diagonal coefficients are strictly positive (which can be ensured thanks to the use of the exponential function) to guarantee that when multiplied by its transposed version we get the correct covariance properties throughout the optimization path.

Unfortunately, we empirically realized that this kind of parametrization is not enough because of the eigenvalues behavior of that square root matrix: they numerically explode or implode resulting in ill-conditioned corresponding covariance matrices and create instability. More precisely, it appeared that when the means are far from the optimal means, the system chooses to modify its covariance matrices first instead of the means which should be prevented because of the spurious maxima. Limiting the eigenvalues range is better than an unbounded exponential function.

We use a distorted Sigmoid 7 function bounded by two constants initially

6 SoftMax(v) k = exp(v k ) ∑ K k ′ =1 exp(v k ′ ) 7 Sigmoid(t) = 1 1+exp(-t)
parametrized by the variable-wise standard deviations σ j of all initial codes variables E 0 (x i ) j :

(∀t ∈ R) λ + (Λ -λ) × Sigmoid(t) ∈]λ, Λ[ (76)
where:

λ = 0.3 × min j∈{1,...,d} σ j and Λ = 3 × max j∈{1,...,d} σ j (77)
to be read with

m j = 1 N ∑ N i=1 E 0 (x i ) j and σ j = 1 N ∑ N i=1 (E 0 (x i ) j -m j ) 2 .
In supervised classification, for a known labeled dataset (but hidden to the trained system), the confusion matrix (a. k. a. matching matrix, contingency matrix or table of confusion) is counting for every pair (actual class k, predicted class ℓ) the number of occurences of points falling into that pair configuration (being in group k and category ℓ). Thus, that table can be an interpretation and visualization tool after a supervised classification or an unsupervised clustering (with ground truth) on data.

From the confusion matrix M between a given clustering (indexed by the rows k) and a revealed ground truth (indexed by the columns ℓ), several accuracy measures can be built afterwards:

H U N G A R I A N A C C U R A C Y ( A C C )
Thanks to the Hungarian Method [START_REF] Kuhn | The hungarian method for the assignment problem[END_REF][START_REF] Stephens | Dealing with label switching in mixture models[END_REF] run on the confusion matrix (while taking care of the sign thanks to a minus sign and adding the maximum entry because the Hungarian Method minimizes a sum and we want to maximize that a sum of occurences), we can measure how good is a clustering

ACC = max σ K ∑ k=1 M k,σ(k) ⇐⇒ min σ K ∑ k=1 C k,σ(k) ( 78
)
where the cost matrix C is adapted with these entries

C k,ℓ = (max k ′ ,ℓ ′ M k ′ ,ℓ ′ ) - M k,ℓ
and an optimal K-permutation σ tells which cluster σ(k) should cluster k be assigned to. In the case where we need to over-cluster (a. k. a. over-segment the dataset) which is predicting more groups than ground truth semantic ones, it is interesting to see if the several small predicted clusters fit the fewer and bigger semantic human ones. In this case, normalized conditional entropy NCE seems appropriate:

NCE = CE H ( 82
)
where the conditional entropy CE is:

CE = ∑ k,ℓ Mk,ℓ log Mk,ℓ Mk,• (83)
which makes NCE a kind of an asymetric version of NMI.

Introductory Examples for Generative Wasserstein Clustering

The probabilistic perspective tackled by Wasserstein-GAN in [START_REF] Arjovsky | Wasserstein generative adversarial networks[END_REF] approach can inspire one to think of the Kullback-Leibler alternative with the famous Expectation-Maximization algorithm for Gaussian mixture models (EM-GMM) [START_REF] Dempster | Maximum likelihood from incomplete data via the EM algorithm[END_REF]. Let's try to solve the same problem as EM-GMM by using the Wasserstein distance with Wasserstein-GAN instead of the Kullback-Leibler divergence (or likelihood up to a sign and a constant term) with EM: take two Gaussian mixtures A and B both of K = 3 components each in a D-dimensional space (D = 2) defined by:

• K proportions (π A k ) k=1,...,K that are positive and sum to one for A and K proportions (π B k ) k=1,...,K with same properties; • K means (µ A k ) k=1,...,K for A and K means (µ B k ) k=1,...,K for B where the means live in the same 2D space R D ;

• K full covariance matrices (Σ A k ) k=1,...,K that are symmetric definite positive in R D×D for A and K full covariances (Σ B k ) k=1,...,K with same properties for B. Now for good convergence at step 4, we decompose the successive minimizations and maximizations by just focusing on the maximization steps of Wasserstein-GANs in a warm-up phase. Indeed, in section 2.3.4, after steps 1, 2 and 3, we believe that the decoder and the mixture are well-initialized but the critic neural network is not. It means that our Wasserstein distance estimator (provided by the critic neural network) is not good and neither its gradient with respect to the decoder and the mixture. Concretely, this means that even though the generator (i.e. mixture and decoder) is well initialized, the first iterations will decrease the initial generation quality because of a poorly initialized critic that computes our Wasserstein distance value and gradients.

To cope with this bad critic initialization problem at step 4, we choose to optimize the critic alone which is estimating the Wasserstein distance between 2 fixed distributions (real and generated data do not change but the critic does) at the beginning of that step 4. This warm-up step concludes a nice initialization for every parameter as you can see in the toy example represented in Fig. 13 and the step 4 finally converges into Fig. 14. Of course, this warm-up step initially developed without generator neural network is also used in our case of WAMiC with the mixture and the decoder as the two components of the generator: first, we do not change the generator and we only change the critic neural network until convergence and then we trigger a second step after this warm-up and the actual min-max Wasserstein GAN optimization scheme takes place.

To illustrate our approach further, we now work on a toy dataset that we call "Three Moons" with 1000 points for each of the 3 groups in 2 dimensions as presented in Fig. 16. With a code space of dimension 1, even though the clusters are not linearly separable in the original data space, our system is able to cluster them in 3 groups with 100% of accuracy.

For that simple toy dataset, most of the clustering work is done by the vanilla auto-encoder. Indeed, once the auto-encoder is trained, we observed that the 1D codes are already separated according to the cluster labels. Here, we just wanted to see if model selection was plausible in a simple scenario.

For the number of clusters, without labels in our unsupervised context, while a classification-score-based cross-validation is not an option, one can still measure an estimate of the Wasserstein distance (our loss) in step 3 of section 3.3 but on a validation set. The intuition behind is that if we did not overfit, our loss function will still be satisfactory on that validation set (that was not seen during training). After running our first 3 steps algorithm out of 4, we can train a new neural network f to measure the Wasserstein distance between a held-out validation dataset and some generated data empirical distributions. More precisely, we find the results presented in Fig. 15 actually selecting 3 clusters which is satisfactory.

Introductory Examples for Discriminative Wasserstein Clustering

To investigate model selection capabilities of our method and for explanatory reasons, we use here 2 synthetic datasets:

T H R E E M O O N S
3 M To illustrate our approach further, we now work on a toy dataset that we call "Three Moons" with 1000 points for each of the 3 First, as described in section 2.3.4, we trained an auto-encoder on the 3M dataset with a code space dimension of d = 1: even though the clusters are not linearly separable in the original data space, our system is able to cluster them in 3 groups with 100% of accuracy (an EM-GMM or a simpler k-Means run in the codes data is enough!). For that simple toy dataset, most of the clustering work is done by the vanilla auto-encoder. Indeed, once the auto-encoder is trained, we observed that the 1D codes are already perfectly separated along one axis with respect to the cluster labels. Here, the point is not to check clustering capabilities but to see if model selection is plausible in this simple scenario. Now we study the case where that number of clusters K is not known in order to check if model selection is possible in these relatively easy settings. The intuition behind our validation procedure for model selection is that if a given hyper-parameter configuration (the number of clusters in this experiment) does not lead to under/over-fitting, our objective function will still be satisfactory on a validation set (i. e. that was not seen during training) compared to other configurations with one held-out validation set.

In the little more challenging (but still synthetic) VD dataset in Fig. 17, we observe that dealing with various kinds of distributions (even unknown by the system and beyond the Gaussian case) is still compatible with our mixture-type model-based clustering algorithm.

Real Data Experiments

For the real-world experiments, we were interested in 4 datasets with 3 different in nature (images, sparse and dense data) but all with several hundreds of dimensions per item:

• MNIST: 70 000 handwritten digits images dataset living in dimension 784 (for 28 × 28 pixels);

• Reuters: English news stories labeled with a category tree [START_REF] Lewis | Rcv1: A new benchmark collection for text categorization research[END_REF].

Following DEC Xie et al. [2015], we used 4 root categories: corporate/industrial, government/social, markets, and economics as labels and discarded all documents with multiple labels. We computed tf-idf features on the 2000 most frequent words to represent all articles;

• Reuters-10k: a random subset of Reuters with only 10 000 examples (selected with precisely the same random generator seed as DEC);

• HHAR: The Heterogeneity Human Activity Recognition (HHAR) dataset [START_REF] Stisen | Smart devices are different: Assessing and mitigatingmobile sensing heterogeneities for activity recognition[END_REF] contains 10,299 sensor records from smart phones and smart watches. All samples are partitioned into 6 categories of human activities and each sample is of 561 dimensions.

On MNIST, in Fig. ( 18), we generated data from our GeWaC further and further in random directions from the centroids: we see that the digits get fancier away from the centroids. The good quality of the generation is comparable to those of regular GANs but in a cluster-wise fashion.

In these experiments, we used the same encoder (with symmetric decoder) MLP 8 structure from Xie et al. the code space) (and ReLU activations) for fair comparisons with others. In fact, our conditions are tougher than the ones of IMSAT, GAR and DEPICT that use more sophisticated convolutions than matrix-vector products while DEC, VaDE and our DiWaC do not which is very encouraging.

k m k +0.5S k u m k + S k u m k +1.5S k u m k +2S k u m k +2.5S k u m k +3S k u m k +3.5S k u
The overall results of our DiWaC approach compare favorably to the deep clustering state-of-the-art in Table . 2. First, we observe that there is an improvement over standard baseline algorithms ("AE + GMM" compared to GMM and even k-Means) when fed with the output of an AE which is coherent with the results of [START_REF] Xie | Unsupervised deep embedding for clustering analysis[END_REF]. We must admit how surprised we are that this "AE + GMM" baseline works so well on various datasets. Undiscovered mathematical theories could tackle this phenomenon: Why would a simple auto-encoder consistently map data in separated Gaussian groups?. Maybe the auto-encoders layers loosely behave like successive cluster-wise deterministic random projections [START_REF] Bingham | Random projection in dimensionality reduction: applications to image and text data[END_REF] and non-linear functions: indeed, the linear layers matrix entries have Gaussian frequencies but this would require further investigations to get thorough scientific interpretations. On a more practical level, we are particularly interested in real-world industrial cases where we highly recommend this "AE + GMM" baseline for its ease of implementation and good results in practice and also for fast prototyping.

There is a supplementary and significant improvement for DiWaC over its "AE + GMM" initialization. These good results are even comparable to those of supervised non-convolutional networks with just a few layers of the 1990s9 . On MNIST, the important performance gap between VaDE and DiWaC can be explained by difficulties (already studied in [START_REF] Jiang | Variational deep embedding: A generative approach to clustering[END_REF]) with data that lie extremely close to low-dimensional manifolds, like images. In that regard, our algorithm inherits the strenghts of WGAN [START_REF] Arjovsky | Wasserstein generative adversarial networks[END_REF] and models data much more faithfully without our DiWaC algorithm being adversarial which means that Optimal Transport alone is a powerful tool as the non-adversarial work of [START_REF] Genevay | Learning Generative Models with Sinkhorn Divergences[END_REF] also tends to proove. Furthermore, as presented previously (section 2.3.5), our strategy has a natural criterion for model selection whereas VaDE has no principled model selection criterion; in theory, the validation likelihood could be used, however, computing the likelihood of for variational techniques remains an open question (see e.g. [START_REF] Cremer | Inference suboptimality in variational autoencoders[END_REF]). In other works, DEC's authors propose to use the ratio between training and validation loss as a criterion, however this ad hoc solution has little statistical foundations. Model selection is definetely the main feature of our approach on top of providing good results at same hyper-parameter configurations. In fact, better neural networks structure can be found by cross-validation with our algorithm which could be done in future work.

F U T U R E W O R K A N D C O N C L U S I O N

This work presents twos ways to use the Wasserstein distances literature with neural networks in order to achieve efficient clustering at the crossroad of Optimal Transport, Neural Networks and Model-based mixture techniques: first, we did generative clustering and second, we did discriminative clustering. In spite of three careful initialization steps, our technique is still end-to-end trainable as theoretically the final step could be done alone but in practice, with an associated higher risk of being stuck in spurious local extrema and much less reproducibility (i. e. higher dependence on the pseudo-random generators seeds which we do not have thanks to our robust and reasonable initialization steps).

In a mockingly accurate description of the generative part of our work, we could say that placing a mixture distribution at the input of new Generative Adversarial Networks seemed worthwhile for clustering purposes. Indeed, we inherit the recent wealth of literature (see adversarial autoencoders by [START_REF] Makhzani | Adversarial autoencoders[END_REF], Wasserstein autoencoders by [START_REF] Tolstikhin | Wasserstein auto-encoders[END_REF] and even adversarially learned inference [START_REF] Dumoulin | Adversarially learned inference[END_REF]) with tunable input mixtures each mode corresponding to a cluster. In fact, it turns out that generating data in a clustered fashion with the desire to minimize the estimated Wasserstein distance with real data is difficult because of a mono-component collapsing effect due to a LASSO-kind of constraint on the mixture proportions which favors degenerated (sparse) proportions. Another way to explain this failure is the richness of the neural networks expression power that makes the Gaussian mixture clustering useless: one Gaussian prior is enough to generate whole datasets which looses any clustering capabilities hope if done carelessly. Eventually, we acknowledge the better and elegant solutions given in the recent work of [START_REF] Mukherjee | ClusterGAN: Latent Space Clustering in Generative Adversarial Networks[END_REF] to cope with these problems and nothing can prevent us from re-using the same techniques with Wasserstein distance for us in lieu of the Jensen-Shannon divergence they use.

For the discriminative part, we observe that although our DiWaC system sees the data one by one because of its linear computation and memory complexity, it handles complex relationships between the points. In a way that is similar to spectral clustering (SC) with pairwise similarities in particular, but our approach has improved results with an objective function that does not full storage of pairwise similarities but only optimal transport that encodes it implicitly into the systems with still pairwises similarities but at a cluster distributions level. Our algorithm detects the space zones where the density mass is occupied thanks to the probabilistic model inherent to the optimal transport theory associated to our optimization. For that matter, we did manage to handle the trade-off between the flexibility of the neural networks functions and the rigidity of the mixture models choice of distributions. A counter-example would consist in a dataset with non-separable classes: the system would not know precisely when to activate the low-density threshold in order to optimally put a class frontier.

Although our work is far from the well-established theory of the Reproducing Kernel Hilbert Spaces (RKHS), we realize that a similar phenomenon occurs when it comes to satisfactory results: from an input space of a given reasonable dimensionality D, it is better to first increase that dimensionality by changing space thanks to a mapping (the first layers of our encoder are increasing dimensionality compared to D) and only after that, shrink it to a small dimensionality d ≪ D (the last layer of our encoder). The main difference between our work and the RKHS theory is that our transformations are learned from the data and not given analytically a priori (e.g. by a Gaussian kernel) which echoes the work of [START_REF] Unser | A representer theorem for deep neural networks[END_REF].

Within the context of the algorithm laid out above, we empirically observe some symbiosis operating between generative or discriminative clustering and non-linear embedding. A great advantage of our approach is model selection especially for the number of clusters that seems to empirically work well. In other tasks different from clustering for itself, encouraged by future over-clustering investigations, one can be inspired by the work of [START_REF] Liao | Learning deep parsimonious representations[END_REF], that demonstrates the benefits of clustering-based regularization for supervised classification and generalization. For example, inside a supervised classification problem, identifying sub-categories might help to specialize classifiers on more homogeneous classes for improved generalization capabilities and ease of data interpretations. Finally, we now have the intriguing possibility of following other successes of the supervised classification and unsupervised clustering communities, such as the automatic selection of discriminative parts in the spirit of what [START_REF] Sun | Learning Dictionary of Discriminative Part Detectors for Image Categorization and Cosegmentation[END_REF] and eventually [START_REF] Doersch | What makes paris look like paris[END_REF] did and the automatic setting of the number of clusters. Beyond all these promising clustering results, our GeWaC algorithm has the ability to generate cluster-wise data. It is interesting to see that such conditional generation, already explored in the supervised setting [START_REF] Mirza | Conditional generative adversarial nets[END_REF] is not fundamentally harder without supervision.

U N S U P E RV I S E D F E AT U R E I M P O R TA N C E

What Makes Paris Look like Paris?

probably Alexei A. Efros in SIGGRAPH, 2012 A B S T R A C T
Machine learning and pattern recognition requires data analysis: a process of inspecting, cleansing, transforming and modeling data with the goal of discovering useful information, informing conclusion and supporting decision-making. Unsupervised learning as a scientific field provide tools for dimensionality reduction, visuzalization procedures, features extraction etc. in order to empower human beings with enough computational power to grasp the environment we live in. This chapter is more of plea for further investigations towards unsupervised feature importance rather than a scientific contribution. Indeed, we revisit notions such as differenciation for distributions, distribution Wasserstein-based metrics and manifold normal in order to call for both further theoretical and practical research work.

K E Y W O R D S
Feature importance, relevant features, manifold, foreground/background segmentation, cosegmentation, hypersurface normal

I N T R O D U C T I O N

This chapter aims at exploring possible statistical and algorithmical tools for unsupervised feature importance extraction. Indeed, the problem of extracting the absolute or relative importance of data coordinates is of high interest for understanding data. In supervised learning, revisiting the principle of Occam's razor gave birth to a considerable literature and we suggest at least one thesis to the reader: Structured Sparsity-Inducing Norms: Statistical and Algorithmic Properties with Applications to Neuroimaging, the Ph.D. manuscript of [START_REF] Jenatton | Structured sparsity-inducing norms: Statistical and algorithmic properties with applications to neuroimaging[END_REF] which combines coordinate selection, weighting and structure in high-dimensional problems such as neuro-imaging. At the same time, the decision trees literature also carefuly studied the problem in several works by [START_REF] Breiman | Random forests[END_REF][START_REF] Breiman | Classification and regression trees[END_REF] and more recently, there is a will to break the black box taboo of neural networks being supposedly non-interpretable with interesting attempts by [START_REF] Knight | Neural networks debunk phrenology[END_REF] and de [START_REF] De Sá | Variance-Based Feature Importance in Neural Networks[END_REF]. The emerging popularity of add-on toolboxes such as Captum [START_REF] Kokhlikyan | Pytorch captum[END_REF] is a solid proof showing a research trend towards input space interpretability in supervised settings at least.

Back in unsupervised learning, the problem of selecting or weighting coordinates by relevance is probably an ill-posed problem because there is no supervision. As usual in pattern recognition, we assume that data live in an instrisically low-dimensional manifold compared to the whole data space dimensionality. If we are able to find a machine learning procedure able to compute an hypersurface normal of that manifold at each point of it, then we can interpret that normal direction coordinates as relevant or not for describing the manifold. More precisely, we can ask ourselves Is it possible to maximally change a manifold of data with an infinitesimally small distortion? and the distortion would be a function of space giving high amplitudes to coordinates that one should not change in order to preserve the manifold consistency. This question is reminiscent to the notion of gradient and we make the hypothesis that perturbating data in an infinitesimally small fashion can be done with gradient of a Wasserstein distance between the real data distribution and a pertubated version of it.

One possible application of this work could be unsupervised foreground / background segmentation. Indeed, from a dataset of images containing the same high-level semantic category of content (e. g. "wolves") in several outdoor / indoor conditions, the revealed coordinates would select foreground pixels from background non-content-manifold-specific pixels (that can intuitively be changed without breaking the semantic meaning of the image category). Generalizing such an automatic tool would be of great interest in many scientific fields beyond computer vision.

I N WA M A D I : I N FI N I T E S I M A L WA S S E R S T E I N M A X I M A L D I S T O R T I O N

In the previous chapter 2, we reviewed some consequences of the impossibility theorem by [START_REF] Kleinberg | An impossibility theorem for clustering[END_REF]. In particular, the metric invariance is an interesting and difficult subject. Indeed, it seems that choosing a particular metric is a heavy commitment. This is especially true in unsupervised learning probably because determining a metric is choosing the algorithms lenses for seeing the data without supervision which is redefining a notion of neighborhood tainted by the curse of dimensionality we mentioned earlier in introduction section 1.5.2.

For a random variable x coming from distribution p living in a space X (say R D ), we can consider the distortion function D = t × F mapping X to X (and t ∈ R + ) which creates a second random variable y defined by:

y = x + D(x) = x + t × F (x) (84)
which defines a new distribution q t,F . For the sake of simplicity, we impose:

(∀x ∈ R D ) F (x) 2 = 1 (85)
so that the length of the distortion is simply t. Our goal is to measure how different a perturbated distribution q t,F can be from the original distribution p with an infinitesimal length t and constrained energy F (x) 2 = 1. Thanks to an already successful probabilistic approach in Machine Learning [START_REF] Murphy | Machine Learning, a Probabilistic Perspective[END_REF], we rephrase our question set out in our introduction: What is the infinitesimal steepest distortion of data? Indeed, this kind of approach would give a function F such that when applied to each data point x i ∈ R D , the computed vector F (x i ) ∈ R D would tell which coordinate (x (j) i ) j=1,...,D is relevant i. e. characteristic in an interpretable way for deeper data analysis especially when the dimensionality D is high.

There is a natural mathematical and geometric tool to measure a distortion for distributions: the Wasserstein distance when the associated data space metric is d. Thus, inspired by the optimization idea of steepest gradient direction, we can first define such a function measuring a quantity corresponding to the discrepancy L(t, F , d) induced by the distortion D = t × F , namely the Wasserstein distance (based on metric d) between the distributions p and q t,F :

L(t, F , d) = min γ∈Γ(p,q t,F ) E (x,y)∼γ [d(x, y)] (86)
and because

y ∼ q t,F ⇐⇒ y = x ′ + t × F (x ′ ) and x ′ ∼ p (87)
we get (without applying a change of variable formula involving a Jacobian term):

L(t, F , d) = min γ ′ ∈Γ(p,p) E (x,x ′ )∼γ ′ d(x, x ′ + t × F (x ′ )) (88)
with Γ(p, q) being the set of coupled distributions with marginals p and q. Please note that without the min operator, this change of variable would not have been valid: between Eq. ( 86) and Eq. ( 88), the optimal transport changed from γ to γ ′ because they are different since Eq. ( 87). Now we revisit the notion of steepest function Wasserstein direction gives an optimization problem:

F * d = arg max F |∇ t L(0, F , d)| (89)
where for a given direction F and metric d, the quantity ∇ t L(0, F , d) is the derivative of function t → L(t, F , d) on t = 0 + which is:

∇ t L(0, F , d) = lim t→0 + L(t, F , d) -L(0, F , d) t -0 (90)
which simplifies in:

∇ t L(0, F , d) = lim t→0 L(t, F , d) t (91)
because p = q 0,F for all direction F and so, we get:

F * d = arg max F lim t→0 + L(t, F , d) t (92)
At this point, we ignored the distance d operating in the data space X = R D but we can build such a distance in a form that parses a large variety of metrics:

d φ : R D × R D → R + (93) (x, y) → φ(y) -φ(x) 2
and we note that

d φ (x, y) = (L 2 • φ)(x, y) = L 2 (φ(x), φ(y)) (L 2 being the eu- clidean distance).
We make sure that φ is a smooth bijection so that we inherit injectivity (and differentiability for optimization reasons we will see later). Indeed, such a bijective φ allows the associated function d φ to verify the distinguishability property of a distance, namely:

(∀(x, y) ∈ R D × R D ) x = y ⇐⇒ φ(x) = φ(y) (94) ⇐⇒ 0 = L 2 (φ(x), φ(y)) = d φ (x, y) (95)
and the other required properties for being a distance: positivity, symmetry and triangular inequality are given for free thanks to the euclidean distance. Moreover, we also want to avoid some equivalence class explosion effect due to the fact that there is no practical difference of interpretation between choosing a given metric d and a proportional one α × d (with α > 0) especially for explosively large coefficient α. Avoiding such annoying properties can be obtained by constraining the variations of the function φs indexing the distances space made of d φ s. Mathematically, the notion of variation for a multivariate bijective function φ is studied thanks to the derivative matrix ∇φ(x) ∈ R D×D called the Jacobi matrix at each point x ∈ R D and the main variation directions are given by its eigen values (λ(x) (j) ) j=1,...,D . Thus we can propose two ways to constrain these variations:

B O U N D I N G A M P L I T U D E ( B A ) We limit the Jacobi eigenvalues amplitude:

(∀x ∈ R D )(∀j ∈ 1, D ) λ(x) (j) ∈ [J min , J max ] ⊂ R * + (96) Z E R O I N G G L O B A L L O G -A M P L I T U D E ( Z G L A )
We enforce a null logarithm of Jacobi determinant mean:

0 = κ = E x∼U (p,q t,F ) log | det ∇φ( x)| = E x∼U (p,q t,F ) D ∑ j=1 log λ( x) (j) (97)
where:

x ∼ U (p, q t,F ) means x = u × x + (1 -u) × y (98)
with u ∼ U R (0, 1) and x ∼ p and y = x ′ + t × F (x ′ ) with x ′ ∼ p Indeed, we need that "zero overall Jacobian property" being true but only over the convex enveloppe of the original and distorted points and not necessarily on the whole space R D (manily because we will not evaluate the functions at hand anywhere else outside that enveloppe). This convex enveloppe sampling technique has been used for maintaining a Lipschitzian constraint by [START_REF] Gulrajani | Improved training of wasserstein gans[END_REF].

These two combined constraints over the set of smooth bijections defines the functions set Φ. For more mathematical details, we highly recommend the reader the thorough academic book by [START_REF] Ambrosio | Gradient flows: in metric spaces and in the space of probability measures[END_REF] to conduct further and more principled studies.

A Simplified Case: Empirical Distributions

Before doing some mathematical proposals, we study in this section a simplified case where only empirical distributions are at stake. For that simplified scenario with euclidean distance, we handle p = 1 B ∑ B b=1 δ x i b an empirical distribution from p, we also have qt,F = 1 B ∑ B b=1 δ x i b +t×F (x i b ) from q t,F . For all t ∈ 0, 1 2 min i,i ′ x i ′x i 2 the optimal transport between p and qt,F is the natu- ral one as the Fig. 19 shows (we will proove the general case later):

W L 2 ( p, qt,F ) = min γ∈Γ( p, qt,F ) E (x,y)∼γ y -x 2 (99) = min γ ′ ∈Γ( p, p) E (x,x ′ )∼γ ′ x ′ + t × F (x ′ ) -x 2
but with empirical distributions, the transport γ ′ is in fact an assignment π of integers i b parsing the points in the p = 1 B ∑ B b=1 δ x i b sum to the integers i b ′ parsing the points in the qt, Unfortunately, this is not useful because it means that all F directions are equally distorting the original distribution with respect to the euclidean Wasserstein distance in a quantity that does not even depend on the distribution p.

F = 1 B ∑ B b ′ =1 δ x i b ′ +t×F (x i b ′ ) with π(i) = i ′ : W L 2 ( p, qt,F ) = min π 1 B B ∑ b=1 x π(i b ) + t × F (x π(i b ) ) -x i b 2 (100) = 1 B B ∑ b=1 x i b + t × F (x i b ) -x i b 2 = t so that lim t→0 + W L 2 ( p, qt,F ) t = 1 for all directions F .
Based on that failed study, we decide to also optimize the distance d φ to get non trivial F directions of distortion because we just saw that the euclidean distance is independent from the manifold at hand and thus too much isotropic in a data-independent fashion. Optimizing the distance, gives a new optimization problem: max

F ,φ lim t→0 + W d φ (p, q t,F ) t (101)
instead of max F lim t→0 W L 2 (p,q t,F ) t = 1 that previously left us with no optimiza- tion hope. This Eq. ( 101) is much more powerful because we do not only get the steepest distortion but also the pair of steepest distortion direction F and optimized associated metric d φ .

Now that we have a better grasp on what is mathematically going on, we tackle the general (smooth) distributions case.

General Distribution Case

The ideas we just briefly tackled are appealing and now we give a rather theoretical result in order to pursue the optimization side in a general case beyond the euclidean distance.

Proposition. For an infinitesimally small distortion D = t × F , the optimal transport between the original distribution p and distorted distribution q t,F is the natural transport for all smooth metric d φ indexed by bijection φ:

(∃M ∈ R * + )(∀t ∈ R * + ) t < M =⇒ W d φ (p, q t,F ) = min γ∈Γ(p,q t,F ) E (x,y)∼γ d φ (x, y) = E x∼p d φ (x, x + t × F (x))
Proof. For all pair (x, y) of points drawn from any transport γ ∈ Γ(p, q t,F ), there exists a point x ′ such that y = x ′ + t × F (x ′ ) because this is how q t,F is built. Applying a Taylor expansion on the function f x,x ′ defined by:

f x,x ′ (t) = d φ (x, x ′ + t × F (x ′ )) (102)
gives:

d φ (x, y) = d φ (x, x ′ + t × F (x ′ )) (103) = f x,x ′ (t) = f x,x ′ (0) + t × ∇ t f x,x ′ (0) + o(t) = d φ (x, x ′ ) + t × ∇ t f x,x ′ (0) + o(t) (104)
Let's focus on the second term ∇ t f x,x ′ (0), we know thanks to a Taylor expansion on φ around x ′ that:

f x,x ′ (t) = φ(x ′ + t × F (x ′ )) -φ(x) 2 (105) = φ(x ′ ) -φ(x) + t × ∇φ(x ′ ) × F (x ′ ) + o(t) 2
and thus the only remaining non-negligible term depending on t gives:

∇ t f x,x ′ (0) = ∇φ(x ′ ) × F (x ′ ) 2 (106)
and in the end, for Eq. ( 103), we obtain:

d φ (x, y) = d φ (x, x ′ ) + t × ∇φ(x ′ ) × F (x ′ ) 2 + o(t) (107)
and from that Eq. ( 107) as for all transport γ ∈ Γ(p, q t,F ) there is an associated coupling γ ′ ∈ Γ(p, p) such that each pair (x, y) ∼ γ corresponds to the pair (x, x ′ ) ∼ γ ′ (as seen earlier) we get:

E (x,y)∼γ d φ (x, y) = E (x,x ′ )∼γ ′ d φ (x, x ′ ) (108) +t × E x ′ ∼p ∇φ(x ′ ) × F (x ′ ) 2 +t × ǫ F ,γ (t)
where lim t→0 ǫ F ,γ (t) = 0 and one can note that the second term decouples the (x, x ′ ) pairing. If we consider the natural transport γ * (i. e. (x, y) ∼ γ * ⇐⇒ (x ∼ p and y = x + t × F (x)), then we can measure the difference D(γ, γ * ) for any other non-natural transport γ:

D(γ, γ * ) = E (x,y)∼γ d φ (x, y) -E (x,y)∼γ * d φ (x, y) (109) = E (x,x ′ )∼γ ′ d φ (x, x ′ ) + t × ǫ F ,γ (t) -ǫ F ,γ * (t) thanks to Eq. ( 108
)
which is positive for a t > 0 sufficiently small because lim t→0 |ǫ F ,γ (t)ǫ F ,γ * (t)| = 0 which proves that the natural transport γ * has the minimal transport cost and thus is the optimal transport.

When we briefly saw the empirical distribution and euclidean case, it appeared that when the distortion is sufficiently low in amplitude (i. e. with a small t), the optimal transport is the natural one (assigning each original point to its distorted version). With that proposition above in mind and in order to build an optimization objective and an algorithm, we propose a procedure to get a value of such a t > 0 satisfying:

max b φ(x i b ) -φ(x i b + t × F (x i b )) 2 ≤ 1 2 min b,b ′ φ(x i b ) -φ(x i b ′ ) 2 (110)
The right term does not depend on t and is easily computed from mini-batches of data. The left term can be approximated to guess the right order of magnitude for t that we can divide later until the inequality Eq. ( 110) is satisfied in a dichotomic fashion. Indeed:

φ(x i b ) -φ(x i b + t × F (x i b )) 2 ≃ t × ∇φ(x i b ) × F (x i b ) 2 (111) ≤ t × J max because F (x i b ) 2 = 1 and thus t k = min b,b ′ φ(x i b ) -φ(x i b ′ ) 2 2 k+1 × J max (112)
with an increasing k ∈ N * is a good strategy until satisfying Eq. ( 110).

O P T I M I Z AT I O N

In this section, we present a draft of an optimization strategy to summarize the ideas we just presented with neural networks implementations for functions. This is about leveraging the research effort for deep learning in general and GANs in particular for our representation learning and data analysis purposes. Indeed, we have the function F implemented by a neural network of pararameters θ F mapping R D to R D . We add some layers: a SoftMax layer followed by an element-wise square-root layer with kept sign layer such that we maintain the norm 1 constraint on F .

Thanks to some work accomplished in a different context by [START_REF] Dinh | Density estimation using real nvp[END_REF], the implementation of bijection φ in a special neural network of parameters θ φ is already done and the constraints needed in section 3.2 are easily applicable. More specifically, we keep on zeroing the overall jacobian mean in a way that is similar to online k-Means [START_REF] Bottou | Convergence Properties of the K-Means Algorithms[END_REF] via an intermediate bijection ψ that we divide by scalar κ for control the variation of φ = exp( -κ D ) × ψ. All these parameters are concatenated in θ that the algorithm 20 optimizes.

Maintaining the constrains that we required Eq. ( 96) and Eq. ( 97) in section 3.2 is facilitated by the structure of our bijective neural network that we took from [START_REF] Dinh | Density estimation using real nvp[END_REF]. Indeed, for bijection ψ implemented by L bijective layers

ψ = Layer L • • • • • Layer ℓ • • • • • Layer 1 ,
we have at each layer ℓ among L:

x [:d] → y [:d ℓ ] = Layer ℓ (x) [:d ℓ ] = x [:d ℓ ]
(113)

x [d ℓ :] → y [d ℓ :] = Layer ℓ (x) [d ℓ :] = x [d ℓ :] × exp(s ℓ (x [:d ℓ ] )) + t ℓ (x [:d ℓ ] )
where the functions s ℓ s and t ℓ s are free neural networks (with pythonic notations for dimensions and without loss of generality in the coordinates order but with an arbitrary pivot d ℓ ) and then we can bound the associated Jacobi matrices:

∇Layer ℓ (x) = I d ℓ 0 d ℓ ×(D-d ℓ ) -diag(exp(s ℓ (x [:d ℓ ] ))) (114) and log | det ∇Layer ℓ (x)| = D-d ℓ ∑ j=1 s (j) ℓ (x [:d ℓ ] )
Thanks to the chain rule applied to such a bijection ψ as a composition of layers from x 0 = x to x ℓ = Layer ℓ (x ℓ-1 ) for ℓ ∈ 1, L , we can collect through the forward computations of the function ψ output and sum the outputs of s ℓ in order to get κ:

κ = E x∼U (p,q t,F ) [log | det ∇ψ(x)|] = E x∼U (p,q t,F ) L ∑ ℓ=1 D-d ℓ ∑ j=1 s (j) ℓ (x ℓ-1 [:d ℓ ] ) (115)
Having an updated estimate of κ allows us to use φ = exp( -κ D ) × ψ instead of ψ directly such that log | det ∇φ(x)| = log | det ∇ψ(x)|κ which as zero mean.

For numerical stability reasons, we can impose on the neural networks s ℓ s to have a final element-wise 3×tanh

∑ L ℓ=1 (D-d ℓ )
layer so that κ ∈ [-3, 3] is bounded and thus for the local stretching amplitude it gives: J ℓ min/max = exp( ±3

D-d ℓ

). Recently, some Generative Adversarial Networks articles by [START_REF] Pan | Recent Progress on Generative Adversarial Networks (GANs): A Survey[END_REF] and [START_REF] Detlefsen | Reliable training and estimation of variance networks[END_REF] report some evidence that the Kantorovich-Rubinstein formulation implemented by neural networks in Lipschitz functions [START_REF] Miyato | Spectral normalization for generative adversarial networks[END_REF] has positive regularization effects on the computation and optimization of the Wasserstein-based losses. To benefit from these advances, let's recall what we wanted: max

F ,φ lim t→0 + W d φ (p, q t,F ) t (116)
Thanks to that Kantorovich-Rubinstein duality (used for example for Wasserstein Generative Adversarial Networks first by [START_REF] Arjovsky | Wasserstein generative adversarial networks[END_REF]), we know that for sufficiently low value of t (with some approximation):

W d φ (p, q t,F ) = W L 2 (p φ , q φ,t,F ) = max C∈Lip 1 E x∼p C φ(x) -C φ(x + t × F (x) (117)
In fine, we can forget about the limit operator in Eq. ( 116) for a low value of t because the optimal transport is the natural transport even in the Kantorovich-Rubinstein formulation and thus we have:

max θ F ,θ φ ,θ C L(θ F , θ φ , θ C ) (118) with L(θ F , θ φ , θ C ) = E x∼p C φ(x) -C φ(x + t × F (x) t (119)
and in a stochastic gradient descent strategy, the only important quantity is a bias-free estimate of the gradient [START_REF] Robbins | A stochastic approximation method[END_REF]:

∇L(θ F , θ φ , θ C ) = ∇ ∑ B b=1 C φ(x i b ) -C φ(x i b + t × F (x i b ) t × B ( 120
)
for a small t and some B random indices

i b ∼ U N (1, N).
In conclusion of this section, the optimization problem finally gets:

max θ F ,θ φ ,θ C E x∼p C φ(x) -C φ(x + t × F (x) t (121)
such that φ is bijective, C is 1-Lipschitz and F Jacobian the BA and ZGLA properties in Eq. ( 96) and Eq. ( 97) for small t given by Eq. ( 112). For the specific computer vision scientific field, pixels are made of one or three coordinates (for gray levels or color images respectively), so showing relevant pixels in images of the same category could lead to an unsupervised foreground/background segmentation where the only supervision is the fact that all images belong to the same semantic category. A similar problem was tackled by Joulin et al.

[2010] a few years ago but with a little more supervision: we have access to several images categories labels and they call the problem cosegmentation. As they stated:

Purely bottom-up, unsupervised segmentation of a single image into foreground and background regions remains a challenging task for computer vision.

This remains true for all kind of data but without different meanings. Fundamentally, if one has a dataset, one could interpret relevance measurements this kind of algorithms could provide. In computer vision, the core idea of co-segmentation is that the availability of multiple images that contain instances of the same "object" classes makes up for the absence of detailed supervisory information. Some research has been efficiently conducted for interactive foreground / background segementation [START_REF] Rother | grabcut": interactive foreground extraction using iterated graph cuts[END_REF]] but here we would want to avoid user interaction and benefit from a whole dataset: a class of data sharing a common pattern that we want to highlight. The only form of supervision is knowing that data share some information of interest without knowing what precisely.

F U T U R E W O R K A N D C O N C L U S I O N

Investigating Wasserstein distances with varying metric seems promising for future work. This sketch of contribution can be a stepstone answer to metric invariance pointed out by [START_REF] Kleinberg | An impossibility theorem for clustering[END_REF] for clustering (which is an unsupervised task like feature importance extraction in this work). This shows that this work can be improved in terms of machine learning and optimization and also engineering on real world data.

Indeed, computer vision in general and foreground/background unsupervised segmentation (in the way we present it) in particular are ways to provide a better understanding between unsupervised metric learning and feature importance extraction thanks to large cardinality datasets.

Les espoirs que suscite l'informatique quantique annoncent des améliorations révolutionnaires en matière d'intelligence artificielle (IA) embarquée dans les produits du quotidien. Bien que la recherche en IA évolue rapidement vers l'informatique quantique avec des personnes clés comme Pr. Julia Kempe (succédant à un pionnier de l'IA, le Pr. Yann LeCun à la tête du Centre de Data Science de l'Université de New York depuis quelques années).

Néanmoins, des changements plus imminents et plus radicaux sont déjà en cours grâce à l'électronique de pointe : par exemple, l'accélération matérielle amplifiée par les "processeurs neuronaux" d'Amazon et d'Apple (c'est-à-dire des processeurs dédiés contrairement aux processeurs traditionnels "CPU" généralistes). En particulier, l'électronique de faible puissance est cruciale, car la nonstockage de l'électricité reste un problème malgré les efforts de recherche d'Apple, de Tesla (et de nombreux autres).

A B S T R A C T Uncertainty estimation with neural networks predictions is not yet well studied (except in the recent and controversial field called Bayesian deep learning). Moreover, in the academic and business worlds, some misinterpretations persist about the probabilities as outputs of a supervised classifier: a probability of belonging to a category or its associated vector across all categories does not fully express information on uncertainty. Indeed, most practitioners extract the maximum entry of a probability vector, which is hacky and not sufficient for a more thorough interpretation of uncertainty. At the same time, in many sensitive applications where security, health or even justice issues are involved, it seems that uncertainty estimation is crucial.

The core idea of our Hypothesis of an Uncertainty Model (HUM) contribution is predicting parameters of an output law instead of predicting an output estimate. On the one hand, putting a parametrized probabilistic law on the output side makes our approach Bayesian. On the other hand, using neural networks to produce these parameters makes it less principled but more pragmatic. By combining the best of these two worlds in a small number of additional lines of code is convenient in practice to estimate uncertainty for automatic or assisted decisions in supervised learning in an engineeringly feasible but yet theoretically principled fashion. There is a large class of supervised learning problems that can be cast in an optimization of this form:

min θ F L(θ F ) (122) L(θ F ) = E (x,y)∼Nature ℓ y, F (x)
where the quantity ℓ(y, F (x)) expresses how much one tolerates error for confusing the prediction F (x) with the groundtruth output y for input x. Our notations abide by the deliberate starting epistemological choices of this dissertation in section 1.2.2 especially for the Nature idealized distribution. The prediction function F is implemented here by a neural network (whose parameters are θ F ) but we insist on the fact that it could be implemented by other tools that we can pick from the enormous machine learning zoo of supervised learning of algorithms [START_REF] Bishop | Pattern recognition[END_REF] while still keeping this current contribution relevant. We will describe how to modify this Eq. ( 122) and how to create the right optimization scheme while still being both mathematically principled and engineeringly friendly in a general case with many applications. But first, we want to debunk some widespread misunderstanding among practitioners in supervised classification about how to interpret an output classification probability vector with respect to uncertainty (which was our engineering starting point). Typically, for an image classification problem (say with an AlexNet convolutional neural network [START_REF] Krizhevsky | Imagenet classification with deep convolutional neural networks[END_REF] trained on ImageNet [START_REF] Fei-Fei | Imagenet: crowdsourcing, benchmarking and other cool things[END_REF]), the output is a 1-sum positive vector of high dimensions (for K ≃ 10 3 categories for Ima-geNet) in the form of conditional probability estimates (P(c = k|x)) k=1,...,K for an image x belonging to which category k ∈ 1, K . Unfortunately, there is a detail that is often overlooked concerning to the "knowing part |x" of the notation P(c = k|x) which states in an implicit (and sometimes forgotten) fashion that the input x should be taken from the same theoretical distribution both at training and testing times. Indeed, at training time, the formulas are by construction true because the images are taken from Nature (better approximated by large cardinality dataset). At test time, the input examples must be coming from the same Nature distribution for the knowing part of the probability being true which is just a reasonable hypothesis in theory but that is difficult (and rather impossible) to guarantee in practice (i. e. in every real-world and practical industrial applications). The best proof of it is that if we give a completely random input x taken from an arbitrary distribution, the system would still answer a 1-sum positive vector output although the input does not belong to any of the categories corresponding to these output vector entries due to its independently random nature.

From these considerations, we draw 2 conclusions for each input x at test time: (i) taking the maximum entry index k * of such outputs (P(c = k|x)) k=1,...,K has indeed some significance in a maximum a posteriori multinomial perspective in order to make a concrete decision but nevertheless (ii) the corresponding maximum entry estimating P(c = k * |x) has no or little uncertainty significance. If uncertainty is a matter of interest, then the output should be considered as a whole K-dimensional point (and not just the maximum of its entries) living on the unit simplex and uncertainty should be a possibly blurry zone around one predicted point in that space. Now that we have described some practical issues related to uncertainty in classification, we present a general solution to tackle supervised learning problems that takes into account some uncertainty information.

The spectacular come back of deep learning techniques can be dated back to 2012, when [START_REF] Krizhevsky | Imagenet classification with deep convolutional neural networks[END_REF] dramatically improved image classification scores. The suprise was considerable because one particular reason (among others): beyond larger experimental scales, there was almost no conceptual difference between the initial LeNet [LeCun and Bengio, 1995] in 1995 and that implementation except: the activation function (from the Sigmoid function to the positive part function) for marginal accuracy improvements, impressive implementation improvements with the original usage of GPU hardware (Graphics Processing Units) for speed (otherwise, training are infeasible for this kind of year-long-type optimization with classic CPUs of traditional computers at that time) and DropOut regularization [START_REF] Srivastava | Dropout: a simple way to prevent neural networks from overfitting[END_REF] to cope with overparametrization of larger neural networks than usual. Indeed, DropOut consists in injecting artificial random noise to learning paramaters during training with the intuition that without that noise, testing accuracy will be more robustly improved for generalization purposes.

The originality of what [START_REF]Uncertainty in Deep Learning[END_REF] proposed consists in continuing DropOut at test time which implicitly provides (dependent) several outputs for each input and then estimating law parameters on those (say mean and variance of a Gaussian law). Averaging correlated outputs (the learned parameters are very close) has unclear theoretical implications for the Monte Carlo estimation but this method has the merit of pointing to a practical need to investigate sound statistical approaches to cope with the black-box deserved reputation of neural networks (outputs without "as-is" possible interpretations).

In the past, [START_REF] Nix | Estimating the mean and variance of the target probability distribution[END_REF] proposed to predict the mean and covariance of each prediction by minimizing the Kullback-Leibler divergence between the output empirical data distribution and the Gaussian output prediction (i. e. maximizing the likelihood). It is based on the assumption that there is a sufficiently large data set, i. e. , that their is no risk of overfitting and that the neural network finds the correct regression.

In our work, we chose to suggest a lighter machinery (than the ones of [START_REF] Nix | Estimating the mean and variance of the target probability distribution[END_REF] and [START_REF]Uncertainty in Deep Learning[END_REF]): instead of changing learning parameters to implicitly produce an ensemble effect on the prediction side, we directly estimate the parameters of a law on the prediction side. We removed the randomness on the weights that the approach of [START_REF]Uncertainty in Deep Learning[END_REF] implies to put it entirely on the prediction side giving a beneficial out of local minimum escape effect to our system while still keeping a valid Kullback-Leibler interpretation like the approach of [START_REF] Nix | Estimating the mean and variance of the target probability distribution[END_REF] and even [START_REF] Bishop | Mixture density networks[END_REF]. We highly recommend careful readers about uncertainty in deep learning to see the video of Zoubin Ghahramni1 to get an hindsighful overview of that kind of literature.

As recently describes by [START_REF] Detlefsen | Reliable training and estimation of variance networks[END_REF], estimating uncertainty is desirable goal in several contexts:

• for time series with Gaussian processes [START_REF] Seeger | Gaussian processes for machine learning[END_REF] where means and variances of predictions are estimated at each timestamps;

• the historical work of [START_REF] Bishop | Mixture density networks[END_REF] and [START_REF] Nix | Estimating the mean and variance of the target probability distribution[END_REF] using the seminal idea of predicting the parameters of a law instead of a value that we build upon;

• meanwhile the probabilistic framework of Bayesian techniques was also useful thanks to [START_REF] Mackay | A practical bayesian framework for backpropagation networks[END_REF] and revisited by [START_REF] Kingma | Auto-encoding variational bayes[END_REF];

• Monte Carlo Dropout of [START_REF]Uncertainty in Deep Learning[END_REF] with unclear theoretical consequences as predictors are correlated;

• active learning while uncertainty prioritizes choice of informative training examples [START_REF] Huang | Active Learning by Querying Informative and Representative Examples[END_REF].

In this work, we aim to generalize the inspiration of [START_REF] Nix | Estimating the mean and variance of the target probability distribution[END_REF] for a larger part of supervised deep learning problems.

R E V I S I T I N G D E E P S U P E R V I S E D L E A R N I N G

As described in the state of the art chapter at section 1.2.3, classification and regression are two important fields applications of supervised learning and are thus both applications of this current uncertainty estimation contribution.

Classification

Neural Networks are a fundamentally continuous technology that progressively train algorithms by taking into account its own errors on mini-batches of data in order to improve thanks to stochastic gradient descent and back-propagation throughout all the parameters [START_REF] Goodfellow | Deep Learning[END_REF].

In fact, even plain classification with neural networks is a special case of regression where a degenerated histograms encode categories (a clean category corresponding to a vector with zero entries except for one entry one, hence the name one-hot encoding) which has the advantage of naturally manipulating probabilities with float real numbers and of letting us access to the continuous regression tools. For potentially K different categories, the supervised classification literature tends to encourage the negative log-likelihood loss corresponding to:

ℓ(y, z) = - 1 K y ⊤ log(z) = - 1 K K ∑ k=1 y (k) × log(z (k) ) (123)
where the label y adopts the one-hot encoding which gives naturally access to its relaxed probabilistic interpretations. Depending on the scientific culture of the reader, Eq. ( 123) can also be seen as the cross-entropy loss or even the Kullback-Leibler divergence between the two discrete densities represented in 1-sum positive vectors z and y. Indeed, the empirical likelihood of given predictions z = F (x) measured on (input, output) pairs (x, y) from training data is

Likelihood = N ∏ i=1 K ∏ k=1 (z (k) i ) y (k) i 1 NK (124) = N ∏ i=1 K ∏ k=1 (F (x i ) (k) ) y (k) i 1 NK
and the empirical negative log-likelihood gives:

-log(Likelihood) = -1 NK N ∑ i=1 K ∑ k=1 y (k) i × log(z (k) i ) (125) = -1 NK N ∑ i=1 K ∑ k=1 y (k) i × log(F (x i ) (k) )
which justifies Eq. ( 123).

At the same time, if we put the Kullback-Leibler divergence KL(y i , z i ) definition between two discrete densities y i and z i = F (x i ):

KL(y i , z i ) = 1 K K ∑ k=1 y (k) i × log y (k) i z (k) i = 1 K K ∑ k=1 y (k) i log(y (k) i ) -log(Likelihood) (126) but the first term 1 K ∑ K k=1 y (k) i log(y (k)
i ) does not depend on the input x i nor the prediction z i which makes it removable towards a predictor optimization scheme which also justifies the minimization of Eq. ( 123) in this classification memory aid:

Maximizing the classification likelihood ⇐⇒ Minimizing the classification cross-entropy loss ⇐⇒ Minimizing the mean Kullback-Leibler divergence between groundtruth labels and predictions discrete densities (127)

Regression

For regression, y is continuous (floats) and possibly multivariate, e. g. the price of rent corresponding to each location in a city, supply chain forecasts, weather forecast maps. A stable and well-studied optimization tool is the quadratic regression:

ℓ(y, z) = y -z 2 2 = K ∑ k=1 y (k) -z (k) 2 (128)
another loss is the Manhattan one:

ℓ(y, z) = y -z 1 = K ∑ k=1 |y (k) -z (k) | (129)
to focus on small errors (as big errors are considerably more punished in quadratic rather then Manhattan loss ).

In a way that is similar to establishing a likelihood for cross-entropy classification, it is possible to draw the links between regression and likelihood. As a matter of fact, the likelihood of a Gaussian model on the output side centered on the predictions F (x) with an unknown but fixed diagonal homothety covariance matrix σI K (for the sake of understanding and lack of any sophisticated hypothesis, it turns out we do not loose generality with σ = 1 in our particular case):

Likelihood =   N ∏ i=1 exp -1 2 y -F (x) 2 2 (2π) K 2   1 N (130)
as the neg-log-likelihood gets:

-

log(Likelihood) ∝ 1 N N ∑ i=1 y i -F (x i ) 2 2 (131)
and likewise for a the Kullback-Leibler divergence between that Gaussian density q centered at F (x) (with identity covariance without loss of generality) and the distribution p of ground-truth labels y | x (when (x, y) ∼ Nature), we eventually have:

KL(p, q) = E y∼q log q(y) p(y) = E y∼q log (q(y)) -E y∼q log (p(y)) (132) (133)
but the first term is once again useless for optimizing F , which justifies the least squares approach because:

-

E y∼q log (p(y)) = -E (x,y)∼Nature   log exp -1 2 y -F (x) 2 2 (2π) K 2   ∝ E (x,y)∼Nature y -F (x) 2 2
In the end we re-establish a regression memory aid:

Maximizing the regression likelihood (134)

⇐⇒ Minimizing the regression quadratic loss ⇐⇒ Minimizing the mean Kullback-Leibler divergence between the groundtruth outputs distribution and the gaussianized predictions distribution Please note that one could do exactly the same for Laplacians and Manhattan L 1 loss instead of Gaussians and quadratic (eucidlean) L 2 criterion.

H U M : H

Y P O T H E S I S F O R A N U N C E R TA I N T Y M O D E L
The aim of this work is to propose a way to introduce some uncertainty modelization in supervised machine learning presented above focusing on deep neural networks for the sake of recent research attention. We want to build a model that provides a probabilistic law prediction instead of a value prediction which is following a Bayesian tradition. At test time, we would use for example the main mode of our output law as the output value (i. e. in lieu of a value prediction) but with some additional uncertainty information depending e. g. on how big the support of the distribution is around that mode.

To put this in simple layman's terms and with some abuse of terminology, traditional approaches try to predict a specific value F (x) hoping to match the label y in order to ideally reach even at test time some sort of conformity y ≃ F (x)

(approximately) as the loss function evaluation ℓ(y, F (x)) role is to avoid too much discrepancy during training between prediction F (x) and groundtruth label y.

In contrast, we propose to focus on uncertainty thanks to outputs corresponding to interpretable parameters of a smooth distribution G(x) hoping that we reach another kind of conformity

y ≈ G(x) (136)
(approximately). Indeed, the graphic slight difference of notation between the symbols "≈" and "∼" inspired us. More mathematically, this work relies on two statements. First, thanks to the definition of Dirac distributions, we easily establish that we always have:

(∀(x, y) ∈ R D × R K ) ℓ y, F (x) = R K ℓ y, z × δ F (x) (z)dz (137) = E z∼δ F (x) ℓ y, z
where δ a is the Dirac distribution located at a that gives that same a value for when integrated throughout space. Second, as recalled in section 1.5.1 with the Robbins-Monro theorem, we do not need to have access at precisely the exact evaluation of the minimization objective or of its gradient in order to minimize it by stochastic gradient descent: only a biased-free estimator of the gradient is required.

The core idea of this chapter is to ask: Why not taking a more convenient and more interpretable distribution than the Dirac distribution? in Eq. ( 137). This intuition is motivated by the fact that using sum Dirac distributions is considering Nature as just a set of scattered points whereas a smooth distribution gives some consistency to our so-called Nature modelization. Now, if we are emancipated from the previous uncertain-free formulation of section 4.2, then in terms of modelling combining Eq. ( 137) and Eq. ( 122) gets:

min G E (x,y)∼Nature E z∼G(x) ℓ y, z (138)
and G(x) is a prediction law whereas F (x) previously was a prediction value.

A discretized version of this approach has been extensively applied to supply chain forecasts by the Lokad company [START_REF] Vermorel | Probabilistic forecasting[END_REF] since few years ago. In this chapter, dealing with a wider range of applications (especially suitable with neural networks as we will see) in a more continuous fashion is our contribution to the original and seminal idea of [START_REF] Vermorel | Probabilistic forecasting[END_REF].

In terms of optimization of the parameters θ G of function G predicting the output law G(x), the good news is that the same Robbins-Monro theorem [START_REF] Robbins | A stochastic approximation method[END_REF] as in section 1.5.1 justifies the usual Hebbian learning rule of stochastic gradient descent. For the sake of clarity we admit that what is true to the standard stochastic gradient descent algorithms is also true for modern stochastic algorithmical routines like Adam [START_REF] Kingma | A method for stochastic optimization[END_REF] which benefits from the rest of regular contemporary deep learning tools [START_REF] Abadi | Ten-sorFlow: Large-scale machine learning on heterogeneous systems[END_REF][START_REF] Paszke | K. Pearson. Liii. on lines and planes of closest fit to systems of points in space[END_REF].

For the usual uncertainty-free settings, we usually had:

θ t+1 F = θ t F -α t ft (139)
where α t is the learning rate and ft = 1 In our proposed uncertainty modelling, we get:

B ∑ B b=1 ∇ θ F ℓ y i b , F (x i b ) with i b ∼ U N (1, N)
θ t+1 G = θ t G -α t ĝt (140)
with same kind of learning rate α t as before. The main difference is in the bias-free gradient estimation that needs a little Monte Carlo estimation (as M stands for Monte Carlo number of iterations):

ĝt = 1 MB M ∑ m=1 B ∑ b=1 ∇ θ G ℓ y i b , z m,i b (141)
where we still have i b ∼ U N (1, N) the same kind of random uniform index as before but we also have M sampled outputs from the same prediction law for each input x i b totalling MB gradient evaluations per mini-batch: z m,i b ∼ G(x i b ) which is M times more computations than before but there is no guarantee confirming or infirming that a big B or a big M could improve convergence or results as the two sources of stochasticity are related (maybe even a simple M = B = 1 scenario could give good resuls).

There is a strong link between our technique and the Reparametrization Trick [START_REF] Kingma | Auto-encoding variational bayes[END_REF] for initially variational auto-encoders problematics. Indeed, for the authors, several probabilistic laws, the law parameters can be represented in the gradient calculus from its normalized version (e. g. without mean nor unusual standard deviation for the Gaussian law) of the same simplified probabilistic law accordingly distorted to fit the desired parametrized law. Thanks to three special cases, we will now precisely see how this is all working thanks to pseudo-random generators of our computers.

Uncertain Logistic Regression for Classification

In logistic regression classification, we recall from section 4.2.1 that for the general optimization min F E (x,y)∼Nature [ℓ(y, F (x))] is specified for logistic classification by:

ℓ(y, z) = - 1 K y ⊤ log(z) (142)
and to maintain the simplex-constraint we propose a uniform law (as segments always remain inside a convex set) for applying Eq. ( 141):

z m,i b ∼ U R K A(x i b ), B(x i b ) (143) ⇐⇒ t m,i b ∼ U R (0, 1) and z m,i b = t m,i b × A(x i b ) + (1 -t m,i b ) × B(x i b )
where functions A and B give class membership probabilities and are implemented with regular neural networks with a last SoftMax layer and it all becomes:

ĝt = 1 MB M ∑ m=1 B ∑ b=1 ∇ θ G ℓ y i b , t m,i b × A(x i b ) + (1 -t m,i b ) × B(x i b ) (144) = 1 MB M ∑ m=1 B ∑ b=1 ∇ θ G -y ⊤ i b log t m,i b × A(x i b ) + (1 -t m,i b ) × B(x i b )
to be compared to Eq. ( 141). In these settings, the prediction uncertainty can be measured in 0 ≤ 1 2 A(x i b ) -B(x i b ) 1 ≤ 1 which justifies our probabilistic model. In terms of Kullback- Leibler divergence, minimizing the mean divergence between the fixed discrete groundtruth density y and the random discrete prediction density z ∼ G(x) this way, is useful even at a category level by superposing both A(x) and B(x) histograms (i. e. positive 1-sum vectors).

Out of curiosity, we can explore what is happening for more than 2 histograms for the sake of algorithms.

z m,i b ∼ U R K A 1 (x i b ), . . . , A c (x i b ), . . . , A C (x i b ) (145) ⇐⇒ t m,i b ∼ ∆ C and z m,i b = C ∑ c=1 t (c) m,i b × A c (x i b )
where ∆ C is the uniform distribution over the C-order simplex (i. e. ∑ C c=1 t (c) m,i b = 1). To uniformly sample inside that C-order simplex, first, we first sample uniformly C -1 values between 0 and 1 u m,i b ,c ∼ U R (0, 1) (and we rearranged the indices corresponding to c such that they are re-ordered in increasing order), second we recursively build:

t (1) m,i b = u m,i b ,1 (146) t (c) m,i b = u m,i b ,c+1 -u m,i b ,c for c ∈ 2, C -1 t (C) m,i b = 1 -u m,i b ,C-1
that cancels out in a sum of 1 and finally we get:

ĝt = 1 MB M ∑ m=1 B ∑ b=1 ∇ θ G ℓ y i b , C ∑ c=1 t (c) m,i b × A c (x i b ) (147)
which implies an output law that fires in a convex set drawn by the polygon made of (A c (x)) c=1,...,C for each input x. Although this formulation is theoretically pleasing, it did not provide, in our experience any additional interpretation ease nor accuracy improvements, yet.

Uncertain Least Square for Regression

Likewise, in least square regression:

ℓ(y, z) = y -z 2 2 ( 148
)
and we choose the Gaussian distribution to obtain:

z m,i b ∼ N µ(x i b ), C(x i b ) × C(x i b ) ⊤ (149) ⇐⇒ ǫ m,i b ∼ N (0, I K ) and z m,i b = µ(x i b ) + C(x i b ) × ǫ m,i b
where function µ is a vanilla neural network and C is a lower-triangular matrix (with possible bounded diagonal values thanks to affine Sigmoids see Table 1 at page 14). For the Hebbian rule we update with:

ĝt = 1 BM B ∑ b=1 M ∑ m=1 ∇ θ G y i -µ(x i b ) + C(x i b ) × ǫ m,i b 2 2
(150)

In these settings, the uncertainty information is contained in the C(x i b ) × C(x i b ) ⊤ covariance matrix as we have fitted a Gaussian to our prediction outputs. Obviously, for the sake of completeness, we have modelled a Gaussian distribution with a full covariance matrix but a diagonal covariance matrix or a simple proportional to the identity covariance matrix is most of the times enough in real-world applications.

In terms of Kullback-Leibler divergence, we still have the same probabilistic interpretations as we minimize the mean over input x of the divergence between the deterministic groundtruth labels y Dirac distribution and the smooth Gaussian distribution N µ(x), C(x) × C(x) ⊤ which means we did manage not to sacrifice interpretation for uncertainty. The mean µ(x) plays the role of the good old-fashioned value prediction that is tainted by incertainty C(x) × C(x) ⊤ .

Uncertain Mixtures for Regression and Classification

Sometimes, beyond mono-modal distribution for regression (e. g. often Gaussian) and uniform law for classification, one can imagine a non-continuous spectrum of outcomes like in supply chain forecasts [START_REF] Vermorel | Probabilistic forecasting[END_REF] to handle several discrete scenarios at the same time. This way, one can benefit from the computational power of machines to process newly predictable outcomes. This is the reason why for each input x we tackle here one finite mixture q(x) of C distributions with proportions π(x) and components (g k (x)) K k=1 (that are themselves distributions) which can written as:

q(x) = C ∑ c=1 π c (x) × g c (x) (151)
Just by the continuous integral definition of the expectation, we obtain:

E z∼q (ℓ(y, z)) = C ∑ c=1 π c (x) × E z c ∼g c (x) [ℓ(y, z c )] (152)
which facilitates our mixture use because we can directly use the formulas of the previous sections 4.3.1 and 4.3.2 depending on the applications we have:

ĝt = 1 MB M ∑ m=1 B ∑ b=1 C ∑ c=1 ∇ θ G π c (x) × ℓ y i b , z m,i b ,c with z m,i b ,c ∼ G c (x i b ) (153)
where the previous technique of reparametrization as previously can re-applied to generate samples from G c (x i b ).

In this mixture scenario, the parameters set θ G contains the parameters θ π of the neural network π that gives the C proportions and each of the parameters θ G c s that parametrize the laws G c (x)s. Two special cases of mixtures attract our attention:

G A U S S I A N M I X T U R E F O R R E G R E S S I O N
Interestingly, this corresponds to the Mixture Density Networks designed by [START_REF] Bishop | Mixture density networks[END_REF] in the 1990s; Interestingly, our straightforward yet efficient approach has some connections with other techniques. Let us take into account that the discrepancy function ℓ in Eq. ( 138), if we recall that in classification and regression cases, we chose ℓ(y, z) = -y ⊤ log(z) and ℓ(y, z) = yz 2 2 respectively. In both cases, for any label y, we know that the function z → ℓ(y, z) is convex and thus, thanks to Jensen inequality for any distribution q of predicted output (e. g. q = G(x)):

D I R A C M I X T U R E F O R C L A S S
E z∼q (ℓ(y, z)) ≥ ℓ y, E z∼q (z) (154)
which implies that:

E (x,y)∼Nature E z∼G(x) ℓ y, z ≥ E (x,y)∼Nature ℓ y, E z∼G(x) z (155)
which means that for a non-Dirac global minimum G * of Eq. ( 138) corresponding to the lower bounded right hand side of Eq. ( 155), we know that the global minimum is also reached at

G d : z → E z∼G * (x) (z) (156)
One possible interpretation of that could be The deterministic optimization and our optimization share some global minima if they exist. What is interesting to mention is that in practice, first, space zones that is very close to the training data manifold have roughly the same result for with or without uncertainty model and second, in-between zoom have better results and interpretation with out uncertainty model. Another way to look at it is to consider our technique like a DropOut [START_REF] Srivastava | Dropout: a simple way to prevent neural networks from overfitting[END_REF] technique consisting in adding some randomness to avoid being stuck in local minima. Our approach also looks like the one of ? but we do linear combination on the output instead of doing it at the input like them. Our method needs further investigation beyond showing good results because of intringuing self-regularized optimization phenomena thanks to randomness. Although, it is rigorous to separate model from optimization, in our case, it is possible that the additional randomness of our optimization technique actual makes our model more robust and it is unclear to attribute improvements between model and optimization as we will see in the next section. Based on the updated Hebbian rules, we previously described in normal and out uncertain contexts, the stochastic gradient descent algorithm optimizes the parameters of our neural networks. Using the (pseudo)-random number generators of our computers, we were able to generalize the Robbins-Monro theorem allowing us to calculate our bias-free estimators of gradients at each stochastic gradient iteration.

Safe Computations

Some numerical overflow instabilities have been studied and avoided thanks to the so-called logsumexp trick. This technique must be revisited in our case.

Let's first briefly describe what it usually is and then how we can adapt it for our implementations.

Even though the logsumexp trick can be used in many scenarios, we choose to describe it in the special example of the previouly presented logistic regression of classification for the sake of pragmatism. In neural networks usual implementations, the cross-entropy loss is applied after the computation of a SoftMax layer on what we call logits c = [c 1 , . . . , c k , . . . , c K ] (see the previously described in Table 1 at page 14). At some point, there is a need to compute the gradient of the quantity

v ∈ R K : v = log(SoftMax(c)) (157)
with respect to the neural network parameters involved in the calculus of c. If we take a closer look to v, we expand it to:

v k = log exp(c k ) ∑ K k ′ =1 exp(c k ′ ) = log(exp(c k )) -log K ∑ k ′ =1 exp(c k ′ ) (158)
The direct computation of exp(c k ) where c k ≥ 30 is already bigger than 10 13 which makes it all unfeasible for a modern computer even though we know that log(exp(c k )) ≃ log(10 13 ) ≃ 30 is a reasonable number. For any value M ∈ R, we can use the fact that:

v k = M + log(exp(c k -M)) -M -log K ∑ k ′ =1 exp(c k ′ -M) = log(exp(c k -M)) -log K ∑ k ′ =1 exp(c k ′ -M) (159)
and for the specific choice M = max k c k , our numerical problems disappear because all exponential arguments c k -M ≤ 0 are more suitable for the computation of exp(c k -M) ≤ 1 without approximation. Nowadays deep learning tools such as PyTorch [START_REF] Paszke | K. Pearson. Liii. on lines and planes of closest fit to systems of points in space[END_REF], Tensorflow [START_REF] Abadi | Ten-sorFlow: Large-scale machine learning on heterogeneous systems[END_REF] and MXNet [START_REF] Chen | Mxnet: A flexible and efficient machine learning library for heterogeneous distributed systems[END_REF] applies that logsumexp trick implicitly behind the scene but this can reused in this present work. We can remark that the idea consisted in avoiding big exponentials by guaranteeing negativity of their arguments thanks to the maximum entry that is substracted for safer computations.

Likewise, in our uncertainty estimation context, we have to compute at some point the quantity w ∈ R K for t ∈ [0, 1]:

w = log (t × SoftMax(a) + (1 -t) × SoftMax(b)) ( 160
)
where a ∈ R K and b ∈ R K come from neural networks. With that negative argument exponentials technique in mind, we understand that the SoftMax function is invariant to uniform shifting:

(∀M ∈ R) SoftMax(c) = SoftMax(c -M) (161)
thus we propose for any P ∈ R and Q ∈ R:

w k = log t × exp(a k ) ∑ K k ′ =1 exp(a k ′ ) + (1 -t) × exp(b k ) ∑ K k ′ =1 exp(b k ′ ) = log t × exp(a k -P) ∑ K k ′ =1 exp(a k ′ -P) + (1 -t) × exp(b k -Q) ∑ K k ′ =1 exp(b k ′ -Q)
and to guarantee the same negative exponentials properties as before, we choose: P = max k a k and Q = max k b k which provides much more numerical stability when it comes to computing w and its gradient with respect to the parameters a fortiori.

In the special case where the number of classes K = 2, practioners use a Sigmoid function applied on a logit value c and the 2 estimated probabilities get [Sigmoid(c), 1 -Sigmoid(c)] but the SoftMax function is a generalization of the Sigmoid for that K = 2 binary case thanks to that relationship:

SoftMax ([0, -c]) = [Sigmoid(c), Sigmoid(-c)] (162) = [Sigmoid(c), 1 -Sigmoid(c)]
which gives us Eq. ( 162) in that K = 2 binary special case for no supplementary effort.

Re-Using Uncertainty-Free Models

In 2018, Apple and Amazon reached $10 12 of market capitalization and in 2019 Microsoft followed in 2020 by Alphabet that owns Google also did reach that financial milestone while Facebook is at roughly $0.6 × 10 12 . So we humbly imagine a way to use all these powerful open-source means to better train our uncertainty models just to stand on the shoulders of giants.

In order to re-use all that wealth of trained models, we use the idea that an uncertainty-free model is a model that has no uncertainty: for example, a Gaussian density with no covariance becomes a Dirac distribution, a uniform density with equal bounds also becomes a Dirac distribution. We can briefly explain this strategy in two supervised instances: A = B ⇐⇒ no uncertainty (164)

In our experiments, approximately 10% more of training epochs is enough to get relevant additional uncertainty information. We must insist on the fact that we did not want to improve accuracy of already impressively sophisticated algorithms (at the very least with an engineering perspective and for means of distributed computational power) but only we only wanted to get some new uncertainty information.

Algorithms

There exists many ways to implement the neural networks from our HUM (Hypothesis of Uncertainty Model). In this section, we provide a binary classification pseudo-code and another pseudo-code for regression.

Neural Networks are often built in layers which correspond to a composition of elementary functions (even in deep residual learning [START_REF] He | Deep residual learning for image recognition[END_REF], we can still separate blocks of layers). This is good news for cheap learning of heavy neural networks with uncertainty because they can be initialized with worldwide companies means but without uncertainty. For example, in classification, a neural network corresponding to a function F can be thus separated it into F = G • H In order to go deeper with our approach, we study a toy dataset called "Two Moons" with 1000 points for each of the 2 groups (blue and red 2D dots) in 2 dimensions as presented in Fig. 21 in two scenarios: an easily separable one and another much noisier one. In this supervised classification exercize, the easy case wihtout noise in not very interesting from an uncertainty perspective which makes our approach unnecessary on the left. On the contrary, the difficult case is more relevant because some space regions imply legitimate doubt levels on the right.

We trained a small multi-layer perceptron on this dataset to investigate what was empirically at stake: for each input [x 1 , x 2 ] ⊤ we predict 2 probability vectors p = [a, b] ⊤ and q = [c, d] ⊤ as bounds of a uniform prediction law ŷ ∼ U R (p, q) following section 4.3.1 (page 99). Thus, it is fair to consider u = p-q 2 ∈ [0, 1] as an evaluation of uncertainty. Our goal is not to provide any classification accuracy improvements but only some more additional information about uncertainty. So we empirically checked on this synthetic problem if our technique was able to preserve good classification results which is the case in Fig. 22 with yellow for the blue class, a purple for the red class and small region of green for the blurry boundary with a thickness that is understandable when considering the mixed classes. Now we can see if it is also providing a useful uncertainty information in Fig. 23.

Here we see an interesting phenomenon: the prediction is highly confident almost everywhere in space (purple low level of uncertainty) except in the boundary which can be split in two: (i) low uncertainty boundary in purple at the center of the figure and high uncertainty far from where the points are (whereas the boundary is green in the previous Fig. 22 wihtout distinction of having a 0.5 probability of being blue.)

These encouraging synthetic results exhibits two ways of not knowing: with uncertainty in yellow in Fig. 23 and with certainty in purple in Fig. 23. The Dirac phenomenon we suspected in section 4.3.4 is occuring: close to the training manifold, the uniform distribution we get is a Dirac but going further from the training manifold grows some thickness away from the Dirac distribution. (the higher, the less certain). Blue became white for readibility reasons.

Dogs and Wolves

To illustrate our approach, we took a image classification task that is difficult even for us as human beings: distinguishing domesticated canide animals (say mainly dogs) and wild canide animals (say wolves, foxes etc.) that we took from the ImageNet dataset [START_REF] Fei-Fei | Imagenet: crowdsourcing, benchmarking and other cool things[END_REF] (and its tree-structured label ontology). We formed a dataset2 of a balanced 800 images per class training set and 200 images per class validation set on which we measure accuracy. We also formed a test set from the ImageNet test dataset (annotated by the ResNet-152 available via PyTorch) with 100 cats, 100 dogs and 100 wolves to evaluate test error and extrapolation uncertainty (i. e. out of distribution) as cats are not wolves nor dogs3 . As [START_REF] Greenspan | Guest editorial deep learning in medical imaging: Overview and future promise of an exciting new technique[END_REF] explains, using a convolution neural networks (CNN) trained from a natural images classification problem with plenty of labeled data (e. g. ImageNet [START_REF] Fei-Fei | Imagenet: crowdsourcing, benchmarking and other cool things[END_REF]) is useful for a medical images classification problem that suffered from a notorious lack of large scale labeled datasets until recently. Indeed, chopping the head off a CNN is relevant as a feature extractor even for a very different problem which is understandable because early computer vision cues such as the relationship between contours and gradients in the first layers seems to be independent from the image domain and thus shared among all computer vision tasks. Meanwhile, some empirical studies tend to show that the first layers of a convolutional network often correspond to Gabor functions for general-purpose computer vision applications and the deeper layers can be interpreted in more sophisticated ways [START_REF] Zeiler | Visualizing and understanding convolutional networks[END_REF].

Inspired by these successful ideas and intuitions, first we build a convolutional neural network like in Fig. 25 close to AlexNet [START_REF] Krizhevsky | Imagenet classification with deep convolutional neural networks[END_REF] but we double the number of ending layers and second we make all the first layers shared for two separate ending dense neural networks. Indeed for our "Dogs and Wolves" images classification problem, we need two outputs with the belief that some features are shared among these virtually two neural networks. We do not entirely double the number of parameters with two separate convolutional neural networks because we want to prevent our learning procedure from some bad overfitting effects of over-parametrization and gain some more speed. Our "Dogs and Wolves" images classification problem looks like our "Two Moons" noisy case with much higher dimensionality and where even human beings might disagree sometimes. In fact, the very notion of being domesticated is human-related for example. Numerically, we observe a slight improvement of less than 1% of validation accuracy from a 76% in normal uncertainty-free training and 77% with our HUM technique. In fact, we initialized as previously stated in section 4.4.2 our HUM convolutional neural network with a classic uncertaintyfree regular approach. The goal here was not to improve the validation score but only to check if we are not loosing some which is valid as Fig. 26 is showing little overlaid differences. The real benefit of our HUM approach is the uncertainty estimation. Thus, it seems that our HUM approach on a difficult computer vision task is efficient.

In Table . 3 we provide some of the most uncertain examples according to our HUM CNN: we took the top-8 uncertain (considering p-q 1 2 as an evaluation of uncertainty) images of our testing set (that includes cats that were not seen during training). Qualitatively, we see that just taking the entropy of a regular CNN as an uncertainty criterion does not perform as well in Table. 4. In this work, we revisited an old problematic in Machine Learning and more generally in Statistics about dealing with estimation error that is inherent to realworld studied phenomena. An interesting avenue of research could be to refine the "far from dataset" kind of uncertainty: revisiting the Wasserstein distance between the query Dirac and the training dataset would pragmatically give a geometric aspect to uncertainty. This field of research has a wide of applications: active learning would choose examples based on predicted uncertainty which is intuitively (at least) reasonable:

to discover what we are unsure of first. There is also another promising application: reinforcement learning without reward to learn how to navigate in a new environement without task which would give some complementary understanding to the notion of curiosity in the scientific effort of [START_REF] Burda | Largescale study of curiosity-driven learning[END_REF] and we could consider that diminishing an uncertainty estimation would be the cornerstone of a reinforcement learning policy. Indeed, even in regular reinforcement learning, training inescapably boils down to distentangling an exploration / exploitation dilemma [START_REF] Sutton | Reinforcement learning: An introduction[END_REF]: on the one hand, exploration is testing unusual (or risky) patterns of actions with the hope to expect possibly higher rewards and on the other hand, exploitation is benefiting from known (or conservative) patterns of actions in order to better guarantee a certain level of rewards. This contribution about uncertainty estimation could help breaking the dilemma. The rationale behind could consist in saying that one should explore where uncertainty is high and exploit where uncertainty is low and value estimation is high.

It seems interesting to also investigate the blurry prediction provided by uncertainty models to help for more robustness to adversarial examples (example designed to fool supervisedly trained algorithms, see the survey of [START_REF] Yuan | Adversarial examples: Attacks and defenses for deep learning[END_REF] for details). At least, this kind of systems would be fooled but the predictions would have larger standard deviation for example for this deceiving points.

C O N C L U S I O N

In this dissertation, learning representations using neural networks and optimal transport provided promising statistical and algorithmic tools for clustering, for an attempt towards unsupervised feature importance and some practical supervised prediction with uncertainty. Three axes have been present throughout these contributions in large scale settings for thrice: clustering for simplifying the cardinality axis, unsupervised feature importance for simplifying the dimensionality axis and finally uncertainty estimation for the output space axis.

For further research work, some interactions between our three contributions can be drawn: (i) optimizing over a set of Wasserstein metrics as we suggest in unsupervised feature importance extraction can be investigated for clustering in order to tackle the metric-invariance, (ii) unsupervised feature importance extraction and clustering can be done in a simultaneous fashion combining both feature selection and clustering which is desirable for high-dimensional data and (iii) out of distribution uncertainty estimation can be revisited with Wasserstein distances.

Away from traditional supervised algorithms, we have shown that neural networks are appealing tools especially when empowered by optimal transport and data analysis purposes. A lot of programming has been done to tackle these issues and interacting between real-world pragmatism and theoretical reasoning back and forth has been key to conduct this work. For dramatic computational speed improvements reasons, some pretraining techniques from 1990s are less popular nowadays in neural networks training. Recently, the rise of Generative Adversarial Networks (GANs) initiated by [START_REF] Goodfellow | Generative adversarial nets[END_REF] is revisiting a new kind of optimization: not an usual minimization problem but a min-max problem hence the adjective adversarial which is thus more difficult to monitor. In practice, for same domain (or medium) data, it is possible to use the same structure and hyper-parameters as other authors e. g. the DCGAN structure by [START_REF] Radford | Unsupervised representation learning with deep convolutional generative adversarial networks[END_REF] for images. Unfortunately, for custom data such as genomics data, there is a need for practical rule of thumbs to begin working for new kinds of data and neural networks and particularly wirh GANs.

Historically, neural networks experts used to train an unsupervised autoencoder, in order to keep the encoder as an initialization of the neural network now trained in a supervised fashion. This pre-training technique seems to save a lot of time during an era of expensive and long computations. The aim of this appendix is to follow this kind of reasoning to initialize GANs that are notoriously painful to train and monitor without further data knowledge.

S I M I L A R O B J E C T I V E S

Let's study the Wasserstein GAN [START_REF] Arjovsky | Wasserstein generative adversarial networks[END_REF] where the Wasserstein distance W(p, q G ) is defined between the distributions of real data p and generator data q G (defined by y ∼ q G ⇐⇒ z ∼ n, and y = G(z))

W(p, q G ) = min γ∈Γ(p,q G ) E (x,y)∼γ [ yx 2 ] = min γ ′ ∈Γ(p,n)

E (x,z)∼γ ′ [ G(z) -x 2 ] (166)
where Γ(p, q) stands for the set of couplings γ based on the two marginal distributions p and q. In the GAN literature, this is reformulated by thanks to the Kantorovich-Rubinstein duality but in this appendix, we will still use the Wasserstein formulation.

Indeed, if we look at that auto-encoder objective function made of an encoder E and decoder D: where the AE operator is simply the mean distance between the original data from p and the same data distorted by an encoder E followed by a decoder D:

AE(p, p E ,D ) = E x∼p [ D • E (x) -x 2 ] (168)
We believe that equations Eq. ( 166) and Eq. ( 168) are similar. Training an autoencoder and training a generative adversarial network are both unsupervised but auto-encoder are easier to train because there is only a one-way minimization which is not the case for GANs that have two-ways type of minimization (minmax) w.r.t. different parameters.

By training an auto-encoder, we get a reasonable reconstruction loss for output functions E encoder and D decoder. Now, it is possible to estimate the mean µ = 1 N ∑ N i=1 E (x i ) and covariance σ = 1 N ∑ N i=1 (E (x i )µ) 2 (with element- wise square root and square functions). The intuition is that for a Gaussian noise ǫ ∼ N (0 d , I d ), then the decoded points D(µ + σ × ǫ) (with element-wise product) has a distribution that is close to decoded codes D • E (x) which is close to the original data x. This way, we made a reasonable initialization for a generator made of a Gaussian noise into a linear layer from empirical mean µ and standard deviations σ followed by the decoder function D.

Early on since the beginning of the GAN literature in 2014, deep learning scientists saw the potential of seeing the min-max optimization scheme as the minimization of a discrepancy between real and generated data where that discrepancy required and adversarial optimization for estimation [START_REF] Goodfellow | NIPS 2016 tutorial: Generative adversarial networks[END_REF] (i. e. a maximization of an objective to be minimized in fine). Back to our work, at this point, the auto-encoder training was a traditional one-way minimization but now the adversarial network (or critic against the generator) is needed to measure the discrepancy we just mentioned between real and generated data. In our Wasserstein discrepancy case, this correponds to only maximization via the Kantorovich-Rubinstein duality of the Wasserstein distance. Once both the generator from the auto-encoder minimization and next the critic from the discrepancy maximization for estimation are done. We actually not only know how far our generator (i. e. decoder) combined with our noise generator is from the real data distribution but we also get a convincing initialization for regular GANs training.

To sum up, we propose three steps:

1. Auto-Encoder Minimization for the euclidean reconstruction error to be low;

2. Critic Maximization until real v. s. generated discrepancy has plateaued out;

3. Regular GAN Min-Max optimization with a generator initialized by a Gaussian fitting the encoded codes followed by the decoder on the one side and the critic previously maximizing the real v. s. generated discrepancy on the other side;

which is much easier for monitoring and prototyping as getting rid off an unpredictable objective variations behavior especially at the first epochs contributes to acceleration of development and robustness w.r.t. uncontrolled objective variations. In these conditions, finding hyper-parameters is much more simplified. Generative Adversarial Networks [START_REF] Goodfellow | Generative adversarial nets[END_REF] really attracted the research community attention when image rendering started to get convincing especially with the spark ignition DCGAN work of [START_REF] Radford | Unsupervised representation learning with deep convolutional generative adversarial networks[END_REF] with a recent research and engineering accomplishment in extremely large scale settings of BigGAN [START_REF] Brock | Large scale GAN training for high fidelity natural image synthesis[END_REF] and mainstream press demos4 . Following this trend, we use the same neural network structure as DCGAN [START_REF] Radford | Unsupervised representation learning with deep convolutional generative adversarial networks[END_REF] shown in Fig. 27 for a Wasserstein GAN trained on CIFAR-10 images dataset [?]. In this section, we use Adam [START_REF] Kingma | A method for stochastic optimization[END_REF] default parameter with 10 -6 for generator learning rate and 10 -5 for critic learning rate and also 10 -5 the one-way minimization learning rate. Moreover, in min-max optimizations, we do 10 max iterations followed by 1 min iteration for each data mini-batch. The first step according to our section 4.6 consists in training an auto-encoder. Fig. 28 is showing the evolution of the minimized loss w.r.t. the number of iterations.

Second step consists in estimating the mean µ and standard deviation σ of the codes coordinates (codes are encoded version of data). Then a critic maximization is taking place so that we can estimate the Wassertein distance between the real data and the decoded Gaussian noise associated with N (µ, σ) as shown in Fig. 29. Here, the Gaussian parameter (µ, σ) and the decoder playing now the role of generator are all remaining frozen (i. e. with zero learning rate) so that the critic neural network can be considered as just a tool to estimate the Wasserstein distance between fixed paramaters noise and generator data distribution and real data distribution in a Kantorovich-Rubinstein maximization procedure. Third and final step is doing the original optimization scheme but now from a relevant starting point (which is crucial in non-convex deep learning optimization notoriously prone to spurious extrema). Fig. 30 is showing a direct minimization in spite of the adversarial min-max learning procedure. Here all parameters are allowed to vary meaning that the first two steps are just initialization steps following the end-to-end recommendations.

All these steps (except the first and easy auto-encoder step) are greatly simplified thanks to the ease of using the Lipschitz property for the critic neural network offered by [START_REF] Miyato | Spectral normalization for generative adversarial networks[END_REF] (this would not be reasonable with clipping [START_REF] Arjovsky | Wasserstein generative adversarial networks[END_REF] nor regularization term [START_REF] Gulrajani | Improved training of wasserstein gans[END_REF]). Compared to direct GAN training shown in Fig. 32, with fewer iterations and less tedious selection of learning rates, our approach gets better results. Indeed, one-way optimization schemes are easier than two-ways min-max optimization ones to handle. With mastered tools described by [START_REF] Lecun | Efficient backprop[END_REF] and decades of shared good practice, one-way optimization (say minimization and not min-max optimization) is well studied nowadays. This is to be put in contrast with the two opposite learning rates required to the direct stochastically alternated method. To be honest, in our case, the final step has the same min-max GAN-like algorithm but the variations are better controlled in practice because the system is better initialized.

Figure 32: Generated samples of our Wasserstein GAN trained on CIFAR-10

  C A R D I N A L I T Y N, the number of data points for train; D I M E N S I O N A L I T Y D, the dimensionality of one data point; O U T P U T D I M E N S I O N K, the dimensionality of the automatic output decision

Figure 1 :

 1 Figure 1: Several input/processing/output scenarios from standard neural networks to recursive ones. Red: Input, Green: Processing, Blue: Output. From left to right: one (input) to one (output) vanilla neural network, one to many like in image captioning (one image to many words in a sentence), many to one like in sentiment analysis (many words in a sentence to a category of mood), and the (delayed or not) many to many case like in langauage translation (many words of a sentence in one language to many words of a new sentence in a different language). Diagrams taken from the pedagogical blog of Andrej Karpathy: http://karpathy.github.io/ 2015/05/21/rnn-effectiveness/

Figure 2 :

 2 Figure 2: Images Datasets Explosion

Figure 3 :

 3 Figure 3: GAN Principle: Imitating Data thanks to some Artificial Low Dimensional Noise in Purple from a Code Space (R d ) transformed thanks to a Generator function (or forger) into Generated Data in Red living in the Data Space ((R D )) so that they are supposed to be close to the Real Data in Blue thanks to the Critic function (or detective) -Adapted from https://optimaltransport.github.io

Figure 4 :

 4 Figure 4: Volume of an euclidean sphere of unit radius in dimension D

  Figure 5: 0.1-peel ratio R(D)f or an euclidean sphere of unit radius in dimension D

Figure 6 FrequencyFigure 7 :

 67 Figure 6: "Two Moons" Toy Dataset

FrequencyFigure 8 :Figure 9 :

 89 Figure 8: Two Moons Projection through kernel-PCA into 1D (with Gaussian Kernel)

S

  C A L E I N VA R I A N C E Clustering output should not change if we multiply data by a constant R I C H N E S S O R C L U S T E R S H A P E S I N VA R I A N C E all separable cluster shapes should be possible (e. g. beyond linear separation or ball-shaped clusters) C O N S I S T E N C Y O R M E T R I C I N VA R I A N C E Clustering output should not change with respect to the choice of distance

  Figure 10: Data and codes spaces

2 Figure 11 :

 211 Figure 11: Graphical Model for Data Generation

Figure 12 :

 12 Figure 12: GeWaC Optimization Scheme

: 2 :

 2 for j = 1, . . . , N critic do Sample some B × K Gaussian noise codes from a pseudorandom generator e i,k ∼ N (0, I d ) i = 1, . . . , B k = 1, . . . , K 3: Compute the critic evaluation on the mini-batch of decoded noise and its original version concatenation

  ) with α = 10 4 , β = 2 × 10 4 , σ = 10 -3 , b = [m, 0] ⊤ and c = [-M, 0] ⊤ (m = 9 and M = 100).

2 .

 2 1. Sample two clusters indices (k, k ′ ) independently from the multinomial distribution of parameters π Observe the Wasserstein distance W(p k , p k ′ ) between the components/clus- ters associated with k and k ′

Algorithm 4

 4 Optimization algorithm (step 4) 1: Input: Data (x i ) i=1,...,N Number of clusters K 2: Initialization:θ E and θ M = (ω k , µ k , S k ) k=1,...,K , # initialized fromsteps 1, 2 and 33: π ← ω # both π and ω are (re)-parametrized thanks to a softmax on free paramaters 4: T ← 0 5: while θ has not converged do -batch of size B from the dataset

N

  O R M A L I Z E D M U T U A L I N F O R M AT I O N ( N M I ) Measuring an estimator of dependence normalized by the entropy information contained in both signals is also a classic tool for clustering evaluation but this time evaluated on the divided confusion matrix M = M ∑ k,ℓ M k,ℓ • = ∑ ℓ Mk,ℓ and M•,ℓ = ∑ k Mk,ℓ N O R M A L I Z E D C O N D I T I O N A L E N T R O P Y ( N C E )

Figure 13 :

 13 Figure 13: Same problem as EM-GMM but with a Wasserstein GAN without a generator neural network but a Gaussian mixture generator instead. A is the red distribution and B the blue one and they get purple when in local superposition -best seen in colors

Figure 14 :

 14 Figure 14: Converged mixtures (purple because of the superposition of the red and blue mixtures) -best seen in colors

  .

  [2015] D-500-500-2000-d (D is the dimensionality of the input space e. g. D = 784 for MNIST and d = 10 is the dimensionality of 8 MLP stands for Multi-Layer Perceptron

  m

Figure 18 :

 18 Figure 18: Generated digits images. From left to right, we have the ten classes found by GeWaC and ordered thanks to the Hungarian algorithm. From top to bottom, we go further and further in random directions from the centroids (the first row being the decoded centroids). More specifically, u is sampled from the uniform random density on the unit hypersphere in the code space.

Figure 19 :

 19 Figure 19: Optimal and Natural Transports are the same in Euclidean Distance Case for close Disttributions

Figure 20 :

 20 Figure 20: Unsupervised Feature Importance Algorithm

  K E Y W O R D S Bayesian deep learning, Bayesian neural networks, Uncertainty, Confidence, Reparametrization Trick, Neural Networks, Distributions 4.1 I N T R O D U C T I O N

  is a random uniform index parsing the N-cardinality training dataset in mini-batches of size B totalling B gradient evaluations per mini-batch.

  4.4 A L G O R I T H M A N D P R A C T I C A L I M P L E M E N TAT I O N D E TA I L S

R

  E G R E S S I O N In our Gaussian regression case, µ can be reasonably well-initialized by its uncertainty-free already trained counterpart with a zero-initialized Cholesky function C for the covariance matrix. C = 0 ⇐⇒ no uncertainty (163) C L A S S I FI C AT I O N In our uniform classification case, both A and B can be initialized by the same weights (while the random generator of t will give different gradients starting from the very first mini-batch gradient computation).

Algorithm 5 2 : 6

 526 and we create 2 functions, namely A = G A • H and B = G B • H where function H is shared whereas both G A and G B are equally initialized by the previous training output G. Classification with Uncertainty 1: Input: Data (x i , y i ) i=1,...,N where x i ∈ R D and y i ∈ {0, 1} K ∩ ∆ K , a mini-batch size B and a Monte Carlo parameter M ∈ N * Initialization: θ = {θ H , θ G A , θ G B } initialized from a regu- lar uncertainty-free classification optimization scheme for respectively H, A and B; 3: while θ has not converged do -batch of size B from the dataset(x i b , y i b ) b=1,...,B where i b ∼ U N (1, N) 6:Sample BM values from a pseudo-random generatort b,m ∼ U R (0, 1) for b = 1, . . . , B and m = 1, . . . , M 7: Compute h b ← H(x i b ), a b ← G A (h b ) and b b ← G B (h b ) for b = 1, . . . , B 8: Apply the logsumexp trick a b ← a b -max k a (k) b and b b ← b b -max k b t b,m × SoftMax(a b ) + (1t b,m ) × SoftMax(b b ) 10: Perform a gradient descent step of L to update θ 11: end while Algorithm Regression with Uncertainty 1: Input: Data (x i , y i ) i=1,...,N where x i ∈ R D and y i ∈ R K , a mini-batch size B and a Monte Carlo parameter M ∈ N * 2: Initialization: θ µ for a mean function µ initialized from a regular uncertainty-free regression optimization scheme; θ C for a lower triangular matrix function C initialized by a diagonal of training standard deviations; θ = {θ µ , θ C } 3: while θ has not converged do 4: Free all gradients accumulators 5: Sample a mini-batch of size B from the dataset (x i b , y i b ) b=1,...,B where i b ∼ U N (1, N) 6: Sample BM values from a pseudo-random generator ǫ b,m ∼ N (0 K , I K ) 7: Compute m b ← µ(x i b ) and S b ← C(x i b )

Figure 21 :

 21 Figure 21: "Two Moons" Dataset in two Scenarios: a clear and easy case (left) and a noisy and difficult case (right)

Figure 22 : 2 Figure 23

 22223 Figure 22: "Two Moons" Map of Prediction Probability of being blue a+c 2

Figure 24 :

 24 Figure 24: Two examples of images from our "Dogs and Wolves" dataset taken from Imagenet

Figure 25 :

 25 Figure 25: Our Convolutional Neural Network in two parts separated by the tick blue line: (i) some shared convolutional and fullyconnected layers on the left and (ii) a few fully-connected layers on the right to finally output the two bounds of our uniform prediction law.

Figure 26 :

 26 Figure 26: Validation score w.r.t. the training epochs comparing our HUM CNN with its initialized regular CNN counterpart (the higher, the better, 1.0 means 100% of good classification on the validation set)

  A P P E N D I X : G E N E R AT I V E A D V E R S A R I A L N E T W O R K S P R E T R A I N I N G I N T R O D U C T I O N

  objective function with a generator G, a real data distribution p living in R D , and noise generator distribution n living in R d min G W(p, q G ) (165)

E

  X P E R I M E N T S T R A I N I N G D C G A N O N C I FA R -1 0

Figure 27 :

 27 Figure 27: DCGAN encoder structure, the critic has same structure(with an additional 100-1 perceptron layer) and the decoder/generator has a symmetric structure

Figure 28 :

 28 Figure 28: Autoencoder mini-batch loss w.r.t. the number of iterations (the lower the better, 50 iterations for 1 epoch)

Figure 30 :

 30 Figure 30: Our pre-trained GAN mini-batch objective w.r.t. the number of iterations. (the lower the better, 550 iterations for 1 epoch as each data mini-batch is seen 10 times for the critic and once for the generator)

  The prediction could be wrong because test data do not come from the same distribution as training data. Out from the training distribution, test data go against one of the most fundamental machine learning hypothesis to make any system work. For example, on-

line self-adapting systems

[START_REF] Bertsekas | Dynamic programming and optimal control[END_REF]

look like reinforcement learning systems and deal with out of distribution data uncertainty purposefully. During training, prediction systems never see anything but data coming from training data. Neural networks have the deserved reputation of quickly over-fitting on training data empirical distribution (if used carelessly) which is both a warning against both interpolation and extrapolation data prediction in worst case scenarios as explained by ? partially avoided by regularization

Table 2 :

 2 Experimental accuracy results (%, the higher, the better) based on the Hungarian method. (the last rows correspond to methods without neural networks )

	Datasets	MNIST Reuters Reuters-10k HHAR
	DiWaC (ours)	98.42	84.24	84.87	92.42
	GeWaC (ours with fixed proportions from AE+GMM )	97.37	82.14	82.27	87.54
	ClusterGAN [Mukherjee et al., 2019]	90.97	-	-	-
	VaDE [Jiang et al., 2016]	94.06	79.38	79.83	84.46
	DEC [Xie et al., 2015]	84.30	75.63	72.17	79.82
	AE + GMM (full covariance)	82.56	70.98	70.12	78.48
	IMSAT [Hu et al., 2017]	98.40	-	71.00	-
	GAR [Kilinc and Uysal, 2018]	98.32	-	-	-
	DEPICT [Dizaji et al., 2017]	96.50	-	-	-
	GMM (diagonal covariance)	53.73	55.81	54.72	60.34
	k-Means	53.47	53.29	54.04	59.98

  I FI C AT I O N For a high number of possible categories (K ≃ 10 3 in ImageNet[START_REF] Fei-Fei | Imagenet: crowdsourcing, benchmarking and other cool things[END_REF]), allowing a mixture of a much reduced number C ∼ 10 1 of components is useful to handle both uncertainty and close categories in a more interpretable way than in the usual supervised deep learning apparatus.

	e
	4.3.4 Links with Other techniques

https://www.deepl.com/press.html

Multi-Layered Perceptron

http://yann.lecun.com/exdb/mnist/

see this scientific blog for details: https://gregorygundersen.com/blog/ 2018/04/29/reparameterization

EM-GMM and k-Means++ provided by scikit-learn[START_REF] Pedregosa | Scikit-learn: Machine learning in Python[END_REF]

The "clustering no-free lunch theorem" is not an official nickname but we choose it because it helps understanding the work of[START_REF] Kleinberg | An impossibility theorem for clustering[END_REF]. The original authors would rather use the clustering "impossibility theorem".

see http://yann.lecun.com/exdb/mnist for a complete supervised MNIST benchmark although they train with supervised labels whereas we do not have this supervised information.

3.4 A L G O R I T H M

https://www.youtube.com/watch?v=FD8l2vPU5FY

http://harchaoui.org/warith/dogs_wolves.zip

http://harchaoui.org/warith/imagenet_test_resnet152_pytorch.zip

https://www.thispersondoesnotexist.com