Aurélien Appriou

Andrzej Cichocki

Fabien Lotte

Aurelien Appriou

Modern machine learning algorithms

Keywords: Brain-Computer Interfaces (BCI), Mental Workload, Emotions, Riemannian geometry, Deep Learning, EEG

come

I. INTRODUCTION

B RAIN-Computer Interfaces (BCIs) enable their users to interact with computers by using brain activity only, usually measured with Electroencephalography (EEG) [START_REF] Clerc | Brain-Computer Interfaces 1[END_REF]. For example, BCIs can enable people with severe motor impairment to control a wheelchair with EEG only, e.g., by imagining left or right hand movements to make the wheelchair turn left or right [START_REF] Millán | Combining Brain-Computer Interfaces and Assistive Technologies: State-of-the-Art and Challenges[END_REF]. Such BCIs are called active BCIs since users are actively sending commands to the system [START_REF] Zander | Towards passive brain-computer interfaces: applying brain-computer interface technology to human-machine systems in general[END_REF]. In contrast, the so-called passive BCIs [START_REF] Zander | Towards passive brain-computer interfaces: applying brain-computer interface technology to human-machine systems in general[END_REF] are not used for direct control, but to monitor users' mental states in real-time, to adapt an application accordingly. For instance, passive BCIs were used to estimate mental workload [START_REF] Mühl | EEG-based workload estimation across affective contexts[END_REF], i.e., the amount of cognitive resources currently engaged by subjects, or affective states [START_REF] Muhl | A survey of affective brain computer interfaces: principles, state-of-the-art, and challenges[END_REF], i.e., the emotions subjects currently feel.

A. Mental workload estimation from EEG

Passive BCIs were used to study mental workload during navigation tasks with different input devices [START_REF] Frey | Framework for electroencephalography-based evaluation of user experience[END_REF], during visualization tasks [START_REF] Peck | Using fNIRS brain sensing to evaluate information visualization interfaces[END_REF] or plane piloting [START_REF] Gateau | Real-time state estimation in a flight simulator using {fNIRS}[END_REF]. Workload estimation was also used to design applications that dynamically adapt to the users' states, e.g., video games with adaptive difficulty [START_REF] Fairclough | BCI and physiological computing for computer games: Differences, similarities & intuitive control[END_REF] or training applications with a sequence of exercises adapted to the cognitive capabilities of each learner [START_REF] Yuksel | Learn piano with BACh: An adaptive learning interface that adjusts task difficulty based on brain state[END_REF]. In general, estimating mental workload from EEG is extensively used in the field of Neuroergonomics [START_REF] Ayaz | Neuroergonomics: The Brain at Work and in Everyday Life[END_REF], [START_REF] Frey | EEG-based neuroergonomics for 3D user interfaces: opportunities and challenges[END_REF] which consists in using neuroscience tools and results, e.g., BCIs, to assess the ergonomic qualities of interactive systems.

However, reliably estimating mental workload from EEG signals, over time, contexts and subjects is difficult [START_REF] Mühl | EEG-based workload estimation across affective contexts[END_REF]. For instance, in [START_REF] Mühl | EEG-based workload estimation across affective contexts[END_REF], discriminating low from high workload in 2s epoch of oscillatory EEG activity was possible with a classification accuracy of only about 69%, using a Filter Bank Common Spatial Patterns (FBCSP) algorithm [START_REF] Ang | Filter bank common spatial pattern algorithm on BCI competition IV datasets 2a and 2b[END_REF] coupled with a Linear Discriminant Analysis (LDA) [START_REF] Blankertz | Optimizing spatial filters for robust EEG single-trial analysis[END_REF]. In [START_REF] Mahmoud | Novel Classification System for Classifying Cognitive Workload Levels under Vague Visual Stimulation[END_REF], authors used a Bilinear CSP and a Linear Probabilistic Support Vector Machine (SVM) to classify 1.1s epochs of oscillatory EEG activity into four levels and obtained 93% classification accuracy. In [START_REF] Brouwer | Estimating workload using {EEG} spectral power and {ERPs} in the n-back task[END_REF], authors applied a SVM to classify two levels of workload, using N-back tasks. They obtained 84% classification accuracy on 0.5 to 1.5s long epochs. In [START_REF] Zarjam | Estimating cognitive workload using wavelet entropy-based features during an arithmetic task[END_REF], authors extracted features using a wavelet entropy, and then applied a Multi-Layer Perceptron (MLP) to classify workload data into 7 levels with 5s long epochs. They obtained 98% classification accuracy on their own data set, and 83% on the data set from [START_REF] Grimes | Feasibility and Pragmatics of Classifying Working Memory Load with an Electroencephalograph[END_REF]. Except the last study above, all classification accuracies have been obtained in offline analysis settings.

B. Affective states estimation from EEG

Passive BCIs estimating affective states are called affective BCIs (aBCIs) [START_REF] Muhl | A survey of affective brain computer interfaces: principles, state-of-the-art, and challenges[END_REF]. Examples of aBCI work include studies that found neurophysiological responses to differentiate between frustration and boredom in e-learning [START_REF] Chanel | Emotion Assessment From Physiological Signals for Adaptation of Game Difficulty[END_REF] and between frustration and normal game play [START_REF] Reuderink | arousal and dominance in the EEG during game play[END_REF]. Moreover, Rani et al. [START_REF] Rani | Maintaining optimal challenge in computer games through real-time physiological feedback[END_REF] showed that both players' enjoyment and skills increased when tasks were adapted to their affective states rather than to their performances. In [START_REF] Fairclough | Fundamentals of physiological computing[END_REF], players could modulate their affective states to influence games' parameters. Affective BCIs can also be used for automatic media recommendation [START_REF] Soleymani | A multimodal database for affect recognition and implicit tagging[END_REF], or real-world emotions detection [START_REF] Zheng | Emotionmeter: A multimodal framework for recognizing human emotions[END_REF].

Despite an increasing number of aBCIs studies, defining and clustering emotion dimensions remains challenging. There are multiple main approaches to define emotion classes [START_REF] Muhl | A survey of affective brain computer interfaces: principles, state-of-the-art, and challenges[END_REF]: the most popular one, and the one used in our study, is the circumplex model of Russell [START_REF] Russell | A circumplex model of affect[END_REF], which assumes that any affective state can be localized on a two-dimensional plane. The first axis of this plane is the valence, ranging from positive feelings to negative ones, and the second axis represents arousal, ranging from calm to excited. In [START_REF] Soleymani | A multimodal database for affect recognition and implicit tagging[END_REF], they used an SVM for a 3-class problem with 15s long epochs, for valence (low vs neutral vs high, 50% accuracy) and arousal (low vs neutral vs high, 62%). [START_REF] Candra | EEG emotion recognition using reduced channel wavelet entropy and average wavelet coefficient features with normal Mutual Information method[END_REF] studied valence and arousal (low vs high) as well, but with 6s EEG epochs. They respectively obtained 77% and 74% classification accuracy by using a combination of wavelet entropy and average wavelet coefficient coupled with an SVM. [START_REF] Al-Nafjan | Recognition of affective states via electroencephalogram analysis and classification[END_REF] used SVMs, Deep Neural Networks or Random Forest to classify data from the DEAP data set [START_REF] Koelstra | DEAP: A database for emotion analysis Using physiological signals[END_REF] (presented in the Methods section of this paper), using a 2-class problem with 60s long epochs. They respectively obtained 79%, 49% and 56% classification accuracy. However, the Deep Learning method was not described in this paper, making it difficult to assess the validity and superiority of this approach. [START_REF] Kothe | Emotion recognition from EEG during self-paced emotional imagery[END_REF] used a Logistic Regression (LR) to discriminate valence levels (low vs high), and obtained 71% accuracy on 6s epochs. In [START_REF] Koelstra | DEAP: A database for emotion analysis Using physiological signals[END_REF], which introduced the DEAP data set we use here, they applied a Naive Bayes Classifier for 2-class discrimination (low vs high) for both valence and arousal dimensions (60 s epochs). They obtained 57% and 62% of accuracy for the valence and arousal dimensions, respectively. Other studies used Deep Learning methods to optimize affective states classification (see supplementary material).

