Bilal Piot
email: bilal.piot@univ-lille3.fr

Matthieu Geist
email: matthieu.geist@centralesupelec.fr

Olivier Pietquin
email: olivier.pietquin@univ-lille1.fr

Imitation Learning Applied to Embodied Conversational Agents

Embodied Conversational Agents (ECAs) are emerging as a key component to allow human interact with machines. Applications are numerous and ECAs can reduce the aversion to interact with a machine by providing user-friendly interfaces. Yet, ECAs are still unable to produce social signals appropriately during their interaction with humans, which tends to make the interaction less instinctive. Especially, very little attention has been paid to the use of laughter in humanavatar interactions despite the crucial role played by laughter in human-human interaction. In this paper, methods for predicting when and how to laugh during an interaction for an ECA are proposed. Different Imitation Learning (also known as Apprenticeship Learning) algorithms are used in this purpose and a regularized classification algorithm is shown to produce good behavior on real data.

Introduction

An important challenge for the future of computer science is to build efficient and user-friendly human-machine interfaces. This will enable a large public to interact with complex systems and reduce considerably the technological gap between people. In the last decade, Embodied Conversational Agents (ECAs), also called avatars, emerged as such interfaces. However, their behavior appear quite unnatural to most users. One possible explanation of this bad perception is the inability of ECAs to make a proper use of social signals, even though there exists some research on this subject [START_REF] Schröder | Building autonomous sensitive artificial listeners[END_REF]. One of these sig- nals is laughter, which is a prominent feature used by human during interactions. Yet, very little attention has been paid to enable ECAs with laughter capabilities until recently [START_REF] Niewiadomski | Laugh-aware virtual agent and its impact on user amusement[END_REF]Piot et al., 2014b).

Enabling ECAs with laughter capabilities is not only about being able to synthesize audio-visual laughter signals [START_REF] Niewiadomski | Visual laughter synthesis: Initial approaches[END_REF][START_REF] Urbain | Evaluation of HMMbased laughter synthesis[END_REF]. It is also concerned by an appropriate management of laughter during the interaction, which is a sequential decisionmaking problem. Thus, there is a need for a laughterenabled interaction manager, able to decide when and how to laugh so that it is appropriate in the conversation. This being said, it remains uneasy to define what is an appropriate moment to laugh and to choose the good type of laugh.

More formally, the task of the laughter-enabled Interaction Manager (IM) is to take decisions about whether to laugh or not and how. These decisions have to be taken according to the interaction context which can be inferred from laughter, speech and smile detection modules (detecting social signals emitted by the users) implemented in the ECA but also by the task context (for example, if the human is playing a game with the ECA, what is the status of the game). This is a sequential decision-making problem. Thus, the IM is a module implementing a mapping between contexts (or states noted s 2 S) and decisions (or actions noted a 2 A). Let's call this mapping a policy, noted ⇡(s) = a. This mapping is quite difficult to learn from real data as the laughs are quite rare and very different from one user to another.

In this paper, we describe the research results for learning such a mapping from data, recorded during some humanhuman interactions, so as to implement, in the IM, a behavior similar to the one of a human. Imitation Learning (IL), also known as Apprenticeship Learning (AL), methods are considered. Indeed, during some human-human interactions, one human can be considered as an expert which actions should be imitated by the IM. This can be framed as a Learning from Demonstrations (LfD) problem as LfD is a paradigm in which an agent (called the apprentice) learns to behave in a dynamical environment from demonstrations of another agent (named the expert). IL is a framework that offers methods that solve the LfD problem. Most IL methods are classification-based and do no take into account the underlying dynamics of the environment (the humanhuman interaction in our case), which can be problematic. Indeed, choosing when and how to laugh may depend on previous states of the interaction and thus the dynamics of the interaction may be an important aspect to take into consideration. However, a recent algorithm called Regularized Classification for Apprenticeship Learning (RCAL) (Piot et al., 2014a) uses the underlying information of dynamics existing in the data. Thus, pure and regularized classification [START_REF] Taskar | Learning structured prediction models: A large margin approach[END_REF]Piot et al., 2014a) methods are compared and regularized classification is shown to efficiently learn a behavior on data sets of real laughs in a natural interaction context in Sec. 4

The remainder of the paper is organized as follows. First, we present the LfD paradigm and a large-margin classification method in Sec. 2. Then, we use the large margin method to derive the RCAL algorithm (Piot et al., 2014a) in Sec. 3. Finally, we compare several IL methods and RCAL on some real human-human interactions data.

Background

Before introducing a regularized classification algorithm for IL (RCAL), it is necessary to recall some notions and definitions relatives to the IL paradigm. As IL is an implementation of the general LfD problem, which can be formalized properly thanks to the concept of Markov Decision Process (MDP), the remaining of the section is organized as follows. First, we provide some definitions relative to the concept of MDP, then we present the general problem of LfD and finally we show how a large-margin Multi-Class Classification (MCC) which is an IL method can solve the LfD problem.

Markov Decision Process