Most studies proved that classifying affective states from EEG remains really challenging since results hardly go over chance level accuracy. Some studies were even unable to obtain better than chance results when reproducing previous works with statistically rigorous evaluation methods [START_REF] Lotte | Turning negative into positives! Exploiting "negative" results in Brain-Machine Interface (BMI) research[END_REF]. Finally, confounding factors due to electromyography (e.g., facial muscles activity during emotion expression) have likely played a role in the performances obtained in many studies. However, other studies have obtained better accuracies when using different EEG patterns. For instance, [START_REF] Frantzidis | Toward emotion aware computing: An integrated approach using multichannel neurophysiological recordings and affective visual stimuli[END_REF] obtained 81% classification accuracy using evoked potentials for a 4-class valence/arousal classification. Note that this paper focuses on oscillatory activity, as this can classify cognitive and affective states from continuous EEG, without the need for stimuluslocked response, which evoked potentials do need.

C. Paper objectives

Thus, the classification accuracies obtained so far -mostly around 70% for workload, and around 60-65% for emotions in oscillatory-based studies -revealed the need for more robust and accurate EEG classification algorithms. Therefore, we propose here to study algorithms that proved efficient either in recent active BCI classification competitions [START_REF] Ang | Filter bank common spatial pattern algorithm on BCI competition IV datasets 2a and 2b[END_REF], [START_REF] Yger | Riemannian approaches in Brain-Computer Interfaces: a review[END_REF], notably Riemannian geometry classifiers, or in other fields of artificial intelligence, such as Deep Learning [START_REF] Schirrmeister | Deep learning with convolutional neural networks for EEG decoding and visualization[END_REF], [START_REF] Lecun | Deep learning[END_REF]. Note that such algorithms have been mostly explored for EEG classification of motor tasks, but not systematically studied and compared for workload/affective states estimation. Here we formally study and compare these various algorithms as well as two new variants we proposed, for both workload, arousal and valence classification from EEG signals 1 . We also propose guidelines about which algorithm to use in which context. As baseline, we use two standard methods for studying workload levels/affective states classification: 1) Common Spatial Pattern (CSP) spatial filters with an LDA classifier and 2) the FBCSP [START_REF] Ang | Filter bank common spatial pattern algorithm on BCI competition IV datasets 2a and 2b[END_REF], which is a CSP extension that won numerous active BCI competitions. Then, we studied two Riemannian approaches, manipulating and classifying EEG signals as covariance matrices: Minimum Distance to the Mean with Fisher geodesic filtering classifier (FgMDM) and Tangent Space Classifier (TSC). Such methods have recently won six international brain signals competitions [START_REF] Yger | Riemannian approaches in Brain-Computer Interfaces: a review[END_REF]. We then propose to improve these Riemannian approaches by working on a bank of band-pass filters such as the ones used for FBCSP, instead of using a unique band-pass filter. We name these new approaches FBFgMDM and FBTSC. Finally, we used a Convolutional Neural Network (CNN), i.e., a Deep Learning algorithm, which recently obtained promising results for many machine learning problems [START_REF] Lecun | Deep learning[END_REF]. We studied the CNN developed in [START_REF] Schirrmeister | Deep learning with convolutional neural networks for EEG decoding and visualization[END_REF], since it obtained promising results for motor imagery-based BCIs.

In this paper, we first present the workload and emotion EEG data sets used, before describing each machine learning algorithm. We perform two evaluation studies: 1) a subjectspecific study, with each algorithm trained on data specific to each subject, and then tested on other data from the same subject. This is the standard way current BCIs are designed, given the large between-subject variability [START_REF] Blankertz | Optimizing spatial filters for robust EEG single-trial analysis[END_REF]; 2) a subjectindependent study, with each algorithm trained on all data recorded from all subjects except that of the target subject, on which algorithms are tested. This is much more challenging, but if successful, would enable BCI-based monitoring without requiring any calibration for new subjects.

II. METHODS

A. Mental workload EEG data set

The data set used comes from [START_REF] Mühl | EEG-based workload estimation across affective contexts[END_REF]. Signals from 28 EEG electrodes (active electrodes in a 10/20 system without T7, T8, Fp1, and Fp2) were recorded from 22 users [START_REF] Mühl | EEG-based workload estimation across affective contexts[END_REF]. To induce mental workload variations, N-back tasks were used: the user had to indicate whether a letter displayed on screen was the same one as the letter displayed N letters before, in a stream of successively displayed letters. Here, 2-sec trials from a 0back task were labeled as "low" workload, while those from a 2-back were labelled as "high" workload. In total, 720 trials were available for each workload level and user. See the supplementary material for more information.

As introduced previously, we studied both subject-specific and subject-independent calibrations. For subject-specific calibration, the first half of each user's trials was used for training and the second half for testing. For the subject-independent calibration, the training set comprised all trials of all users except the current user used for testing, i.e., around 21*1440 = 30240 training trials. To allow for comparisons between calibration types, the testing set of each user was the same testing set as with subject-specific calibration, i.e., the second half of the trials (720 testing trials) from this user.

B. Emotion EEG data set

The data set used for studying emotions was the "DEAP" database [START_REF] Koelstra | DEAP: A database for emotion analysis Using physiological signals[END_REF]. It used music-video clips to influence two types of emotion dimensions -valence and arousal, according to the circumplex model of Russell [START_REF] Russell | A circumplex model of affect[END_REF]. The data set contains 40 trials, corresponding to 40 music video-clips, recorded on 32 participants. EEG were recorded using 32 electrodes (placed according to the international 10-20 system). Valence and arousal levels were measured using Russell's valencearousal scale directly after each videos, by clicking on a 1-9 continuous scale. This self-assessment system on a continuous scale makes the classes definition more complex: in DEAP [START_REF] Koelstra | DEAP: A database for emotion analysis Using physiological signals[END_REF] as well was in our study, 5 was kept as a threshold to split trials into two classes -low and high -for both "emotion-arousal" and "emotion-valence" data sets, making classes unbalanced. All the classifiers used were able to deal with unbalanced classes, except the CNN for which we upsampled the minority class by randomly duplicating trials from this class in order to obtain balanced classes.