A finite MDP [START_REF] Puterman | Markov decision processes: Discrete stochastic dynamic programming[END_REF] models the interactions of an agent evolving in a dynamic environment and is represented by a tuple M = {S, A, R, P, } where S = {s i } 1iN S is the state space, A = {a i } 1iN A is the action space, R 2 R S⇥A is the reward function (the local representation of the benefit of doing action a in state s), 2]0, 1[ is a discount factor and P 2 S⇥A S is the Markovian dynamics represented as a function from S ⇥ A to S where S is the set of probability distributions over S. The Markovian dynamics P gives the probability, P (s 0 |s, a), to reach s 0 by choosing the action a in the state s. A Markovian stationary and deterministic policy ⇡ is an element of A S and defines the behavior of an agent. In order to quantify the quality of a policy ⇡ relatively to the reward R, we define the value function. For a given MDP M = {S, A, R, P, } and a given policy ⇡ 2 A S , the value function

V ⇡ R 2 R S is defined as V ⇡ R (s) = E ⇡ s [ P +1 t=0 t R(s t , a t )]
, where E ⇡ s is the expectation over the distribution of the admissible trajectories (s 0 , a 0 , s 1 , . . . ) obtained by executing the policy ⇡ starting from s 0 = s. Moreover, the function

V ⇤ R 2 R S , defined as V ⇤ R = sup ⇡2A S V ⇡ R , is called the optimal value func- tion. A useful tool is, for a given ⇡ 2 A S , the action-value function Q ⇡ R 2 R S⇥A : Q ⇡ R (s, a) = R(s, a) + E P (.|s,a) [V ⇡ R ].
It represents the quality of the agent's behavior if it chooses the action a in the state s and then follows the policy ⇡.

Moreover, the function

Q ⇤ R 2 R S⇥A defined as: Q ⇤ R = sup ⇡2A S Q ⇡
R is called the optimal action-value function. In addition, we have that [START_REF] Puterman | Markov decision processes: Discrete stochastic dynamic programming[END_REF]. Thus, we have 8s 2 S, 8a 2 A:

Q ⇡ R (s, ⇡(s)) = V ⇡ R (s) and that max a2A Q ⇤ R (s, a) = V ⇤ R (s)
R(s, a) = Q ⇤ R (s, a) X s 0 2S P (s 0 |s, a) max a2A Q ⇤ R (s 0 , a). (1)
Eq. ( 1) links the reward R to the optimal action-value function Q ⇤ R . This equation will be useful in the sequel, as one could obtain R from Q ⇤ R by a simple calculus.

Learning from Demonstrations

LfD is a paradigm in which an agent (called the apprentice) learns to behave in a dynamical environment from demonstrations of another agent (named the expert). To address this problem, we place ourselves in the MDP framework which is used to describe dynamical systems as a set of states, actions and transitions. In this framework, the learnt behavior takes the form of a policy. More precisely, using MDPs, solving the LfD problem consists in finding the policy of the expert agent in states unvisited by the expert.

The expert agent is supposed to act optimally (with respect to the unknown reward function) and the apprentice can only observe the expert policy ⇡ E via sampled transitions of ⇡ E . Moreover, we suppose that the apprentice has some information about the dynamics which he could have collected by previous interactions. The aim of the apprentice is of course to find a policy ⇡ A which is as good as the expert policy with respect to the unknown reward.

More precisely, we suppose that we have a fixed data-set of expert sampled transitions

D E = (s i , ⇡ E (s i ), s 0 i ) {1iN E } where s i ⇠ ⌫ E 2 S and
s 0 i ⇠ P (.|s i , ⇡ E (s i )).
In addition, we suppose that the apprentice has some information about the dynamics via a fixed data-set of sampled transitions D P = (s j , a j , s 0 j ) {1jN P } where s j ⇠ ⌫ P 2 S and s 0 j ⇠ P (.|s j , a j ). We have D E ⇢ D P and no particular assumptions are made considering the choice of the action a j or the distributions ⌫ E and ⌫ P which can be considered unknown. Those requirements (used for example by [START_REF] Klein | Inverse reinforcement learning through structured classification[END_REF]; [START_REF] Boularias | Relative entropy inverse reinforcement learning[END_REF]) are not strong and can be fulfilled by real-life applications.

One can argue that having a data-set of sampled transitions D P is a strong assumption. However, the presented algorithms can be run with D P = D E which is the case shown in the experiments (see Sec. 4).

The large margin approach for Classification

To tackle the problem of LfD, it is possible to reduce it to an MCC problem [START_REF] Pomerleau | Alvinn: An autonomous land vehicle in a neural network[END_REF][START_REF] Ratliff | Imitation learning for locomotion and manipulation[END_REF][START_REF] Ross | Efficient reductions for imitation learning[END_REF][START_REF] Syed | A reduction from apprenticeship learning to classification[END_REF]. The goal of MCC is, given a training set D = (x i 2 X, y i 2 Y ) {1iN } where X is a compact set of inputs and Y a finite set of labels, to find a decision rule g 2 Y X that generalizes the relation between inputs and labels. [START_REF] Ratliff | Imitation learning for locomotion and manipulation[END_REF] use a large margin approach which is a scorebased MCC where the decision rule g 2 Y X is obtained via a score function q 2 R X⇥Y such that 8x 2 X, g(x) 2 argmax y2Y q(x, y). The large margin approach consists, given the training set D, in solving the following optimization problem:

q ⇤ = argmin q2R X⇥Y J(q), (2)
J(q) = 1 N N X i=1 max y2Y {q(xi, y) + l(xi, yi, y)} q(xi, yi),
where l 2 R X⇥Y ⇥Y + is called the margin function. If this function is zero, minimizing J(q) attempts to find a score function q ⇤ for which the example labels are scored higher than all other labels. Choosing a nonzero margin function improves generalization [START_REF] Ratliff | Imitation learning for locomotion and manipulation[END_REF]. Instead of requiring only that the example label is scored higher than all other labels, one requires it to be better than each label y by an amount given by the margin function. Another way to improve generalization is to restrain the search of q ⇤ to an hypothesis space H ⇢ R X⇥Y . However, it introduces a bias.

Applying the large margin approach to the LfD problem is straightforward. From the set of expert trajectories D E , we extract the set of expert state-action couples DE = (s i , ⇡ E (s i )) {1iN E 2N ⇤ } and we try to solve:

q ⇤ = argmin q2H⇢R S⇥A J(q) (3) J(q) = 1 NE N E X i=1 max a2A {q(si, a) + l(si, ⇡E(si), a)} q(si, ⇡E(si)).
The policy outputted by this algorithm would be ⇡ A (s) 2 argmax a2A q(s, a) where q is the output of the minimization. The advantages of this method are its simplicity and the possibility to use a boosting technique [START_REF] Ratliff | Imitation learning for locomotion and manipulation[END_REF] to solve the optimization problem given by Eq. (3). However, this is a pure classification technique which does not take into account the dynamics information contained in the sets D E and D P . In the following section, we present an algorithm, RCAL (Piot et al., 2014a), based on the large margin approach which uses the dynamics information by adding an original regularization term to J(q).

Regularized Classification

In this section, RCAL , a non-parametric AL algorithm using the information contains in the dynamics, is introduced.

First, it is important to remark that the large-margin classification problem described by Eq. ( 2) tries to find a function minimizing an empirical criterion obtained from sparse data in a given set of functions called the hypothesis space. This is, in general, an ill-posed problem (infinite number of solutions when the hypothesis space is rich enough) and a way to solve it is the regularization theory of [START_REF] Tikhonov | Methods for solving ill-posed problems[END_REF], which adds a regularization term that can be interpreted as a constraint on the hypothesis space. [START_REF] Evgeniou | Regularization networks and support vector machines[END_REF] show how the work of [START_REF] Vapnik | Statistical learning theory[END_REF] set the foundations for a general theory which justifies regularization in order to learn from sparse data. Indeed, the basic idea of Vapnik's theory is that the search for the best function must be constrained to a small (in terms of complexity) hypothesis space. If the hypothesis space is too large, a function that exactly fits the data could be found but with a poor generalization capability (this phenomenon is known as over-fitting). [START_REF] Evgeniou | Regularization networks and support vector machines[END_REF] show that the choice of the regularization parameter corresponds to the choice of an hypothesis space: if is small the hypothesis space is large and vice-versa. In the large margin framework, a general and natural way to introduce regularization is to add to J(q) a regularization term W (q) and to consider the following optimization problem:

q ⇤ = argmin q2H⇢R S⇥A J W (q) = argmin q2R S⇥A (J(q) + W (q)) .
where 2 R ⇤ + and W is a continuous function from R S⇥A to R + .

In order to introduce the dynamics information contained in the data-set D P , Piot et al. (2014a) assume that the unknown reward function for which the expert is optimal is sparse. This assumption helps to choose an appropriate regularization term which constrains the hypothesis space, hence reduces the variance of the method. To do so, Piot et al. (2014a) remark that a good score function must verify:

8s 2 S, ⇡ E (s) 2 argmax a2A q(s, a), which means that there exists a reward function R q 2 R S⇥A for which ⇡ E is optimal and such that q(s, a) =

Q ⇤ R (s, a) (see Sec. 2.1
). This reward is given via the inverse Bellman equation (Eq. ( 1)) 8s 2 S, 8a 2 A: R q (s, a) = q(s, a)

X s 0 2S P (s 0 |s, a) max a2A q(s 0 , a).
As the reward R q is assumed to be sparse, a natural choice for W (q) is kR q k 1,⌫ P where ⌫ P 2 S⇥A is the distribution from where the data are generated. However, as P (.|s, a) is unknown, it is not possible to compute R q (s, a). Thus, Piot et al. (2014a) rather consider, for each transition in D P , the unbiased estimate of R q (s j , a j ) noted Rq (j) which is can be obtained from data: Rq (j) = q(s j , a j ) max a2A q(s 0 j , a).

Therefore, in order to introduce the dynamics information contained in D P , the choice of the regularization term W is:

W (q) = 1 NP N P X j=1 | Rq(j)| = 1 NP N P X j=1 |q(sj, aj) max a2A q(s 0 j , a)|,
Even if W (q) is not an unbiased estimate of kR q k ⌫ P ,1 , it is shown by Piot et al. (2014a) that the constraints imposed by the regularization term W (q) are even stronger that the ones imposed by kR q k ⌫ P ,1 .

Thus, the algorithm RCAL consists in solving the following optimization problem:

q ⇤ = argmin q2H⇢R S⇥A J W (q) J W (q) = J(q) + N P N P X j=1 | Rq (j)|.
Then, the policy outputted by RCAL is ⇡ A (s) 2 argmax a2A q(s, a) where q is the output of the minimization of J W . In order to minimize the criterion J W (q), a boosting technique, initialized by a minimizer of J q for instance, is used (see the work of Piot et al. (2014a) for more details) in our experiments. Such a technique is interesting as it avoids to choose features which is often problem dependent. However, there is no guarantee of convergence.

4 An ECA deciding when and how to laugh

In order to build a laugh-aware agent, we searched for an interaction scenario which is realistic and simple enough to be set up. Thus, we opted for the one used by Niewiadomski et al. ( 2013) that implies telepresence. This scenario involves two subjects watching a funny stimulus. These two persons do not share the same physical space: they watch the same content simultaneously on two separate displays. However, they can see each other reactions in a small window placed on the top of the displayed content. This scenario corresponds to very common situations in real life when someone wants to share interesting content over the web. Such interactions are recorded thanks to cameras and microphones in order to build an ECA. More precisely, features such as laugh intensity, speech probability, smile probability extracted by analysis components are computed each 400 ms from the first person's recordings and form the state s (more details are provided by [START_REF] Niewiadomski | Laugh-aware virtual agent and its impact on user amusement[END_REF]). In addition, from the recordings of the second person (playing the role of the expert) is extracted the type of laugh which corresponds to the action to be imitated by the ECA and form the expert action ⇡ E (s). They are 4 types of laugh (thus 4 actions): strong laugh, normal laugh, quiet laugh and silence. Our data base named D E is composed of several expert trajectories with a total of 2378 state-action couples.

Thus, we have a data-set D E on which we can apply any IL method and more particularly RCAL (with D E = D P and = 0.1). [START_REF] Niewiadomski | Laugh-aware virtual agent and its impact on user amusement[END_REF] use a k-Nearest Neighbors (k-NN with k = 1) algorithm [START_REF] Cover | Nearest neighbor pattern classification[END_REF]. Here, we compare the results obtained via an averaged K-fold cross validation (K = 5) between RCAL, Classif (RCAL with = 0, that is the method of [START_REF] Ratliff | Maximum margin planning[END_REF] presented in Sec. 2.3), k-NN and a classification tree [START_REF] Breiman | Classification and regression trees[END_REF] In our data-set, the laughs and silences are unbalanced (approximately 80% of silence for 20% of laughs). That is why, it is a difficult task to classify correctly the laughs. We observe that the introduction of the regularization improves drastically the performance between RCAL with = 0.1 and the Classif algorithm (RCAL with = 0). Thus, taking into account the underlying dynamics helps to improve the performance. RCAL has also a better performance compared to the k-NN algorithm and the classification tree.

Conclusion and Perspectives

In this paper, a method for learning when and how an avatar (ECA) should laugh during an interaction with humans was presented. It is based on data-driven IL algorithms and especially on structured and regularized classification methods. It is shown, in a telepresence scenario, that RCAL outperformed other classification methods. Compared to previous experimentations [START_REF] Niewiadomski | Laugh-aware virtual agent and its impact on user amusement[END_REF]Piot et al., 2014b), this method provides better results as it takes into account the underlying dynamics of the interaction.

Here, LfD is reduced to an MCC problem. Yet, LfD can also be solved by other methods such as Inverse Reinforcement Learning (IRL) [START_REF] Russell | Learning agents for uncertain environments[END_REF][START_REF] Klein | Inverse reinforcement learning through structured classification[END_REF]. Actually, IRL has been shown to work better for some types of problems [START_REF] Piot | Learning from demonstrations: Is it worth estimating a reward function?[END_REF] and has already been used to imitate human users in the case of spoken dialogue systems [START_REF] Chandramohan | User simulation in dialogue systems using inverse reinforcement learning[END_REF]. Therefore, it seems natural to extend this work to IRL in the near future. Also, this method could be used to generate new simulation techniques for optimizing human machine interaction managers in other applications such as spoken dialogue systems [START_REF] Pietquin | A Probabilistic Framework for Dialog Simulation and Optimal Strategy Learning[END_REF][START_REF] Pietquin | A survey on metrics for the evaluation of user simulations[END_REF].

  Appearing in Proceedings of the 4 th Workshop on Machine Learning for Interactive Systems (MLIS) 2015, Lille, France. JMLR: W&CP volume 40. Copyright 2015 by the authors.

  on the data set D E .

	Algorithms	Global perfomance	Good laugh chosen	Performance on Silence
	RCAL	0.7722	0.3147	0.8872
	Classif	0.7671	0.2316	0.8934
	k-NN	0.7381	0.2681	0.8521
	Tree	0.7533	0.2914	0.8683


Timothé Collet
email: timothe.collet@centralesupelec.fr

Olivier Pietquin
email: olivier.pietquin@univ-lille1.fr

Optimism in Active Learning with Gaussian Processes

Keywords: Active Learning, Classifiaction, Multi-Armed Bandits

In the context of Active Learning for classification, the classification error depends on the joint distribution of samples and their labels which is initially unknown. The minimization of this error requires estimating this distribution. Online estimation of this distribution involves a trade-off between exploration and exploitation. This is a common problem in machine learning for which multi-armed bandit theory, building upon Optimism in the Face of Uncertainty, has been proven very efficient these last years. We introduce two novel algorithms that use Optimism in the Face of Uncertainty along with Gaussian Processes for the Active Learning problem. The evaluation lead on real world datasets shows that these new algorithms compare positively to state-of-the-art methods.