For the subject-specific study, given the low number of trials, we performed a "leave-one-out" cross-validation. Thus, we used 40 models for each subject, each model being trained on 39 trials and tested on 1 trial. For the subject-independent study, we kept all trials of all subjects to compose the training set, except the current subject used for testing (i.e., 31*40 = 1240 trials for the training set). The testing set of each subject was composed of all trials of this subject, i.e., 40 trials.

C. Machine learning algorithms explored

The existing algorithms we evaluate here were all studied on EEG-based motor imagery classification, a widely used BCI design, and obtained impressive results. Since both motor imagery, workload and emotions lead to change in EEG oscillatory activity, it is likely that methods that proved effective for motor imagery can prove effective for workload or emotion classification as well. However, to the best of our knowledge, such methods have never been tested and compared together neither on workload nor on emotions data sets nor with subject-independent calibration. We thus propose this evaluation in this paper. We also propose some new variants of some of these algorithms. Altogether, we studied 7 algorithms. First, CSP and LDA were used as a baseline since they are widely used by the BCI community [START_REF] Blankertz | Optimizing spatial filters for robust EEG single-trial analysis[END_REF]. We then explored the FFBCSP and LDA [START_REF] Ang | Filter bank common spatial pattern algorithm on BCI competition IV datasets 2a and 2b[END_REF], a CNN [START_REF] Schirrmeister | Deep learning with convolutional neural networks for EEG decoding and visualization[END_REF], and four different methods based on Riemannian geometry: two existing ones, namely the Fisher geodesic Minimum Distance to the Mean classifier (FgMDM) and the Tangent Space Classifier (TSC) [START_REF] Yger | Riemannian approaches in Brain-Computer Interfaces: a review[END_REF], and two new extensions we propose here to better exploit the spectral information, namely the Filter Bank FgMDM and Filter Bank TSC. For the workload data set, we assess performances using classification accuracy, i.e., the percentage of test trials correctly classified. For the emotion data set, we used balanced accuracy, i.e. the average of recall obtained on each class, since the classes were unbalanced.

1) Common Spatial Patterns (CSP)

CSP is a widely used algorithm for binary EEG classification, for oscillatory activity-based BCI. It has been shown that changes in both workload [START_REF] Brouwer | Estimating workload using {EEG} spectral power and {ERPs} in the n-back task[END_REF] and emotions [START_REF] Muhl | A survey of affective brain computer interfaces: principles, state-of-the-art, and challenges[END_REF] induce changes in EEG oscillatory activity.The CSP algorithm optimizes spatial filters, i.e., a linear combination of the original EEG signals. It is done such that the variance of a spatially filtered signal, i.e. the band power of this signal, is maximized for one class and minimized for the other class. Formally, CSP optimizes spatial filter w by either maximizing or minimizing:

J CSP (w) = wX 1 X T 1 w T wX 2 X T 2 w T = wC 1 w T wC 2 w T (1)
where T denotes transpose, X i is the band-pass filtered training signal matrix for class i (with the samples as columns and the channels as rows) and C i the spatial covariance matrix from class i. In practice, the covariance matrix C i is defined as the average covariance matrix of each trial from class i [START_REF] Blankertz | Optimizing spatial filters for robust EEG single-trial analysis[END_REF]. The spatial filters w that maximize or minimize J CSP (w) are the eigenvectors corresponding to the largest and lowest eigenvalues, respectively, of the Generalized Eigen Value Decomposition of matrices C 1 and C 2 . In this study, we used six filters, corresponding to the three largest and three lowest eigenvalues, as recommended in [START_REF] Blankertz | Optimizing spatial filters for robust EEG single-trial analysis[END_REF]. Once these filters are obtained, we use as CSP features f = log(wXX T w T ), i.e., the band power of the spatially filtered signals. We used these features as input to an LDA classifier. The CSP requires EEG signals to be band-pass filtered in a specific narrow frequency band. The Alpha rhythm (8-12Hz) being known to vary according to both workload [START_REF] Mühl | EEG-based workload estimation across affective contexts[END_REF] and emotions [START_REF] Koelstra | DEAP: A database for emotion analysis Using physiological signals[END_REF], we applied CSP after band-pass filtering in 8-12 Hz.

2) Filter Bank Common Spatial Patterns (FBCSP)

The FBCSP is an algorithm that optimizes both spatial and spectral filters. To do so, FBCSP first filters EEG signals into multiple frequency bands using a filter bank. Here we used nine band-pass filters in 4Hz-wide bands (in 4-8 Hz, 8-12 Hz, . . . , 36-40 Hz) as in [START_REF] Ang | Filter bank common spatial pattern algorithm on BCI competition IV datasets 2a and 2b[END_REF]. Then, for each bandpassed signals, CSP is used to optimize two spatial filter pairs. From the resulting 36 features (9 bands × 4 CSP filters per band), the four most relevant ones were selected using minimal Redundancy Maximal Relevance (mRMR) [START_REF] Peng | Feature Selection Based On Mutual Information: Criteria of Max-Dependency,Max-Relevance, and Min-Redundancy[END_REF], and used as input to an LDA. The FBCSP algorithm proved its efficiency when winning the Fifth International BCI competition [START_REF] Ang | Filter bank common spatial pattern algorithm on BCI competition IV datasets 2a and 2b[END_REF].

3) Riemannian Geometry Riemannian approaches represent EEG trials as covariance matrices, which are symmetric positive definite (SPD) matrices, and manipulate them with an appropriate geometry, the Riemannian geometry [START_REF] Yger | Riemannian approaches in Brain-Computer Interfaces: a review[END_REF], [START_REF] Congedo | Riemannian geometry for eeg-based brain-computer interfaces; a primer and a review[END_REF]. Classifiers based on such geometry are called Riemannian Geometry Classifiers (RGC).

First, in a Riemannian manif old we can estimate intrinsic non-Euclidean distances between two SPD matrices, i.e. two points (here C 1 and C 2 ), using the Riemannian distance:

δ 2 (C 1 , C 2 ) = n log 2 λ n (C -1 1 C 2 ), (2)
where λ n (M) is the n th eigenvalue of matrix M. The set of tangent vectors to point G on the manif old defines the manifold tangent space at G. More generally, any SPD matrix C i can be projected onto the tangent space at point G using:

S i = Log G (C i ) = G 1/2 logm(G -1/2 C i G -1/2 )G 1/2 , (3)
S i being the projection of C i onto the tangent plane, and logm(•) denotes the logarithm of a matrix.

In this paper, we studied two existing RGCs -the Mean classifier (FgMDM) and the Tangent Space classifier (TSC)and introduced two new ones -the Filter Bank TSC (FBTSC) and the Filter Bank FgMDM (FBFgMDM).

Existing methods: FgMDM [START_REF] Barachant | Riemannian geometry applied to BCI classification[END_REF]: FgMDM projects training matrices C i onto the tangent space at point G (the mean of all training data) using Eq. ( 3), to obtain matrices S i . Then, a Fisher geodesic filter is obtained by optimizing an LDA classifier on S vec i , the vectorized upper-triangular elements of S i , to discriminate classes using such vectors. This results in a matrix of weights W = LDA(S vec i ). The projected SPD matrices S i from both the training & the testing sets are then filtered with weights W, using Ŝi = W(W T W) -1 W T S vec i , where Ŝi denotes the geodesic filtered SPD matrices from S i . Then, these filtered matrices Ŝi are projected back onto the manifold using equation:

Ĉi = Exp C ( Ŝi ) = G 1/2 expm(G -1/2 Ŝi G -1/2 )G 1/2 , (4)
where Ĉi are the filtered SPD matrices projected onto the manifold and expm(M) denotes the exponential of matrix M. Finally, this approach uses a Minimum Distance to the Mean classifier to classify testing geodesic filtered matrices Ĉi . To do so, during the training step, the class centroids G k of each class k are computed by averaging the geodesic filtered covariance matrices Ĉk i from each class k. During testing, the Riemannian distances between the testing geodesic filtered matrix Ĉj and each class centroid G k are first calculated, using Eq. ( 2). The matrix Ĉj is assigned class label k for which the centroid G k is the nearest. In our study, FgMDM was applied on EEG band-pass filtered in 8-12Hz, as for the CSP.

TSC: TSC first projects all training SPD matrices C i onto the tangent space at point G (the mean of all training matrices). Then, it uses any classifier such as LDA, SVM or Logistic Regression (LR) on the vectorized upper-triangular elements of the projected matrices [START_REF] Barachant | Multiclass brain-computer interface classification by Riemannian geometry[END_REF]. We used LR with L 2 regularization (with the default C = 1.0 in scikitlearn [START_REF] Pedregosa | Scikit-learn: Machine Learning in Python[END_REF]), As for FgMDM, TSC used data filtered in 8-12Hz.

New methods:

Filter Bank FgMDM (FBFgMDM): Contrary to FgMDM which exploits EEG signals in a single frequency band, this method applies FgMDM in multiple bands separately, and combines the resulting distances to exploit additional spectral information. This should possibly improve classification performances, as FBCSP did to improve CSP. To do so, FBFgMDM first filters EEG signals in multiple bands using a filter bank, as for FBCSP. Here we used the same bands as the FBCSP, i.e., 4-8 Hz, 8-12 Hz, . . . , 36-40 Hz. Then for the EEG signals in each frequency band j, this method first uses a regular FgMDM, i.e., it computes the Riemannian distances δ 2 (G kj , Ĉij ) between a geodesic filtered SPD matrix Ĉij and each class centroid G kj . Then, from all nine bands j, the four most useful ones for classification are selected with mRMR feature selection [START_REF] Peng | Feature Selection Based On Mutual Information: Criteria of Max-Dependency,Max-Relevance, and Min-Redundancy[END_REF] on the Riemannian distances δ 2 (G kj , Ĉij ) used as features, on the training set. For testing, we compute the squared Riemannian distances for the four bands selected using mRMR only and sum them:

γ 2 (G k , Ĉi ) = j∈Ω δ 2 (G kj , Ĉij ), ( 5
)
where Ω is the set of frequency bands selected with mRMR.

We thus obtain k new distances γ 2 (G k , Ĉi ) to each class k for each trial i. The classification prediction results in choosing the class y i for which the summed squared distance to the centroid is the smallest, i.e.,

y i = argmin k (γ 2 (G k , Ĉi )).
Filter Bank TSC (FBTSC): FBTSC also exploits more spectral information than TSC, by using a filter bank. FBTSC indeed projects matrices C ij , band-pass filtered in bands 4-8 Hz, 8-12 Hz, . . . , 36-40 Hz, to the tangent space using Eq. ( 3). Then, the probabilities that the vectorized upper-triangular elements of the projected SPD matrix S ij belongs to class k is calculated using standard classification algorithms with probabilistic outputs, such as LDA or LR. Here we used LR that directly provides such probability with its softmax function. Since we did so for nine frequency bands, in two classes k, we ended up with nine pairs of probabilities. From these pairs of probabilities, the four most relevant are selected using mRMR on the training set. Finally, we multiplied the probabilities associated to each class k, for the selected bands only, to end up with two probabilities, using P ki = j∈Ω P kij , where P ki is the probability of trial i to be part of class k, and P kij the probability of a projected SPD matrix S ij , band-pass filtered in frequency band j, to be part of class k. The classification prediction results in choosing the class y i for which P ki is the highest, i.e., y i = argmax k (P ki ).

4) Convolutional Neural Networks (CNN)

Shortly, a CNN is a feedforward neural network with at least one convolutional layer. This type of network flows information uni-directionally from the input to the hidden layers and finally to the output. A recent study presented a new type of CNN dedicated to motor task classification in EEG: the Shallow ConvNet [START_REF] Schirrmeister | Deep learning with convolutional neural networks for EEG decoding and visualization[END_REF]. The shallow ConvNet architecture consists in a 3-layer CNN with parameters that have been experimentally tested and validated by their authors [START_REF] Schirrmeister | Deep learning with convolutional neural networks for EEG decoding and visualization[END_REF]. The first layer is a convolutional layer along the temporal dimension, while the subsequent one is a convolutional layer along the spatial dimension, i.e., over EEG electrodes. The first temporal convolution aims at optimizing band-pass filters, and the spatial convolution aims at optimizing spatial filters. Then, signals are squared, a mean pooling is performed (to compute signals band power) and the CNN ends by a fully connected linear classification layer. Overall this CNN thus processes EEG data similarly to the FBCSP and LDA. In contrast to FBCSP, all these filters are optimized simultaneously though, which made it outperform the FBCSP on motor EEG signals [START_REF] Schirrmeister | Deep learning with convolutional neural networks for EEG decoding and visualization[END_REF]. The Shallow ConvNet uses minimally preprocessed EEG signals as input, so we filtered them in 4-40 Hz.

III. RESULTS

Figure 1 summarizes the mean performance obtained on each data set. As a reference, the statistical chance levels [START_REF] Combrisson | Exceeding chance level by chance: The caveat of theoretical chance levels in brain signal classification and statistical assessment of decoding accuracy[END_REF] were estimated at 50.47% for the mental workload study (1440 trials and 22 subjects) and 52.27% for the affective state study (40 trials and 32 subjects). Note that for statistical tests (ANOVA), we checked the data sphericity, and used Greenhouse-Geisser (GG) correction in ANOVA if needed. Post-hocs analyses -Student t-test for paired samples -with Bonferroni adjustments showed no significant differences between algorithms in the subject-specific or subject-independent studies. However, performances obtained suggested better (but non-significantly so) results with the CNN compared to other algorithms, in both subject-specific and subject-independent studies. Riemannian geometry classifiers (RGC), in particular the newly proposed ones (FBfgMDM and FBTSC) provided the second best performances, just after the CNN. On the other hand, the baseline CSP+LDA obtained the worst results.

B. Valence