Introduction

Classification is a supervised learning framework consisting in the association of instances to labels, existing in finite number. It uses a set of instances already associated to labels called the training set. In order to get a high prediction accuracy, the training set should contain a high number of instances. Recently, it has become easy to collect, store and process large data sets. However, the association of each instance of the training set to labels requires the manual annotation by an expert. This task is time consuming, and may involve other costs.

In Active Learning, the goal is to minimize the number of requests to the expert needed to achieve a targeted performance of the classifier. Equivalently, to maximize its performance under a fixed budget of requests. Hence, an Active Learning algorithm dynamically constructs the training set by sequentially deciding which instance to present to the expert. The instance is chosen considering all the previously received labels, such that its inclusion in the training set lead to the best performance of the classifier. We work under the pool-based sampling scheme in which a set of unlabeled instances is available. The Active Learning algorithm successively selects instances from the pool to be labeled.

Among the many ways that exist to formulate this problem, which are reviewed in [START_REF] Settles | Active learning literature survey[END_REF], this work focus on the Error Minimization framework. In this framework, the selection criterion is derived from the true risk, or misclassification rate. It is thus directly linked to the (mostly used) measure of performance for classifiers. We distinguish two strategies in this framework: the function to be minimized is either the maximum or the average misclassification rate among the instance space. In the first case, the strategy is to estimate the true risk related to each unlabeled instance and to sample the one for which it is maximum. In the second case, the strategy is to simulate the sample of each instance and to effectively sample the one resulting in the maximum decrease of the true risk. This second strategy has the advantage of representativeness, which means that an instance with a lot of unlabeled data around it is preferred to an isolated one.

The true risk can not be known perfectly, since it requires the true label distribution which need to be learned. But it can be estimated, along with a measure of uncertainty, using the labels received so far. Among the methods that exist in Error Minimization, the difference lies in how to manage the uncertainty about the true risk. In [START_REF] Zhu | Combining active learning and semisupervised learning using gaussian fields and harmonic functions[END_REF], the true risk is directly replaced by its estimation. In [START_REF] Kapoor | Active learning with gaussian processes for object categorization[END_REF] and [START_REF] Kapoor | Selective supervision: Guiding supervised learning with decision-theoretic active learning[END_REF], the binary loss is expected over the possible values of the true conditional density of the class label. In [START_REF] Huang | Active learning by querying informative and representative examples[END_REF], a min-max approach is used to ensure a minimum decrease of the risk.

In this paper, we propose a new solution to manage the uncertainty in the Error Minimization problem based on the Optimism in the Face of Uncertainty approach. Indeed, querying an instance may be used in order to increase the performance of the classifier directly, but also to improve future decisions. This is thus a case of exploration/exploitation dilemma for which the Optimism in the Face of Uncertainty is a well-renown approach, particularly for finite budget problems. The method is thus to establish first a probability bound on the value of the ideal criterion. Then to select the instance for which it is maximum.

Optimism in the Face of Uncertainty has been successfully used for Active Learning. In [START_REF] Ganti | Building bridges: Viewing active learning from the multiarmed bandit lens[END_REF], the authors focused on the reduction of the version space. In [START_REF] Collet | Active learning for classification: An optimistic approach[END_REF], the instance space is partitioned and an allocation strategy is defined.

We apply our solution to Gaussian Processes. Indeed, they provide a simple Bayesian learning framework which outputs a probability distribution on the conditional density of the class labels. It is thus well fitted for Active Learning, and particularly for Optimistic methods.

We briefly review Gaussian Processes in Section 2. Then, we present our proposed algorithms for both strategies in Sections 3 and 4. We present the experiments end their results in Sections 5 and 6.

Gaussian Processes

A stochastic process is a generalization of a probability distribution to functions, each input is associated with a probability distribution on the output. In the case of Gaussian Processes (GP), this distribution is Normal and can therefore be characterized only by its mean and variance. It thus provides a simple Bayesian framework for the supervised learning problem.

Let us consider the instance space X and the label set Y . In binary classification, the label set is composed of two elements, here Y = {-1, 1}. In the case of noisy classification, the expert can be represented by a Bernoulli distribution:

P (y ∈ Y = 1|x ∈ X) = µ(x) + 1 2 ,
where µ(x) is the mean label for instance x. We use a GP to estimate the values of µ(x) from S.

Let us denote D p the set containing all the instances from the pool, it is split into D l and D u , respectively the sets of labeled and unlabeled examples of size n l and n u . Suppose that we are given a training set S = {(x i , y i )} i∈ 1,n l containing all the labeled instances and the associated label given by the expert.

A GP make use of a kernel function k(., .) to specify the covariance between outputs. Let K l,l = [k(x i , x j )] n l ×n l be the covariance matrix of labeled instances. Let K x,l = [k(x, x i )] xi∈D l be the covariance of an instance x ∈ X and the labeled instances.

Conditioning the joint Gaussian distribution of µ(x) on the observed labels gives the following joint posterior distribution:

µ(x)|S ∼ N (m S (x), σ 2 S (x)). with m S (x) = K x,l (K l,l + σ 2 I) -1 y l and σ 2 S (x) = k(x, x)-K x,l (K l,l + σ 2 I) -1 K x,l .
The GP classifier predicts labels according to the sign of the posterior distribution's mean.

l(x) = sign(m S (x)) (1)
Note that this mean is the same as the Regularized Least Square Regression (RLSR). The interest of GP is in Active Learning, where the confidence measure is able to tell in which region of the instance space the model has to be refined.

Local Risk Minimization

In this section, we propose an algorithm based on Optimism in the Face of Uncertainty that minimizes the maximum local risk of the classifier. We first define the loss function then derive an criterion that tend to minimize it. Suppose that we are given a budget n. At each time step t ∈ 1, n , the Active Learning algorithm selects an unlabeled instance x t ∈ D u , submits it to the expert, receives a label y t and adds the pair to the training set.

At any time step t, the current training set S t is used to train the GP. This one is then used to predict the label l t (x) for any instance x ∈ X from Eq. ( 1). The expected prediction error, or local risk, is the probability that the expert would give a different label:

r t (x) = 1+µ(x) 2 if l t (x) = -1 1-µ(x) 2 if l t (x) = +1 . ( 2
)
We want to minimize the maximum risk one would encounter by presenting an unknown instance x ∈ X. We therefore define the following loss:

L t = arg max x∈X r t (x). (3)
In order to minimize this kind of loss, the solution is to select the instance for which the local risk is maximum.

x t+1 = arg max x∈Du r t (x).
Indeed, sampling an instance will necessarily lower the risk at its location, sampling where the maximum risk is attained guarantees to decrease this maximum: However, this loss can not be used as such. Indeed, this suppose to know the true value of µ(x) at least for every unlabeled instances, which is unrealistic. If it was known, the prediction would be direct. We propose to use Optimism in the Face of Uncertainty in order to approach this solution knowing only a distribution of µ(x), given by the GP. In this paradigm, a confidence bound on the ideal criterion is establish relatively to a probability δ, then its upper bound is used as a heuristic. We show how to use the distribution of µ(x) to derive the distribution of the local risk.