The balanced classification accuracies obtained are reported on Figure 3. We ran a 2-ways ANOVA for repeated measures to evaluate the impact of Algorithm on the emotionvalence data set, regarding the Calibration Type. The results showed significant differences in Algorithm [GG [START_REF] Clerc | Brain-Computer Interfaces 1[END_REF][START_REF] Yger | Riemannian approaches in Brain-Computer Interfaces: a review[END_REF] Post-hoc analyses -Student t-test for paired sampleswith Bonferroni corrections showed a significant difference between FBTSC and CNN for subject-specific calibration [perf F BT SC = 61.09%, perf CN N = 46.32%; p ≤ 0.05]. No algorithm showed better results than others with the subjectindependent calibration. Overall, FBFgMDM and FBTSC obtained the best accuracy (both about 61%) for subject-specific calibration, while FBCSP obtained the best performances for the subject-independent one (55.2%).

C. Arousal

The balanced classification accuracies for the emotionarousal data set are reported on Fig. Post-hoc analyses -Student t-test for paired samples-with Bonferroni corrections showed significant differences with the subject-specific calibration between CNN and all other classifiers (see results in the supplementary material). No algorithm showed better results than others with the subject-independent calibration. Overall the best results were all obtained by RGCs, FBFgMDM and FBTSC for the subject-specific calibration, and FgMDM for the subject-independent one. The first results to highlight are the CNN classification performances we obtained across the different conditions and data sets. Indeed, this algorithm has a higher mean accuracy (although non-significantly so) than the original authors' results, the baseline CSP+LDA, and more importantly than both FBCSP and Riemannian methods, with both subjectspecific and subject-independent calibrations on the workload data set. Moreover, obtaining reasonable performances in a subject-independent calibration from only two seconds of EEG data and only 21 users for calibration, makes the CNN particularly interesting to design calibration-free neuroadaptive technologies in the future. By contrast, this algorithm significantly under-performed with both subject-specific and subjectindependent calibrations on both the valence and arousal data sets. All algorithms indeed outperformed this CNN in all conditions on the emotion data sets.

Multiple factors could explain the observed algorithm performances. First, the number of trials that are used for training models is important. In [START_REF] Schirrmeister | Deep learning with convolutional neural networks for EEG decoding and visualization[END_REF], authors tested the Shallow Con-vNet on multiple motor-imagery data sets (from 288 to 1168 trials), and often obtained significantly better performances with the CNN than with FBCSP. In our study, the workload data set contained 720 training trials whereas both valence and arousal data sets contained 39 training trials only (with cross validation calibration). This might suggest that the CNN could be useful for mental state classification, but only when large amount of training trials are available (around 700 in our study), which is not always possible. However, other factors also differ between both data sets studied and could also explain differences in CNN performances, including the EEG epochs length (2s epochs for workload and 60s epochs for emotions), and the nature of the mental states studied (workload vs emotions). Indeed, emotions are thought to originate from deep brain areas [START_REF] Muhl | A survey of affective brain computer interfaces: principles, state-of-the-art, and challenges[END_REF] and are thus known as being difficult to estimate from EEG. In the future, deeper analyses would thus be needed to fully disentangle these factors, by systematically varying the mental states studied, the EEG epoch length and the number of training trials.

Another relevant result is the promising classification performances of the proposed RGCs. Indeed, FBTSC and FBFg-MDM outperformed the results from the data sets' authors in most conditions/data sets. Moreover, FBFgMDM with subjectspecific calibration, and FBFgMDM and FBTSC with subjectindependent calibration, reached higher mean accuracies than all other algorithms, except the CNN on the workload data set. More interestingly, the low number of trials in the emotion data sets did not seem to affect their performances since they also reached the highest mean accuracies on both the emotion-valence and emotion-arousal data set, both with subject-specific calibration. These promising results compared to standard RGCs (TSC and FgMDM), are probably due to the extra spectral information extracted with the filter bank, and our study enabled us to quantify this gain.

Finally, FBCSP+LDA obtained a higher mean accuracy than CSP+LDA, although not significantly so, in all conditions/data sets, and the higher overall mean accuracy for valence classification with subject-independent calibration. However, it did not obtain higher mean accuracies than others in any other condition. It should be noted that such results reflect the performances obtained in offline evaluation. As such they are likely to be similar to performances obtained in offline or open-loop mental state monitoring, e.g., for noeuroergonomics (ex: mental workload monitoring) or neuromarketing (ex: emotion monitoring). The performances are likely to change in closed-loop applications, with neuroadaptive technologies, and will thus need to be evaluated in this context as well.

Such results enable us to suggest guidelines about which algorithm to use for mental states classification from EEG. First, the CNN is recommended for mental workload classification with both subject-specific and subject-independent calibration, but seems to need a large amount of training trials (at least several hundreds). It should thus probably be avoided for data sets with little training data (i.e., a few dozens). Second, Filter Bank RGCs (FBTSC and FBFgMFM) should also be recommended to obtain good classification performances notably with subject-specific calibration, for both workload and emotion classification, whatever the amount of training data. However, such methods do not seem suitable for subject-independent classification with little training data and/or for emotion classification. They seem suitable for subject-independent classification of workload with large amount of training data though. Our results also confirmed that passive BCIs with subject-independent calibration is possible but very challenging and with much lower accuracies. Similarly, affective state classification in EEG is possible but much more challenging than workload estimation. However, those suggestions imply computational costs that will differ from an algorithm to another. Indeed, using the FB RGCs or the CNN will require a long calibration time, when the testing phase might also be time consuming and has to be considered before to go towards online uses. See the supplementary material for more information.

For the emotion data set, we labelled trials as in the original paper to allow comparisons, i.e., with a global subjectindependent partition between low/high valence/arousal trials, based on the SAM ratings. Note that better methods for partitioning low/high trials in a per-subject basis can also be used [START_REF] Clerico | Electroencephalography amplitude modulation analysis for automated affective tagging of music video clips[END_REF] in the future, to limit the class imbalance. Still in the future, other deep learning architectures, notably Recurrent Neural Networks (RNN) [START_REF] Lecun | Deep learning[END_REF] may prove promising for EEG classification and passive BCIs as well. It would also be interesting to study whether CNN and RGCs can be used to estimate robustly other cognitive states such as fatigue, curiosity or engagement, and how well the proposed RGCs perform on motor imagery data for active BCIs. Altogether, our results suggested that CNN and the proposed filter bank RGCs are valuable machine learning tools for scientists aiming at decoding cognitive and affective states from EEG signals.

Fig. 1 .

 1 Fig. 1. Mean classification accuracy for each algorithm. The best performance of each study is in green, the worst in red.

Fig. 2 .

 2 Fig. 2. Classification accuracy of each algorithm on the workload data set.

  =6.918, p=0.002], Calibration Type [F(1,32)=21.732, p¡0.0001] and Calibration Type*Algorithm [GG(1,32)=5.374, p=0.003].

Fig. 3 .

 3 Fig. 3. Balanced classification accuracy on the emotion-valence data set.

4 .

 4 We then performed a 2way ANOVA with repeated measures, with factor Algorithms and Calibration Type. Results revealed significant effects for Algorithms [GG(1,32)=9.177, p ≤ 0.0001], Calibration Type [F(1,32)=4.262, p=0.048] and Algorithms*Calibration Type [GG(1,32)=3.894, p=0.008].

Fig. 4 .

 4 Fig. 4. Balanced classification accuracies on the emotion-arousal data set.

Preliminary results on mental workload data only, and with a few existing algorithms only, was published as an extended conference abstract in[START_REF] Appriou | Towards robust neuroadaptive HCI: exploring modern machine learning methods to estimate mental workload from EEG signals[END_REF]

Acknowledgements: This work was supported by the European Research Council (grant ERC-2016-STG-714567) and the Japanese Society for the Promotion of Science.


Sonia Boscolo

Christophe Finot
email: christophe.finot@u-bourgogne.fr

Artificial neural networks for nonlinear pulse shaping in optical fibers

Keywords: pulse shaping, neural networks, nonlinear propagation

Artificial neural networks for nonlinear pulse shaping in optical fibers

Introduction

Pulse shaping based on nonlinear phenomena in optical fibers has become a remarkable tool to tailor the spectral and temporal content of light signals [START_REF] Boscolo | Shaping Light in Nonlinear Optical Fibers[END_REF][START_REF] Hammani | Nonlinear spectral shaping and optical rogue events in fiber-based systems[END_REF], leading to the generation of a large variety of optical waveforms such as ultra-short compressed pulses [START_REF] Mollenauer | Experimental observation of picosecond pulse narrowing and solitons in optical fibers[END_REF], parabolic- [START_REF] Finot | Optical parabolic pulse generation and applications[END_REF], triangular- [START_REF] Boscolo | Passive nonlinear pulse shaping in normally dispersive fiber systems[END_REF] and rectangular- [START_REF] Kutuzyan | Dispersive regime of spectral compression[END_REF] profiled pulses. Very different features of the nonlinear pulse evolution can be observed depending on the sign of the fiber group-velocity dispersion (GVD), which result in specific changes of the pulse temporal shape, spectrum and phase profile.

Specifically, because the nonlinear dynamics of pulses propagating in fibers with normal GVD are generally sensitive to the initial pulse conditions, it is possible to nonlinearly shape conventional laser pulses into various specialized waveforms through control of the initial pulse temporal intensity and/or phase profile [START_REF] Boscolo | Shaping Light in Nonlinear Optical Fibers[END_REF]. Conversely, the temporal and spectral features of pulses propagating in fibers with anomalous GVD are typically governed by soliton dynamics [START_REF]Nonlinear Fiber Optics, Fourth Edition[END_REF]. Yet, due to the typically wide range of degrees of freedom involved, predicting the behavior of nonlinear pulse shaping can be computationally demanding, especially when dealing with inverse problems.

Owing to its power of extracting essential information from large amounts of data, machine learning is bringing a revolutionary reform to research in the physical sciences [START_REF] Carleo | Machine learning and the physical sciences[END_REF]. In the field of photonics, a number of studies have been recently reported in laser design and optimization [START_REF] Boscolo | Performance analysis of dual-pump nonlinear amplifying loop mirror mode-locked all-fibre laser[END_REF][START_REF] Woodward | Towards 'smart lasers': Self-optimisation of an ultrafast pulse source using a genetic algorithm[END_REF][START_REF] Baumeister | Deep learning and model predictive control for selftuning mode-locked lasers[END_REF], complex nonlinear dynamics [START_REF] Närhi | Machine learning analysis of extreme events in optical fibre modulation instability[END_REF], design of photonic crystal fibers and optical components [START_REF] Chugh | Machine learning approach for computing optical properties of a photonic crystal fiber[END_REF][START_REF] Pilozzi | Machine learning inverse problem for topological photonics[END_REF], pulse characterization [START_REF] Zahavy | Deep learning reconstruction of ultrashort pulses[END_REF], and optical communications [START_REF] Giacoumidis | Harnessing machine learning for fiber-induced nonlinearity mitigation in long-haul coherent optical OFDM[END_REF][START_REF] Kopae | Aritificial neural network-based equalizer in the nonlinear Fourier domain for fibre-optic communication applications[END_REF]. In [START_REF] Finot | Nonlinear sculpturing of optical pulses with normally dispersive fiber-based devices[END_REF], we have shown that the combination of a graphical approach with the machine-learning method of neural networks (NNs) can provide a rapid and precise identification of the parameters of nonlinear pulse shaping systems based on pulse propagation in a normally dispersive fiber that are required to generate pulses with preset temporal features. In this paper, we use a supervised learning model based on a NN to predict the temporal and spectral intensity profiles of the pulses that form upon nonlinear propagation in fibers with both normal and anomalous dispersion. We also assess the ability of the NN to solve the inverse problem of determining the nonlinear propagation properties from the pulses observed at the fiber output and to classify the output pulses according to the initial pulse shape. Furthermore, we show how our model can handle the nonlinear shaping of initially chirped pulses.

II. Problem under study and numerical tools

The general scheme for nonlinear shaping that we consider in this paper comprises a pre-chirping stage followed by a nonlinear propagation stage. Within this scheme, an initial pulse 0(t) with a peak power P0 is first propagated through a dispersive medium, such as a pair of diffraction gratings, a prism pair [START_REF] Kalashyan | Generation of transformlimited rectangular pulses in a spectral compressor[END_REF], a segment of hollow core or standard fiber with very low nonlinearity [START_REF] Wang | Generation of triangular-shaped optical pulses in normally dispersive fibre[END_REF][START_REF] Verscheure | Pulse doubling and wavelength conversion through triangular nonlinear pulse reshaping[END_REF]. This linear propagation imprints a parabolic spectral phase onto the pulse:

C0 ω 2 / 2,
where the chirp coefficient C0 equals the cumulative GVD of the medium ( being the angular frequency). The so obtained chirped pulse is then propagated through a fiber that reshapes both its temporal and spectral intensity profiles. Our main interest here is in these physical quantities of the pulse that can be directly recorded in experiments, rather than in the complex envelope of the electric field. Different initial pulse shapes are studied: a Gaussian pulse   Here, T0 is a characteristic temporal value of the pulse and θ(x) is the Heaviside function.

To generate data that characterizes the nonlinear shaping process, we use the nonlinear Schrödinger equation (NLSE) [START_REF]Nonlinear Fiber Optics, Fourth Edition[END_REF], which can describe well a variety of nonlinear phenomena associated with pulse propagation in fibers in spite of the fact that it only includes two physical effects, namely, linear GVD and nonlinear self-phase modulation (SPM). The NLSE for the complex electric field envelope, ψ(z,t), is written as

2 2 2 2 | | 0 2 i z t            ( 1
)
where z is the propagation coordinate, t is the retarded time, and β2 and γ are the respective GVD and Kerr nonlinearity coefficients of the fiber. Note that the effects of linear loss can be neglected given the very low loss of silica fibers in the telecommunication wavelength window.