Suppose that the label l t (x) = -1 is predicted. Then,

P(r t (x) ≥ t (x)) = P(µ(x) ≥ 2 t (x) -1).
The GP trained with labeled instances outputs the following distribution:

µ(x)|S t ∼ N (m St (x), σ 2 St (x)).
Thus,

P(r t (x) ≥ t (x) = δ ⇐⇒ 2 t (x) -1 = Φ -1 (δ, m St (x), σ 2 St (x)),
with Φ the cumulative distribution function of the Normal distribution. This applies symmetrically for l t (x) = 1. Let us remind that l t (x) = sign(m St (x)) from Eq (1), then:

t (x) = 1 + Φ -1 (δ, |m St (x)|, σ 2 St (x)) 2 . (4)
Thus, t bounds r t with probability 1 -δ. Following an optimistic approach our algorithm selects at each time step the unlabeled instance for which t is maximum:

x t+1 = arg max Du t (x p ). (5)
Minimizing the local risk is a solution that has shown good results in stateof-the-art. However, it does not consider the density of instances around the considered instances. Indeed, the criterion is computed using only their own prediction from the GP. In the next section, we change the loss to consider this.

Global Risk Minimization

In this section, we propose an algorithm based on Optimism in the Face of Uncertainty that minimizes the global risk of the classifier. The main improvement of this new algorithm compared to the one of the previous section is that it considers the representativeness (as defined in [?]) of an instance. We first the loss function being used and then show a criterion that tend to minimize it.

We now consider the expected error one would encounter by presenting an instance x ∈ X drawn from the instance distribution which is the global risk of the classifier:

R t (S t ) = X r t (x)dP (x).
It is impractical because of the integral and because we do not known exactly the instance distribution. However, we have access to a pool instances following this distribution. We can thus estimate the global risk from it. We therefore define the loss:

Rt (S t ) = x∈Dp r t (x).
Defining a strategy that optimally minimizes this loss is NP-hard. However, it is common to make a myopic approximation. Thus the ideal strategy is to simulate the sampling of all the instances in the pool of unlabeled instances. Then, to effectively sample the one which results in the highest decrease of the risk:

x t+1 = arg max x∈Du ∆ xs Rt (S t ), (6)
where

∆ xs R t (S t ) = Rt (S t ) -Rt (S t ∪ {x s , y s }),
where y s is the true label of x s . This suppose to know the true values of µ(x) and y s , which is unrealistic. We propose to use the Optimism in the Face of Uncertainty approach in order sample as close as possible to the ideal sampling but with only a distribution of µ(x), given by the GP. We show how to compute a confidence bound on the criterion from Eq. ( 6) relatively to a probability δ.

Let us consider the simulation of the sample of x s ∈ D u and that the label given by the expert is y s . Let l + t (x) = sign(m St∪{xs,ys} (x)). Then,

∀x ∈ D p , r t (x) = 1 lt(x) =sign(µ(x)) -1 l + t (x) =sign(µ(x)) , Thus, ∀x ∈ D p , r t (x) = 1 lt(x) =-1 -1 l + t (x) =-1 with prob. P(µ(x) < 0|S t ∪ {x s , y s }) 1 lt(x) =1 -1 l + t (x) =1
with prob. P(µ(x) ≥ 0|S t ∪ {x s , y s }). We can then combine the cases for every instance in the pool to deduct the probability: P(∆ xs R t (S t )|S t ∪ {x s , y s }). We can also combine the cases for y s , given that P(y s = +1) = P(µ(x s ) ≥ 0|S t ), to deduct: P(∆ xs R t (S t )|S t ). We then compute the probability bound e t relatively to a fixed probability δ as:

P(∆ xs R t (S t ) ≤ e t (x s )|S t ) = 1 -δ. ( 7
)
Following an optimistic approach, our algorithm selects at each time step the unlabeled instance for which e t is maximum:

x t+1 = arg max xs∈Du e t (x s ). ( 8
)
5 Experiments

We evaluated our algorithms on several datasets from the UCI repository [START_REF] Lichman | UCI machine learning repository[END_REF].

Each experiment consisted in a series of 1,000 runs. For each of these runs, the dataset was randomly divided in two equal parts: one was used as the pool of instances in which the algorithms were allowed to pick, the other one was used as a test set, hidden from the algorithm. Our algorithm is initialized with one labeled instance randomly drawn from the pool. At each time step, the accuracy of the classifier is recorded as the proportion of well classified instances on the test set. The global performance displayed on figures results from averaging the accuracy of every run.

We compare OLRM and OGRM with the following methods: (1) Random Sampling: each instance to be labeled is randomly drawn from the pool, (2) Uncertainty Sampling (US) [START_REF] Lewis | A sequential algorithm for training text classifiers[END_REF] samples instance closest to the boundary, (3) QUIRE [START_REF] Huang | Active learning by querying informative and representative examples[END_REF] min-max strategy for minimization of the global risk, (4) VOI [START_REF] Kapoor | Selective supervision: Guiding supervised learning with decision-theoretic active learning[END_REF] minimization of the expected global risk, and (5) GPAL (Uncertainty) [START_REF] Kapoor | Active learning with gaussian processes for object categorization[END_REF] minimization of the expected local risk.

All these algorithms use the GP or RLSR classifier. Indeed, different classifier have different intrinsic performances. Thus, we think it is not meaningful to compare Active Learning algorithms which are based on different classifiers. Here, the Gaussian kernel is used. The optimal parameter σ has been found using grid-search.

Results

Fig. 1[a-d] displays the classification accuracy for each evaluated methods for with different budget of labeled instances. Notice that even though the marker only appears ten times on each curve, the accuracy is plotted for every integer number of labeled instances. As some curves are hard to differentiate, Fig. 1[e-h] shows a zoomed version of, respectively, Fig. 1[a-d], focusing on those methods.

The parameter δ in OLRM and OGRM controls the exploration/exploitation trade-off. In the case of OGRM, this parameter was tuned using a gridsearch. It appeared that any value lower than 0.5 lead to almost exactly the same performance, with a slight tendency for values closer to 0.5. In the case of OLRM, no value of the parameter had an outstanding performance. The optimal parameter always depended on the number of labeled instances considered. However, optimistic algorithm are often studied through a finite budget framework. This suppose that the budget is known from the beginning and the parameter may be chosen accordingly. Though, we do not provide any theoretical analysis to chose its value. The displayed curve for this algorithm shows the maximum performance OLRM would have for each budget if the parameter was chosen accordingly.

We can see on Fig. 1[a-d] that the three algorithms based on minimization of the local risk have very close performance. By looking at the zoomed version on Fig. 1[e-h], we see that they are not exactly equivalent. First, notice that US is strictly equivalent with OLRM with δ = 0.5. While GPAL is designed as an improvement of US that considers the variance outputted by the GP, it is most of the time outperformed by it. We can also see that, for a majority budget, there exist a parameter of ORLM that can achieve better performance than US.

We observe on Fig. 1[a-d] that OGRM outperforms all the methods it is compared to. Using the global risk loss function inherently make the performance increase sooner. Indeed, it considers the density of instances, and has an a pri-ori about where to sample with no observations. It is therefore surprising that QUIRE behaves worse than local risk methods. On Fig. 1[a-b], VOI starts with the same performance as OGRM but does not follow. This is because it considers the decrease of the estimation of the risk and not the estimation of the decrease of the risk.

Conclusion

In this paper, we introduce a new way to deal with uncertainty in the Error Minimization problem that uses Optimism in the Face of Uncertainty. We present two algorithms that work with two different losses. One can be seen as a generalization of Uncertainty Sampling, while the second extend the former to consider the instances distribution. Evaluations on real-world datasets were made that showed that both algorithms compared positively to state-of-the-art methods. Future work will focus on theoretical analysis and on the adaptation to the multiclass case. We will also consider improving the noise distribution being used to stick to the expert's distribution.

Fig. 1 .

 1 Fig. 1. Evaluation of algorithms


Hadrien Glaude

Cyrille Enderli

Olivier Pietquin

Non-negative Spectral Learning for Linear Sequential Systems

Method of moments (MoM) has recently become an appealing alternative to standard iterative approaches like Expectation Maximization (EM) to learn latent variable models. In addition, MoM-based algorithms come with global convergence guarantees in the form of finite sample bounds. However, given enough computation time, by using restarts and heuristics to avoid local optima, iterative approaches often achieve better performance. We believe that this performance gap is in part due to the fact that MoM-based algorithms can output negative probabilities. By constraining the search space, we propose a non-negative spectral algorithm (NNSpectral) avoiding computing negative probabilities by design. NNSpectral is compared to other MoM-based algorithms and EM on synthetic problems of the PAutomaC challenge. Not only, NNSpectral outperforms other MoM-based algorithms, but also, achieves very competitive results in comparison to EM.

Introduction

Traditionally, complex probability distributions over structured data are learnt through generative models with latent variables (e.g Hidden Markov Models (HMM)). Maximizing the likelihood is a widely used approach to fit a model to the gathered observations. However, for many complex models, the likelihood is not convex. Thus, algorithms such as Expectation Maximization (EM) and gradient descent, which are iterative procedures, converge to local minima. In addition to being prone to get stuck into local optima, these algorithms are computationally expensive, to the point where obtaining good solutions for large models becomes intractable. A recent alternative line of work consists in designing learning algorithms for latent variable models exploiting the so-called Method of Moments (MoM). The MoM leverages the fact that low order moments of distributions contain most of the distribution information and are typically easy to estimate. The MoM have several pros over iterative methods. It can provide extremely fast learning algorithms as estimated moments can be computed in a time linear in the number of samples. MoM-based algorithms are often consistent with theoretical guarantees in form of finite-sample bounds. In addition, these algorithms are able to learn a large variety of models [1], [START_REF] Bailly | A spectral approach for probabilistic grammatical inference on trees[END_REF], [START_REF] Thon | Links between multiplicity automata, observable operator models and predictive state representationsa unified learning framework[END_REF], [START_REF] Glaude | Subspace identification for predictive state representation by nuclear norm minimization[END_REF]. In a recent work, [START_REF] Thon | Links between multiplicity automata, observable operator models and predictive state representationsa unified learning framework[END_REF] showed that numerous models are encompassed into the common framework of linear Sequential Systems (SSs) or equivalently Multiplicity Automata (MA). Linear SSs represent real functions over a set of words f M : Σ → IR, where Σ is an alphabet. In particular, they can be used to represent probability distributions over sequences of symbols. Although, for the sake of clarity, we focus on learning stochastic rational languages (defined in Section 2), our work naturally extends to other equivalent or encompassed models: Predictive State Representations (PSRs), Observable Operator Models (OOMs), Partially Observable Markov Decision Processes (POMDPs) and HMMs.

Beyond all the appealing traits of MoM-based algorithms, a well-known concern is the negative-probabilities problem. Actually, most of MoM-based algorithms do not constrain the learnt function to be a distribution. For some applications requiring to learn actual probabilities, this is a critical issue. For example, a negative and unnormalized measure does not make sense when computing expectations. Sometimes, an approximation of a probability distribution is enough. For instance, computing the Maximum a Posteriori (MAP) only requires the maximum to be correctly identified. But even for these applications, we will show that constraining the learned function to output non-negative values helps the learning process and improves the model accuracy. A second concern is the observed performance gap with iterative algorithms. Although they usually need restarts and other heuristics to avoid local minima, given enough time to explore the space of parameters, they yield to very competitive models in practice. Recently, an empirical comparison [START_REF] Balle | Methods of moments for learning stochastic languages: Unified presentation and empirical comparison[END_REF] have shown that MoM-based algorithms still perform poorly in comparison to EM.

In this paper, we propose a new MoM-based algorithm, called non-negative spectral learning (NNSpectral) that constraints the output function to be non-negative. Inspired by theoretical results on MA defined on non-negative commutative semi-rings, NNSpectral uses Non-negative Matrix Factorization (NMF) and Non-Negative Least Squares (NNLS). An empirical evaluation on the same twelve synthetic problems of the PAutomaC challenge used by [START_REF] Balle | Methods of moments for learning stochastic languages: Unified presentation and empirical comparison[END_REF] allows us fairly comparing NNSpectral to three other MoM-based algorithms and EM. Not only, NNSpectral outperforms previous MoMbased algorithms, but also, it is the first time to our knowledge that a MoM-based algorithm achieves competitive results in comparison to EM.

Multiplicity Automata

Definition

Let Σ be a set of symbols, also called an alphabet. We denote by Σ , the set of all finite words made of symbols of Σ, including the empty word ε. Words of length k form the set Σ k . Let u and v ∈ Σ , uv is the concatenation of the two words and uΣ is the set of finite words starting by u. In the sequel, K is a commutative semi-ring, in particular IR or IR + . We are interested in mapping of Σ into K called formal power series. Let f be a formal power series, for a set of words S, we define f (S) = u∈S f (u) ∈ K. Some formal power series can be represented by compact models, called MA. Definition 1 (Multiplicity Automaton). Let K be a semi-ring, a K-multiplicity automaton (K-MA) with n states is a structure Σ, Q, {A o } o∈Σ , α 0 , α ∞ , where Σ is an alphabet and Q is a finite set of states. Matrices A o ∈ K n×n contain the transition weights. The vectors α ∞ ∈ K n and α 0 ∈ K n contain respectively the terminal and initial weights. A K-MA M defines a rational language r M : Σ → IR,

r M (u) = r M (o 1 . . . o k ) = α 0 A u α 0 = α 0 A o1 . . . A o k α ∞ . A function f is said realized by a K-MA M, if r M = f .
Definition 2 (Rational language). A formal power series r is rational over K iff it is realized by a K-MA.

Two K-MA that define the same rational language are equivalent. A K-MA is minimal if there is not an equivalent K-MA with strictly fewer states. An important operation on K-MA is the conjugation by an invertible matrix. Let R ∈ K n×n be invertible and

M = Σ, Q, {A o } o∈Σ , α 0 , α ∞ be a K-MA of dimension n, then M = Σ, Q, R -1 A o R o∈Σ , R α 0 , R -1 α ∞ ,
defines a conjugated K-MA. In addition, M and M are equivalent. In this paper, we are interested in learning languages representing distributions over finite-length words that can be compactly represented by a K-MA. When K = IR, this forms the class of stochastic rational language.

Definition 3 (Stochastic rational language). A stochastic rational language p is a rational language with values in IR + such that u∈Σ p(u) = 1.

Stochastic rational languages are associated to the following subclass of IR-MA.

Definition 4 (Stochastic Multiplicity Automaton).

A stochastic multiplicity automaton (SMA) M is a IR-MA realizing a rational stochastic language.

Hankel Matrix Representation

Let f : Σ → K be a formal power series, we define H f ∈ K Σ ×Σ the bi-infinite Hankel matrix whose rows and columns are indexed by

Σ such that H f [u, v] = f (uv), H f =        ε a b aa . . . ε f (ε) f (a) f (b) f (aa) . . . a f (a) f (aa) f (ab) f (aaa) . . . b f (b) f (ba) f (bb) f (baa) . . . aa f (aa) f (aaa) f (aab) f (aaaa) . . . . . . . . . . . . . . . . . . . . .       
When f is a stochastic language, H f contains occurring probabilities that can be estimated from samples by empirical counts. Details about matrices defined over semirings can be found in [START_REF] Guterman | Rank and determinant functions for matrices over semirings[END_REF]. We note H the Hankel matrix when its formal power series can be inferred from the context. Let for all o ∈ Σ,

H o ∈ K Σ ×Σ , h S ∈ K Σ and h P ∈ K Σ be such that H o (u, v) = f (uov), h S (u) = h P (u) = f (u)
. These vectors and matrices can be extracted from H. The Hankel representation of formal series lies in the heart of all MoM-based learning algorithms, because of the following fundamental theorem.

Theorem 1 (See [START_REF] Carlyle | Realizations by stochastic finite automata[END_REF]). Let r be a rational language over K and M a K-MA with n states that realizes it, then rank(H r ) ≤ n. Conversely, if the Hankel matrix H r of a formal power series r has a finite rank n, then r is a rational language over K and can be realized by a minimal K-MA with exactly n states.

Note that, the original proof assumes that K is a field but remains true when K is a commutative semi-ring, as ranks, determinants and inverses are well defined in semimodules. In addition, the proof gives also the construction of H from a K-MA and vice-versa. For a K-MA with n states, observe that H[u, v] = (α 0 A u )(A v α ∞ ). Let P ∈ K Σ ×n and S ∈ K n×Σ be matrices defined as follows,

P = ((α 0 A u ) ) u∈Σ , S = (A v α ∞ ) v∈Σ ,
then H = P S. Moreover, we have that,

H o = P A o S, h S = α 0 S, h P = P α ∞ . ( 1
)
So the K-MA parameters can be recovered by solving eq. ( 1). Hopefully, we do not need to consider the bi-infinite Hankel matrix to recover the underlying K-MA. Given a basis B = (P, S) of prefixes and suffixes, we denote by H B the sub-block of H.

A basis B is complete if H B has the same rank than H. A basis is suffix-closed if ∀u ∈ Σ , ∀o ∈ Σ, ou ∈ S ⇒ u ∈ S.
In [START_REF] Balle | Learning finite-state machines: algorithmic and statistical aspects[END_REF], the author shows that if B = (P, S) is a suffix-closed complete basis, by defining P over P, S over S and H over B, we can recover a MA using eq. ( 1).

Spectral Learning

This section reviews the Spectral Learning algorithm to learn IR-MA from samples generated by a stochastic rational language. In the literature, severals methods are used to build a suffix-closed basis from data. For example, one can use all prefixes and suffixes that appear in the training set. In addition, we require that sets of prefixes and suffixes contain the empty word ε. Once a basis is chosen, the Spectral algorithm first estimates the probabilities in H B by empirical counts. Then, it recovers a factorized form of H B = U DV through a truncated Singular Value Decomposition (SVD). Finally, setting P = U D and S = V , the algorithm solves eq. ( 1). More precisely, let 1 S ε and 1 P ε be vectors filled with 0s with a single 1 at the index of the empty word in the basis, we have that h

S = (1 P ε ) H B , h P = H B 1 S ε . Let, for all o ∈ Σ , T o ∈ IR |S|×|S| be matrices such that, T o [u, v] = δ u=ov then we have that H o = H B T o .
Using these identities in eq. ( 1), one obtains Algorithm 1. In the experiments, following the advices of [START_REF] Cohen | Experiments with spectral learning of latent-variable pcfgs[END_REF], we normalized the feature-variance of the coefficients of the Hankel matrix by independently scaling each row and column by a factor c u = |S| /(#u + 5), where #u is the number of occurrences of u. In addition, depending on the problem, it can be better to work with other series derived from p. For example, the substring-based series p substring (u) = w,v∈Σ p(wuv) is related to p string . According to [START_REF] Balle | Learning finite-state machines: algorithmic and statistical aspects[END_REF], if p string is realized by a SMA, then p substring is too. In addition, he provides an explicit conversion between string-based SMA and substring-based SMA preserving the number of states. For all algorithms compared in Section 4, we used the series leading to the best results for each problem.

Algorithm 1 Spectral algorithm for IR-MA.

1: Choose a set of prefixes P ⊂ Σ and suffixes S ⊂ Σ both containing ε 2: Using S, build the matrices To for all o ∈ Σ such that To[u, v] = δu=ov 3: Estimate HB by empirical counts. 4: U, D, V = SVDn(HB) using the truncated SVD, where n is a parameter of the algorithm 5: For all o ∈ Σ do Ao = V ToV 6:

α 0 = (1 P ε ) HBV 7: α∞ = V 1 S ε 3 Non-negative Spectral Learning
As mentioned in the introduction, Algorithm 1 is designed to learn IR-MA and is very unlikely to return a SMA. This is a major drawback when a probability distribution is required. In this case, one has to rely on heuristics like thresholding the values to be contained in [0, 1], which introduces errors in predictions. A natural enhancement of the Spectral algorithm would be to constraint the return model to be a SMA. Unfortunately, this is not likely to be feasible due to the underlying complexity between IR-MA and SMA. Indeed, although IR-MA are strictly more general than SMA, checking whether a IR-MA is stochastic is undecidable [START_REF] Denis | On rational stochastic languages[END_REF]. In terms of algorithms, constraining the return model to be a SMA would require adding an infinite number of constraints. As a matter of fact, a IR-MA realizing a language r requires the non-negativity (∀u ∈ Σ , r(x) ≥ 0) and the convergence to 1 (lim L→+∞ L l=0 u∈Σ l r(u) = 1) of the series to be a SMA. The undecidability comes only from the non-negativeness. Note that only the existence of the limit is really required as a convergent series can always be normalized. Thus, we propose to restrict learning to IR + -MA, which by definition produces non-negative series. Although, SMA are not included in IR + -MA, there are IR + -MA realizing probability distributions and called Probabilistic Finite Automata (PFA). Relations between all these classes of MA are summed up in Figure 1.

Definition 5 (Probabilistic (Deterministic) Finite Automaton). A probabilistic finite automaton (PFA

) M = Σ, Q, {A o } o∈Σ , α 0 , α ∞ is a SMA with non-negative weights verifying 1 α 0 = 1, α ∞ + o∈Σ A o 1 = 1.
The weights of PFA are in [0, 1] and can be viewed as probabilities over transitions, initial states and terminal states. A Probabilistic Deterministic Infinite Automata (PDFA) if a PFA with deterministic transitions.

Thus, in contrast to the Spectral algorithm which returns a IR-MA, NNSpectral returns a IR + -MA avoiding the non-negativity probability problem. The NNSpectral algorithm, given in Algorithm 2, also uses the decomposition in Theorem 1, but applied to MA defined on the semi-ring IR + instead of the field IR. Thus, there exists a low-rank factorization with non-negative factors of the Hankel matrix representing a IR + -MA. Finding such a decomposition is a well-known problem, called NMF. It aims at finding a low rank decomposition of a given non-negative data matrix H B ∈ IR n×m such that H B -P S F where P ∈ IR n×r and S ∈ IR r×m are component-wise non-negative and • F is the Frobenius norm. Unfortunately, NMF is NP-Hard in general [START_REF] Vavasis | On the complexity of nonnegative matrix factorization[END_REF] and the decomposition is not unique which makes the problem ill-posed. Hence, in practice, heuristic algorithms are used which have only the guarantee to converge to a stationary point. Most of them run in O(nmr). We refer the reader to [START_REF] Gillis | The Why and How of Nonnegative Matrix Factorization[END_REF] for a comprehensive survey of these algorithms. In the experiments, the alternating non-negative least squares algorithm has shown to be a good trade-off between convergence speed and quality of the approximation. This algorithm iteratively optimizes P and S by solving NNLS problems. NNLS is a non-negative constrained version of the least squares method.

Being equivalent to a quadratic programming problem under linear constraints, it is convex and algorithms converge to the optimum. In the experiments, NNLS problems were solved using the projected gradient algorithm of [START_REF] Lin | Projected gradient methods for nonnegative matrix factorization[END_REF]. So, given a non-negative factorization of H B , solving eq. ( 1) is done by NNLS to ensure the non-negativity of the weights. Although NNSpectral cannot learn SMA (or the equivalent OOMs, PSRs), a SMA can be arbitrarily well approximated by a PFA (or the equivalent HMMs, POMDPs) with a growing number of states [START_REF] Gybels | Some improvements of the spectral learning approach for probabilistic grammatical inference[END_REF]. The Probabilistic Automata learning Competition (PAutomaC) is dealing with the problem of learning probabilistic distributions from strings drawn from finite-state automata. From the 48 problems available, we have selected the same twelve problems than in [START_REF] Balle | Methods of moments for learning stochastic languages: Unified presentation and empirical comparison[END_REF], to provide a fair comparison with other algorithms. The generating model can be of three kinds: PFA, HMMs or PDFA. Four models have been selected from each class. A detailed description of each problem can be found in [START_REF] Verwer | Results of the pautomac probabilistic automaton learning competition[END_REF]. Table 1 compares the best results of NNSpectral between learning from strings or substrings to EM and the best results among the following MoM-based algorithms : CO [START_REF] Balle | Local loss optimization in operator models: A new insight into spectral learning[END_REF] using strings, Tensor [1] using strings and Spectral using strings and substrings. A description and comparison of these algorithms can be found in [START_REF] Balle | Methods of moments for learning stochastic languages: Unified presentation and empirical comparison[END_REF].

The quality of a model can be measured by the quality of the probability distribution it realizes. The objective is to learn a MA realizing a series p close to the distribution p , which generated the training set T . The quality of p is measured by the perplexity that corresponds to the average number of bits needed to represent a word using the optimal code given by p .

Perplexity(M) = 2 -u∈T p (u) log(p M (u))
The quality of model can also be evaluated by the Word Error Rate (WER) metric. It measures the fraction of incorrectly one-step-ahead predicted symbols.

In the simulations, we perform a grid search to find the optimal rank for each of the performance metric. For each problem, the best algorithm is indicated by a bold number and the best score between NNSpectral and other MoM-based algorithms is underlined. The score of the true model is reported for comparison. For the perplexity, NNSpectral outperforms other MoM-based algorithm on the 12 problems and does better that EM for 7 problems. For the other 5 problems, NNSpectral achieves performances very close to EM and to the true model. For the WER metric, NNSpectral beats other MoM-based algorithms on 10 problems and scores at top on 7 problems. Also, for all problems NNSpectral produces perplexity close to the optimal ones, whereas, on 5 problems EM fails to produce an acceptable solution.

Conclusion

In this paper, we proposed a new algorithm inspired by the theoretical developments of MA defined on semi-rings. NNSpectral works by constraining the search space to ensure non-negativity of the rational language. Like other MoM-based algorithms, NNSpectral is able to handle large-scale problems much faster than EM. Here, the consistency is lost due to the use of heuristics to solve the NMF problem. Experimentally, this does not seem to be a major problem because NNSpectral outperforms other MoMbased algorithms and is competitive with EM that sometimes produces aberrant solutions. Thus, NNSpectral provides a good alternative to the EM algorithm with a much lower computational cost. In further works, we would like to investigate how using NNSpectral to initialize an EM algorithm and how adding constraints in NNSpectral to ensure the convergence of the series. In addition, relations with NMF could be further exploit to produce tensor [1] or kernel-based algorithms.

Fig. 1 .

 1 Fig. 1. Solid lines are polynomially decidable. Dashed line are undecidable. See [2][9].

Algorithm 2

 2 NNSpectral algorithm for IR + -MA 1: Choose a set of prefixes P ⊂ Σ and suffixes S ⊂ Σ both containing ε 2: Using S, build the matrices To for all o ∈ Σ such that To[u, v] = δu=ov 3: Estimate HB by empirical counts. 4: P, S ← argmin P S HB -P S F s.t. P ∈ IR |P|×n , S ∈ IR n×|S| ≥ 0 by NMF 5: For all o ∈ Σ do Ao ← argmin A≥0 AS -STo F by NNLS 6: α0 ← argmin α≥0 α S -(1 P ε ) HB F by NNLS 7: α∞ ← S1 S ε 4 Numerical Experiments