Here we also neglect higher-order linear and nonlinear effects as the leading-order behavior is well approximated by Eq. (1

u u i u u N          (2)
where LD = T0 2 /|2| and LNL = 1/(P0) are the respective dispersion length and nonlinear length associated with the initial pulse, and the parameter N ('soliton-order' number) is introduced as N 2

= LD/LNL. This way, the nonlinear shaping problem, which depends on the six physical parameters (T0, P0, C0, 2, L) where L is the fiber length, is mapped onto a problem in the three-dimensional space of (, N,

2 00 / C C T 
). This dimensionality reduction relaxes the complexity of the problem, and for a specific selected set of normalized parameters there are many groups of physical parameters suitable the defining equations of ξ, N and C. Equation ( 2) is solved with a standard split-step Fourier propagation algorithm [START_REF]Nonlinear Fiber Optics, Fourth Edition[END_REF], using a uniform grid of 2 13 points on a time window of length 80 T0. In order to generate the temporal and spectral properties for the NN, we perform anamorphic sampling of the output intensity profile provided by Eq. ( 2). This enables us to well represent both the details of the short pulses that can be encountered in the anomalous dispersion regime of the fiber and the longer pulses that are observed in the normal dispersion regime. Given the symmetry of the problem (symmetry of the NLSE and the initial temporal condition), we can restrict our sampling to the positive times and frequencies only. This way we can maintain the number of useful points moderate. Hence, we consider A = 65 points on the interval 0 to 25 T0 to represent the temporal intensity profiles of the pulses, and B = 35 points on the interval 0 to 3.4/T0 for the spectral intensities. The sampled data and the initial conditions are used to train a NN and validate its predictions. We employ a feed-forward NN relying on the Bayesian regularization back propagation algorithm and including three hidden layers with fourteen neurons each, as shown in Fig. 1. This NN is programmed in Matlab using the neural network toolbox.

III. Nonlinear shaping of initially transform-limited pulses

In this section, we study the problem of the nonlinear shaping of initially-transform limited pulses (i.e., with C0 = 0), which can be handled in the parameter space of normalized propagation length and soliton-order number, (, N).

A/ Prediction of the output pulse properties

The training step here involves feeding the NN with an ensemble of 4500 data (both temporal and spectral) generated from numerical simulations of the NLSE and associated with the propagation of a Gaussian pulse with  ranging from 0.025 to 2.5. 1500 randomly chosen data points are used to cover the parameter space in the normal dispersion regime of the fiber with N ranging from 0.025 to 4, and 3000 randomly chosen points for the anomalous dispersion regime with N ranging from 0.025 to 3. After training, the NN is tested on a distinct ensemble of 25000 data not used in the training step. Figure 2 shows the temporal and spectral profiles of the pulse obtained from the network for  = 2 and N = 4 when the fiber has normal GVD, and for  = 1.8 and N = 2.5 at anomalous dispersion. The predictions from the NN algorithm show excellent agreement with the results of the NLSE propagation model over most part of the pulse shape and spectrum. The network is able to reproduce the large temporal and spectral broadening experienced by the pulse upon propagation in the normal dispersion regime. With anamorphic sampling of its output, the network is also able to resolve the details of the temporally compressed pulse and the concomitant more complex structure of the spectrum that are observed after propagation at anomalous dispersion. Nevertheless, some discrepancies with the expected results are visible, and are more pronounced in the anomalous regime in which the propagation dynamics are more complex.

    2 2 , , , 11 2 2 ,
, , 11

AM T NN i NLS i NLS i ii BM F NN i NLS i NLS i ii M I I I M S S S       (3)
In Eq. ( 3), the expected profiles are interpolated to the same time or frequency points used for sampling the network output. The results obtained for 25000 randomly chosen combinations of input parameters  and N in the normal and anomalous dispersion regimes of the fiber are summarized in Fig. 3(a), and confirm that the nonlinear pulse shaping occurring in a fiber with normal GVD is easier to predict than that occurring in the anomalous dispersion regime. This can be physically explained by the more complex nonlinear pulse dynamics that take place in the presence of anomalous dispersion, involving compression and/or splitting stages over short propagation distances. The distributions of values of the temporal and spectral misfit parameters

B/ Retrieval of the propagation characteristics

The inverse problem at hand is much more complex: from a pulse shape and spectrum that are generated after propagation in the fiber, the properly trained network should be able to retrieve the parameters , N as well as the regime of dispersion of the fiber. One could tackle this inverse problem by using reverse propagation, that is, solving the NLSE to the backward direction [START_REF] Ozeki | Clean similariton generation from an initial pulse optimized by the backward propagation method[END_REF][START_REF] Pawlowska | Parametrically shaped femtosecond pulses in the nonlinear regime obtained by reverse propagation in an optical fiber[END_REF]. However, this method would require knowledge of the peak power (or energy) of the pulse at the fiber output, whereas we only consider here the pulse shape characteristics.

For the inverse NN training phase, we use the same data set as that used for the direct problem.

The results of testing the trained NN on 25000 randomly chosen new simulated output pulses are shown in Fig. 4. We can see from Fig. 4(a) that the estimation error on  (and N (N) (defined as the difference between the retrieved parameter value and the target value extracted from the NLSE simulation data) is very close to zero for all test realizations except those corresponding to low input powers or short propagation lengths for which the changes in the temporal and spectral shapes of the propagating pulse are negligibly small, thus leading to similar shapes for different  values. If we limit the data points to the ranges  ≥0.25 and N≥0.25, then, interestingly, each temporal and spectral shape can be unambiguously associated with a single parameter set (, N). The network is able to work out the sign of the fiber dispersion perfectly. The distribution densities of the estimation errors shown in Fig. 4(b) confirm the remarkable accuracy of the results obtained with the NN algorithm: the root-mean square errors on  and N are below 2 10 -3 and 2.7 10 -3 , respectively. We also study the effect of the initial pulse shape on the network's ability to solve the inverse problem. Maps of values of the estimation error on the normalized propagation length are provided in Fig. 5 for 25000 test realizations for each of the initial shapes: hyperbolic secant, parabolic and super-Gaussian pulse. A difficulty arises with the hyperbolic secant pulse propagating in the anomalous regime of the fiber: the estimation error is non-negligible for values of the soliton-order number around 1. This originates in the evolution of the fundamental soliton (N=1), which propagates without changing its shape, thereby making the propagation length difficult to evaluate. Furthermore, higher-order solitons (characterized by integers N>1) follow a periodic evolution pattern, which entails that the same pulse shape may be associated with different propagation lengths [START_REF] Stolen | Observation of pulse restoration at the soliton period in optical fibers[END_REF]. We expect that the effects of energy dissipation or higher-order linear /nonlinear propagation effects, can alleviate this ambiguity.. Initially parabolic pulses can be better handled by the inverse problem NN. The model still works with initially super-Gaussian pulses, whose dynamics can be rather different especially in the anomalous dispersion regime in which the waveforms undergo splitting followed by compression [START_REF] Finot | Optical pulse doublet resulting from the nonlinear splitting of a super-Gaussian pulse[END_REF]. We can also see from Fig. 5 that, again, the anomalous dispersion regime of the fiber is more complex to deal with.

III. Nonlinear shaping of initially chirped pulses

In this section, we include an additional degree of freedom in our analysis by considering a possible pre-chirping of the pulses that are transmitted through the fiber. Our focus is here on pulse propagation in the normal regime of dispersion. Depending on the chirp of the input pulse, the initial stage of nonlinear dynamics in the fiber, where Kerr-induced SPM dominates over GVD, may be very different. Indeed, input pulses with a negative chirp coefficient C0 will experience spectral compression as a result of SPM [START_REF] Oberthaler | Spectral narrowing of ultrashort laser pulses by self-phase modulation in optical fibers[END_REF][START_REF] Finot | Design rules for nonlinear spectral compression in optical fibers[END_REF], whereas for initially positively chirped (or Fourier transform-limited) pulses, spectral broadening will drive the nonlinear dynamics and eventually lead to optical wave-breaking [START_REF] Anderson | Wave-breaking in nonlinear optical fibers[END_REF]. In [START_REF] Finot | Nonlinear sculpturing of optical pulses with normally dispersive fiber-based devices[END_REF], we used a graphical method to find the combinations of values for the parameters N and C that support the formation of pulses with specified temporal features in the fiber.

A/ Prediction of the output pulse properties

The NN learns the NLSE model from an ensemble of 2 10 4 simulations of the propagation of an initial Gaussian pulse with randomly chosen combinations of values for ,N and C over the ranges 0.025 to 3, 1.5 to 5, and -3.4 to 3.4, respectively. For this configuration, the data is sampled over a temporal interval of 30 T0 and a spectral window of 2.3/T0. We then test the trained NN on 7 10 4 simulations from a distinct ensemble of random initial conditions. Figure 7 shows typical examples of the nonlinear pulse shaping occurring in the fiber: formation of parabolic-like pulses with rectangular-like spectra when the chirp of the input pulses is positive

B/ Inverse problem

The results of testing the inverse problem NN on 7 10 4 new simulated output pulses are summarized in Fig. 8. We note that the fiber nonlinearity is essential to distinguishing the combinations of ,N and C values that are associated with each test pulse. Indeed, pulse prechirping entails a drop of the pulse peak power, which in turn results in a frequent occurrence of estimation errors on ,N and C for the test realizations that relate to low values of N. Therefore, we restrict the statistical error analysis to the test data associated with N >1.5. We can see from

IV. Conclusion

We have successfully used a supervised machine-learning model based on a NN to solve both the direct and inverse problems relating to the shaping of optical pulses that occur upon nonlinear propagation in optical fibers. Remarkably, within the range of system parameters considered, any temporal and spectral shape generated at the fiber output from the propagation of an initially Fourier-transform limited Gaussian pulse can be unambiguously associated with a single set of normalized propagation length and soliton-order number values, (, N). Our results

show that a properly trained network can greatly help the design and characterization of fiberbased shaping systems by providing immediate and sufficiently accurate solutions. We have limited our discussion here to passive non-dissipative propagation. An interesting extension of the method described will be to the analysis of the attraction of parabolic pulses towards a selfsimilar state in normally dispersive nonlinear fibers with linear gain [START_REF] Finot | Optical parabolic pulse generation and applications[END_REF]. Furthermore, although demonstrated here in a fiber optics context, the principle of using NN architectures to solve wave equation-based inverse problems is expected to apply to many physical systems [START_REF] Khoo | SwitchNet: A neural network model for forward and inverse scattering problems[END_REF][START_REF] Sanghvi | Embedding deep learning in inverse scattering problems[END_REF].

Figure 1 :

 1 Figure 1: Neural network model used to: (a) predict the output pulse intensity profiles from the fiber, and (b) to retrieve the nonlinear propagation properties.

Figure 2 :

 2 Figure 2: Temporal (panels 1) and spectral (panels 2) intensity profiles of an initial Gaussian pulse after propagation in (a) a normally dispersive fiber with  = 2 and N = 4, and (b) an anomalously dispersive fiber with  = 1.8 and N = 2.5. The predictions from the NN algorithm (blue circles) are compared with the results of NLSE numerical simulations (black curves). Also shown are the input intensity profiles (red dotted curves

(Fig. 3 (

 3 b)) show that, remarkably, more than 90% of the error realizations are well confined to values below 0.04, but some error values are spread out over a wider range. This deviation mostly occurs in the anomalous dispersion regime.

Figure 3 :

 3 Figure 3: (a) Maps of misfit parameter values between the NN predictions and the NLSE simulation results for the output pulse shape (subplot a1) and optical spectrum (subplot a2) in the two-dimensional space (N,), relating to the propagation of a Gaussian pulse in the normal and anomalous dispersion regimes of the fiber with randomly chosen combinations of  and N. (b) Distribution densities of the temporal (subplot b1) and spectral (subplot b2) misfit parameter values.

Figure 4 :

 4 Figure 4: (a) Maps of estimation error values on the soliton number N (subplot a1) and normalized propagation length  (subplot a2) in the two-dimensional space (N, ) for both normal and anomalous dispersion, when the NN is interrogated with randomly chosen new simulated output pulses from the fiber. The red dashed lines delimit the data domain that is used for the statistical error analysis. (b) Distribution densities of the estimation errors on N (subplot b1) and  (subplot b2).

Figure 5 :

 5 Figure 5: Maps of estimation error values on the normalized propagation length  in the twodimensional space (N, ) for both normal and anomalous dispersion, when the NN is interrogated with randomly chosen new simulated output pulses from the fiber corresponding to input (a) hyperbolic secant, (b) parabolic and (c) super-Gaussian pulses. 

Figure 6 :

 6 Figure 6: Points where errors in the detection of the initial pulse shape occur (subplot a1) and map of estimation error values on the normalized propagation length  (subplot a2) in the twodimensional space (N, ) for both normal and anomalous dispersion, when the NN is interrogated with randomly chosen new simulated output pulses from the fiber corresponding to an unlabeled mix of input Gaussian, hyperbolic secant, parabolic and super-Gaussian pulses. In subplot a1, the color of the circle indicates the class label that is not predicted correctly. The red dashed lines delimit the data domain that is used for the statistical error analysis. (b) Distribution densities of the estimation errors on N (subplot b1) and  (subplot b2).

[

  subplot (a)] [29], and narrowing of the pulse spectrum for negative input chirp [subplot (b)] [27]. It is clear that the NN algorithm performs impressively in reproducing the output pulse shapes.

Figure 7 :

 7 Figure 7: Temporal (panels 1) and spectral (panels 2) intensity profiles of an initial Gaussian pulse after propagation in a normally dispersive fiber with (a)  = 2.3, N = 3, and C= -2.42, and (b)  = 2, N = 3 and C= 2.42. The predictions from the NN algorithm (blue circles) are compared with the results of NLSE numerical simulations (black curves). Also shown are the input intensity profiles (red dotted curves).

Fig. 8

 8 Fig.8that the values of the input parameters obtained from the network algorithm are in agreement with the known values from the simulation data, but the estimation errors are higher than those made in the case of initially transform-limited pulses: the rms deviations are 0.36, 0.09 and 0.55 for  , N and C respectively.

Figure 8 :

 8 Figure 8: (a) Map of estimation error values on the normalized propagation length  in the threedimensional space (N, , C), when the NN is interrogated with randomly chosen simulated output pulses from the fiber. (b) Regressions between the predicted values of N,  and C from the NN algorithm and the exact target values from the simulated data. (c) Distribution densities of the estimation errors on N,  and C.

Acknowledgements

We acknowledge the support of the Institut Universitaire de France (IUF).