Etienne Vareille
email: etienne.vareille@ensea.fr

Adel Abbas
email: adel.abbas@ensea.fr

Michele Linardi
email: michele.linardi@cyu.fr

Vassilis Christophides
email: vassilis.christophides@ensea.fr

Evaluating Explanation Methods of Multivariate Time Series Classification through Causal Lenses

Keywords: Multivariate Time Series, Explainability, Classification, Causal Explanations, Consistency, Benchmarking

Explainable machine learning techniques (XAI) aim to provide a solid descriptive approach to Deep Neural Networks (NN). In Multi-Variate Time Series (MTS) analysis, the most recurrent techniques use relevance attribution, where importance scores are assigned to each TS variable over time according to their importance in classification or forecasting. Despite their popularity, post-hoc explanation methods do not account for causal relationships between the model outcome and its predictors. In our work, we conduct a thorough empirical evaluation of model-agnostic and model-specific relevance attribution methods proposed for TCNN, LSTM, and Transformers classification models of MTS. The contribution of our empirical study is threefold: (i) evaluate the capability of existing post-hoc methods to provide consistent explanations for high-dimensional MTS (ii) quantify how post-hoc explanations are related to sufficient explanations (i.e., the direct causes of the target TS variable) underlying the datasets, and (iii) rank the performance of surrogate models built over post-hoc and causal explanations w.r.t. the full MTS models. To the best of our knowledge, this is the first work that evaluates the reliability and effectiveness of existing xAI methods from a temporal causal model perspective.

I. INTRODUCTION

Multivariate Time Series (MTS) are omnipresent in many science and engineering domains, including health-care, sustainable energy, geoscience, and high-performance computing. An MTS is composed of more than one time-dependent variables that may depend not only on its past values but also on other variables. Neural Network (NN) architectures (e.g., TCNN [START_REF] Boniol | dcam: Dimensionwise class activation map for explaining multivariate data series classification[END_REF], LSTM [START_REF] Hochreiter | Long short-term memory[END_REF] Transformer [START_REF] Vaswani | Attention is all you need[END_REF]) are today state-of-theart solutions in order to implement MTS Classification [START_REF] Nemer | A review of research on industrial time series classification for machinery based on deep learning[END_REF]- [START_REF] Ismail | Benchmarking deep learning interpretability in time series predictions[END_REF].

Deep NN models are essentially uninterpretable black boxes, where one feeds an input and obtains an output without understanding the motivations behind that decision. We have recently witnessed a consistent effort to enhance NN models with explanation capabilities at various levels allowing users to track the time-dependent variables or the training timestamped samples that drive a NN model toward a certain decision. Specifically, eXplainable Artificial Intelligence (xAI) techniques such as dCAM [START_REF] Boniol | dcam: Dimensionwise class activation map for explaining multivariate data series classification[END_REF], XCM [START_REF] Fauvel | Xcm: An explainable convolutional neural network for multivariate time series classification[END_REF], Dynamask [START_REF] Crabbé | Explaining time series predictions with dynamic masks[END_REF], TimeSHAP [START_REF] Bento | Timeshap: Explaining recurrent models through sequence perturbations[END_REF] aim to reveal which subsets of MTS variables are mostly involved in deciding a particular class label.

Post-hoc explanation methods like feature attribution [START_REF] Ribeiro | why should i trust you?": Explaining the predictions of any classifier[END_REF], [START_REF] Lundberg | A unified approach to interpreting model predictions[END_REF] or saliency maps [START_REF] Simonyan | Deep inside convolutional networks: Visualising image classification models and saliency maps[END_REF] can help indicate which features or pixels of input data are most 'relevant' to the output decision of a NN model being 'explained'. However, they are designed with the assumption that a single sample is selfexplanatory and analysts can easily distinguish between two samples, which is not easy in complex MTS with thousands of time-dependent variables. Additionally, the intrinsic measures underlying xAI methods are purely associational and may expose potentially spurious and misleading correlations in input data used to train a NN model. Clearly, correlation-based explanations are not informative enough in order to produce surrogate models of reduced complexity (i.e., by projecting input data on the relevant features) or to debug models (e.g., by identifying data instances and/or model components that cause incorrect predictions).

These kinds of tasks call for more expressive forms of explanations, which can reveal causal-effect relationships between TS variables (input or target) along with adequate lag information. Causal explanations are usually distinguished between sufficient explanations and counterfactual explanations [START_REF] Beckers | Causal explanations and XAI[END_REF]. The former permits users to understand the conditions in which a particular action (e.g., selecting a subset of TS variables as predictors) will produce a desired model outcome. The latter instead identifies actions (e.g., changing the values of TS variables) that can alter an observed input, ensuring a change in a previously observed output.

In our work, we are experimentally evaluating the consistency and effectiveness of model-specific (i.e, dCAM [START_REF] Boniol | dcam: Dimensionwise class activation map for explaining multivariate data series classification[END_REF], Dynamask [START_REF] Crabbé | Explaining time series predictions with dynamic masks[END_REF], XCM [START_REF] Fauvel | Xcm: An explainable convolutional neural network for multivariate time series classification[END_REF]), and model-agnostic (i.e, Feature Ablation [START_REF] Freeborough | Investigating explainability methods in recurrent neural network architectures for financial time series data[END_REF], [START_REF] Suresh | Clinical intervention prediction and understanding using deep networks[END_REF], Feature Permutation [START_REF] Zeiler | Visualizing and understanding convolutional networks[END_REF], Feature Occclusion [START_REF] Zeiler | Visualizing and understanding convolutional networks[END_REF], Integrated Gradients [START_REF] Sundararajan | Axiomatic attribution for deep networks[END_REF], TimeSHAP [START_REF] Bento | Timeshap: Explaining recurrent models through sequence perturbations[END_REF], Gradient Shap [START_REF] Lundberg | A unified approach to interpreting model predictions[END_REF]) post-hoc explanations of MTS models built using different NN architectures (i.e, TCNN [START_REF] Boniol | dcam: Dimensionwise class activation map for explaining multivariate data series classification[END_REF], LSTM [START_REF] Hochreiter | Long short-term memory[END_REF] and Transformer [START_REF] Vaswani | Attention is all you need[END_REF]). More precisely, the contribution of our empirical study is three-fold: (i) evaluate the ability of existing post-hoc methods to provide consistent explanations for highdimensional MTS (ii) quantify how post-hoc explanations are related to sufficient explanations (i.e., the direct causes of the target TS variable) underlying the datasets, and (iii) rank the performance of surrogate models built over post-hoc and causal explanations w.r.t. the full MTS models.

None of the previous empirical studies ( [START_REF] Ismail | Benchmarking deep learning interpretability in time series predictions[END_REF], [START_REF] Jacob | Exathlon: A benchmark for explainable anomaly detection over time series[END_REF]) compare both model-agnostic and model-specific post-hoc explanation methods for a variety of MTS classification models (LSTM, CNN, and Transformer) using as a cause-effect ground truth underlying the benchmark datasets as validated by the domain experts. More precisely, [START_REF] Ismail | Benchmarking deep learning interpretability in time series predictions[END_REF] focuses only on MTS agnostic xAI methods (not tailored for MTS classification models) evaluated using AUPR (Area under the precision-recall curve) given as ground truth feature relevance constructed independently of the underlying data generation processes (i.e., by adding or subtracting a fixed constant from the data). The Exathlon benchmark [START_REF] Jacob | Exathlon: A benchmark for explainable anomaly detection over time series[END_REF] compares two anomaly explanation methods (Macrobase, Xstream) with LIME explanations for three NNbased anomaly detectors (LSTM, Auto Encoder, and BiGAN). Exathlon provides ground truth only for the range-based anomalies occurring in Spark execution traces and explanation conciseness and consistency are evaluated independently of the true anomaly causes. As a matter of fact, these metrics report only the number of features used in the explanations and how often anomalies of the same type receive the same explanation. Although the need for experimentally evaluating sufficient explanations of target variables over time has been raised in previous works [START_REF] Crabbé | Explaining time series predictions with dynamic masks[END_REF], [START_REF] Beckers | Causal explanations and XAI[END_REF], this is the first empirical study demonstrating that surrogate models built over causal relationships outperform not only models built over relevant features detected by xAI methods but also advanced NN models trained over the whole data feature space.

The rest of the paper is structured as follows: Section II introduces the core notions of post-hoc and sufficient explanations. Section III presents the main aspects of the analyzed xAI methods. Section IV describes our empirical evaluation setting and metrics. Section V details the results of our experimental evaluation and highlights the main insights. Finally, Section VI summarizes our findings and presents plans for future works.

II. NOTATION AND PRELIMINARIES

We denote a multivariate time series (MTS) as X ∈ R N ×T , where T is the number of time steps and N is the number of dimensions. Each MTS dimension is a univariate TS variable denoted as X i , where 1 ≤ i ≤ N . The length of a TS variable X i is denoted as |X i |. Since all univariate TS variables of an MTS have the same length, we use the notation |X| to also denote the number of time steps T .

We define a sub-sequence of X i of length ℓ, starting at time-step t, with X i t:t+ℓ-1 ∈ R N ×ℓ , and the value at time t of the i th TS variable as X i t ∈ R. Given X, we define a multivariate sub-sequence (set of sub-sequences of all TS variables) of length ℓ and starting at time-step t as X

(1..N ) t:t+ℓ-1 = {X 1 t:t+ℓ-1 , .., X N t:t+ℓ-1 }. Note that a multivariate sub-sequence is a MTS itself.
To consider a scale-free representation of multivariate subsequences we apply Z-normalization [START_REF]Searching and mining trillions of time series subsequences under dynamic time warping[END_REF], [START_REF] Keogh | On the need for time series data mining benchmarks: A survey and empirical demonstration[END_REF]. Hence, given a MTS X, each TS variable is always in the form

X i = { X i 1 -µ σ , ..., X i |X i | -µ σ
}, where µ and σ are the mean and standard of X i respectively.

We define MTS classification as a function F c that given a multivariate sequence outputs the probability of the sequence to belong to a given class. Formally, we have:

F c (X) = [S 1 (X), ..., S C (X)] ∈ R C ,
where C is the total number of classes and S i (X) represents the probability of the sequence X to belong to the i th class. Hence, C i=1 S i (X) = 1. We consider the explainability of an MTS classification model in terms of relevance attribution, which assigns a score to each TS variable at each time step. This score typically reports the feature's importance with respect to a class assignment. Fig. 1 illustrates an example of a Relevance Attribution Matrix (Fig. 1(b)) computed on a MTS X (Fig. 1(a)). Definition 1 (Relevance Attribution Matrix [START_REF] Ismail | Benchmarking deep learning interpretability in time series predictions[END_REF], [START_REF] Zeiler | Visualizing and understanding convolutional networks[END_REF], [START_REF] Ancona | Towards better understanding of gradient-based attribution methods for deep neural networks[END_REF]): Given an MTS X, for which a model F c provides a probability to belong to the class a, the Relevance Attribution Matrix is defined as R(X,

F c (X)[a]) ∈ R N ×ℓ . Each R i,j ∈ R(X, F c (X)[a]
) reports the importance score of a TS value X i j for the model F c in determining the prediction for the class a. We denote by R k (X, F c (X)[a]) the set containing the k ∈ N highest values in R(X, F c (X)[a]).

In our work, we also consider causal explanations [START_REF] Pearl | The seven tools of causal inference, with reflections on machine learning[END_REF] in terms of cause-effect relationships among MTS variables. Causal relationships are typically stated through a temporal graphical model that we introduce hereafter.

Definition 2 (Window Causal Graph [START_REF] Assaad | Survey and evaluation of causal discovery methods for time series[END_REF]): Given an MTS X ∈ R N ×T a window causal graph is defined by G t = (V, E) at time t and maximal lag ℓ, where V = {X i t-s |1 ≤ i ≤ N ∧ 0 ≤ s ≤ ℓ} is the set of vertices each one representing a feature value. A directed edge (X a t-τ , X b t ) ∈ E exists if X a causes the effect X b at time t with lag τ ≤ ℓ. Note that if a = b =⇒ τ > 0 (self-causation).

A window causal graph is consistent throughout time when all the causal relationships remain constant (Consistency Throughout Time [START_REF] Assaad | Survey and evaluation of causal discovery methods for time series[END_REF] i.e., Stationary), and always precede their effects (Temporal Priority [START_REF] Assaad | Survey and evaluation of causal discovery methods for time series[END_REF]). In this case, we have that G t (V, E) = G t ′ (V, E)∀t, t ′ .1 ≤ t, t ′ ≤ |X| and thus we denote a single window causal graph by G(V, E). In that case, we simply refer to the MTS causal graph. In fact, finding a causal explanation of a model outcome is equivalent to giving the actual causes of that outcome [START_REF] Woodward | Making Things Happen: A Theory of Causal Explanation[END_REF], [START_REF] Halpern | Causes and explanations: A structural-model approach. part i: Causes[END_REF]. Hence, causal explanations contain the direct causes of a target TS variable as defined in a causal graph.

Definition 3 (Causal Explanation [START_REF] Runge | Causal network reconstruction from time series: From theoretical assumptions to practical estimation[END_REF]): Given a MTS X faithful to a causal graph G(V, E) with maximal lag ℓ, the set of direct causes of an outcome X i t is defined as follows:

C i t (X) = {(a, t ′ )|t -ℓ ≤ t ′ < t ∧ (X a t ′ , X i t ) ∈ E}.
Under the assumptions of (I) faithfulness (only the conditional independence relations true in the data are entailed by the Causal Markov condition applied to the graph) (II) causal sufficiency (no unobserved latent factor influences two observed vertices), (III) correct independence tests (a.k.a. Causal Markov condition [START_REF] Assaad | Survey and evaluation of causal discovery methods for time series[END_REF]: nodes in a causal graph are conditional independent of its nondescendant, given their parents), (IV) temporal priority (causes occur before their effects), (V) minimality (the causal graph does not contain dependencies not present in the data), and (VI) absence of selection bias, the set of direct causes contains all the observed TS variables that are parents of a variable X i t in G(V, E), and it is sufficient to predict a target variable X i t [START_REF] Runge | Causal network reconstruction from time series: From theoretical assumptions to practical estimation[END_REF]. Formally, variable X i t is independent of X a t ′ conditionally on C i t (X) where t ′ < t.

III. EXPLANATION METHODS

In this work, we focus on model-specific and model-agnostic post-hoc explanations of TS classification models produced by Perturbation-based and Gradient-based methods.

A. Perturbation based methods

This family of xAI methods assesses feature relevance by perturbing their values (masking, removing, altering, permuting) and measuring the impact on the model outcome. Dimension-wise Class Activation Map (dCAM) [START_REF] Boniol | dcam: Dimensionwise class activation map for explaining multivariate data series classification[END_REF] is a model-specific method that adapts a Class Activation Mapping (CAM) [START_REF] Selvaraju | Grad-cam: Visual explanations from deep networks via gradient-based localization[END_REF] to MTS classification. CAM represents one of the earliest explainability techniques natively applied to image classification in Convolutional Neural Networks (CNN). CAM aims to compute the discriminating data region that induces the model to assign a particular class to a given instance.

We denote by A ∈ R Z the activation map generated by an MTS X ∈ R N ×T in the last convolutional filter of the network (before Global Average Pooling layer). Note that A i,j denotes the activation at time step i (1 ≤ i ≤ T ) of the j th kernel. Given w j c , the weight connecting the activation to the neuron of a class c to the kernel j, the CAM score for class c at time step i is computed as follows:

CAM c i = j w j c A i,j (1)
Each CAM score captures the activation conveying to a single class neuron at timestamp i. Note that such scores would not permit filtering out important features for predicting each class. As a matter of fact, state-of-the-art xAI methods for images (e.g., CAM) are not specifically tailored to MTS, where discovering the relevance of features over different time intervals becomes a crucial requirement. In [START_REF] Ismail | Benchmarking deep learning interpretability in time series predictions[END_REF] is empirically shown that such limitation is principally due to the conflation of time and feature representations in the NN layers.

To overcome this limitation, dCAM proposes a new CNN architecture in which the first convolutional filter applies to a TS cube containing all the features (a.k.a TS variables) permutations (of the input instance) along time. To compute relevance over time, dCAM does not require investigating every possible permutation but summarizes feature importance at each timestamp, by applying a heuristic over the permutations cube to obtain the final relevance attribution matrix. Dynamask (DM) [START_REF] Crabbé | Explaining time series predictions with dynamic masks[END_REF] relies on a perturbation operator to create a modified version of the model input in order to generate relevance scores. The operator is defined by considering neighboring values of each feature at different timestamps. By comparing the perturbed prediction to the original prediction, the produced errors can be backpropagated to adjust the matrix scores in the final output. Unlike DM and dCAM tailored to explain DL-based MTS classifiers, the following perturbationbased methods are model agnostic and aim to assess the impact of individual features in terms of model outcome difference, when the features are masked or permuted. SHapley Additive exPlanations (SHAP) [START_REF] Lundberg | A unified approach to interpreting model predictions[END_REF] adopt a fundamental concept from cooperative game theory, which consists of assigning importance value (a.k.a SHAP value) according to the effect of a given feature on the model prediction. SHAP values are regression coefficients whose computation requires retraining the model to explain each feature, and thus, it rapidly becomes a computationally expensive solution. Hence, many approximations have been proposed in this respect, for non-sequential and sequential models (a.k.a time series). Gra-dientSHAP [START_REF] Lundberg | A unified approach to interpreting model predictions[END_REF] (GS) is a hybrid approach (gradient-based + perturbation) that approximates SHAP values by adding white noise to each input sample multiple times. It works by selecting a random baseline and a random point along the path between the baseline and the input, computing the gradient of outputs w.r.t. the random points. The final SHAP values represent the expected values of gradients multiplied by the difference between the input and the baselines. TimeSHAP [START_REF] Bento | Timeshap: Explaining recurrent models through sequence perturbations[END_REF] (TS) works by creating a parallel linear model to the one to explain, training it with slightly perturbed data (over time). Since the space of different feature sets to perturb grows exponentially, TimeSHAP proposes a sampling strategy, that merges sets (a.k.a. coalitions) of semantically equivalent MTS dimensions and time-steps. The perturbation must change the output of the model, but at the same time, it must be large enough to identify the most important features for predictions. Feature Ablation (FA) [START_REF] Freeborough | Investigating explainability methods in recurrent neural network architectures for financial time series data[END_REF], [START_REF] Suresh | Clinical intervention prediction and understanding using deep networks[END_REF] performs the classification task by replacing each timestamp of an observed TS variable with a given baseline value, and then computing the difference from the outcome generated by the original instance. Such a difference quantifies the importance of the masked feature. Feature Occlusion(FO) [START_REF] Zeiler | Visualizing and understanding convolutional networks[END_REF] masks all attributes over a time interval (multivariate sub-sequence), obtaining thus an importance score per time window. The final attribution of each value becomes the mean value of all the window scores in which the value appears. When Z-Normalization is applied, it is common to consider a zero value and a zero-valued TS baseline in FA and FO, respectively. Feature Permutation(FP) [START_REF] Zeiler | Visualizing and understanding convolutional networks[END_REF] produces relevance scores by permuting each feature individually and drawing a new feature value randomly in batches.

B. Gradient-based Approaches

This family of xAI methods relies on gradient computation of neural network outputs to quantify feature relevance. Gradient CAM (Grad-CAM) [START_REF] Selvaraju | Grad-cam: Visual explanations from deep networks via gradient-based localization[END_REF] is a generalization of CAM aiming to output relevance scores based on all the network layers. Such scores are calculated by backpropagating the gradient of each class with respect to the activation generated by each observed TS variable at each timestamp. XCM [START_REF] Fauvel | Xcm: An explainable convolutional neural network for multivariate time series classification[END_REF] was the first model-specific method adapting Grad-CAM to obtain feature relevance scores for MTS classification. Formally, the importance in terms of the neuron activation of a feature f for a given class y c is denoted by:

α c f = 1 N × T i j δy c δA f i,j (2)
The coefficient α c f reports the global average of the (backpropagated) gradients of a class score y c (output of the softmax layer) w.r.t. the activation map A f obtained by the application of the first 2-dimensional convolution filter on a TS variable instance X i t ∈ R N ×T . To compute the relevance attribution, Grad-CAM involves a combination of the relevance coefficient α c f and all the feature maps. Hence, the final Grad-CAM score attribution for a class c is given according to the following formulae:

L c 2D = ReLU ( f α c f A f ) (3)
Note that here, a linear rectifier (ReLU) filter keeps only positive attributions, while L c 2D reports a score for each observed TS variable at each timestamp. Integrated Gradients (IG) is a model-agnostic method [START_REF] Sundararajan | Axiomatic attribution for deep networks[END_REF] that considers the integral of gradients of the model output with respect to the inputs along the straight line path (in R N ×T ) from a baseline instance X ′ to the input X.

IV. METHODOLOGY AND EXPERIMENTAL FRAMEWORK

Testbed Datasets Our experimental evaluation considers families of datasets that dispose of casual graph validated by domain experts. Other available time series benchmarks like Exathlon [START_REF] Jacob | Exathlon: A benchmark for explainable anomaly detection over time series[END_REF] require to first run and evaluate causal discovery algorithms [START_REF] Assaad | Survey and evaluation of causal discovery methods for time series[END_REF], [START_REF] Nauta | Causal discovery with attentionbased convolutional neural networks[END_REF]- [START_REF] Huang | Fast and accurate causal inference from time series data[END_REF] that is outside the scope of our work. Each dataset family is a collection of instances (MTS) with similar properties. For the sake of reproducibility, the complete framework containing the classification models and the explanation algorithms is available online [START_REF] Vareille | Time Series Causal Lens[END_REF]. We have integrated into our testbed the source code of the model agnostic methods provided by the library Captum 1 , as well as, the original implementations of model-specific methods (dCAM, XCM, TimeShap and Dynamask) provided by their authors.

Table I summarizes the main characteristics of the datasets of the three families included in our testbed. Each dataset family contains several MTS instances (number reported in the previous table) of diverse length and dimensionality (i.e., TS variables) committing to causal graphs (see Fig. 2) of varying structure (w.r.t. node number and in-degree) and temporal dependencies (i.e., linear vs nonlinear). Instances of the CLIM family [START_REF] Runge | The causality for climate competition[END_REF] encode linear temporal relationships estimated from a Vector Auto-Regressive model fit on real-world data and processed by domain experts of climate modeling. Instances of the fMRI family [START_REF] Smith | Network modelling methods for fmri[END_REF] simulate fMRI-BOLD levels in a brain network using a domain-specific model and includes nonlinear relations as it is generated from a linear ODE with constant coefficients. Finally, instances of 7ts2h [START_REF] Assaad | 7ts2h structure[END_REF] are fully synthetically generated, using a structural model that includes sines, cosines, or absolute values.

A. Datasets

7ts2h and fMRI exhibit different distributions of MTS instances, as the generative stochastic process relies on different weights and coefficients from one instance to another. On the other hand, their causal graphs differ by two edges in the worst case. 7ts2h is the only family where the data-generating process contains hidden confounders (i.e., non-observed TS variables), which are part of true causes in the causal graphs.

The causal graph (with lag values information) is available for each instance in all families except for CLIM. For this latter, we obtain causal information running SLARAC [START_REF] Weichwald | Causal structure learning from time series: Large regression coefficients may predict causal links better in practice than small p-values[END_REF] algorithm, as suggested by Runge et al. in their benchmark [START_REF] Runge | The causality for climate competition[END_REF]. In the obtained causal model, we observe a significant variance in the number of edges across the instances. Moreover, we consider that all the cause-effect pairs have a lag, either one or two, as indicated in the context of the causality4climate competition [START_REF] Runge | The causality for climate competition[END_REF]. Problem Settings In our study, we consider the problem of predicting the value of a target TS variable at time t+1, namely X i t+1 , given the MTS subsequence of the preceding window (of length ℓ) denoted by X 1,..,N (t-ℓ+1):t (a.k.a predictor). Since the evaluated explanation methods apply to (deep learning) classification models, we compute binary labels from continuous target values. In this respect, we rely on the KBinsDiscretizer available in Scikit-learn. Hence, the probability of X i t+1 to belong to one class or another is given by F c (X 1,..,N (t-ℓ+1):t ). The total number of timestamps for which we predict the value of the target TS variable is |X| -ℓ. NN Architectures for Classification To implement MTS classification we rely on state-of-the-art Neural Network (NN) architectures such as Temporal Convolutional Neural Network (TCNN)(dCAM [START_REF] Boniol | dcam: Dimensionwise class activation map for explaining multivariate data series classification[END_REF], XCM [START_REF] Fauvel | Xcm: An explainable convolutional neural network for multivariate time series classification[END_REF]), Long-short Term Memory Network (LSTM) [START_REF] Hochreiter | Long short-term memory[END_REF] and Transformer [START_REF] Vaswani | Attention is all you need[END_REF] Neural Networks.

Unlike the standard deployment of TCNN and LSTM architectures [START_REF] Ismail | Benchmarking deep learning interpretability in time series predictions[END_REF], we have made adjustments to the Transformers NNs architecture commonly adopted in the literature [START_REF] Zerveas | A transformer-based framework for multivariate time series representation learning[END_REF]. In particular, we use only a transformer encoder replacing ReLU with GELU activations (in the Position-Wise Feed-Forward Layer). To obtain the final outcome, a single feature vector (per batch) resulting from max-pooling the temporal features is fed to a two-layer multi-layer perceptron, which uses GELU activation prior to producing the classification results (in the softmax layer). Our implementation relies on [START_REF] Lippe | UvA Deep Learning Tutorials[END_REF] that learns positional encoding instead of the default sine-cosine scheme in order to obtain a better classification performance overall.

Model training We train a model of each DL architecture type on the prediction task, for each target variable and each MTS instance (in total 1155 models). To obtain meaningful insights from explanation methods, for each NN architecture we consider models built with MTS instances that have an AUCROC above 0.7. The final number of retained models per dataset family and DL architecture is reported in Table II, for a total of 537 models. We perform MTS classification over windows (i.e., sub-sequences) of eight observations. Such window length permits to cover maximal lags between causes and effects in the causal graph, along with several previous timestamps whose influence might be correlated to the target TS variable. We use a 70:30 splitting ratio for training/validation and test respectively, adopting Forward Chaining Cross Validation [START_REF] Bergmeir | On the use of cross-validation for time series predictor evaluation[END_REF] (FCCV). During training, we oversampled MTS instances to balance classes. We rely on a stratified K-Fold data split to tune NN hyperparameters using the optuna tool. During this process, we apply K-Fold to a subset of instances (up to three, depending on the dataset size and time complexity limit) and uniquely on training data.

B. Explanation evaluation metrics

xAI methods are usually evaluated along two axes: (i) to what extent the generated explanations meet user expectations (Plausibility); (ii) how accurately they reflect the predictive model to explain (Fidelity). In this paper, we focus on the latter criterion, and specifically, we consider the following metrics: (a) temporal consistency of explanations produced in different windows; (b) the relevance of explaining features to the underlying causal graph; (c) the predictive performance of the surrogate models trained only on the explaining features.

1) Temporal Explanation Consistency: In explainability, the notion of consistency typically measures to what extent two samples of the same class are explained in the same manner. In the case of the datasets we use, where each MTS value originates from a unique generation process, we expect that explanations are locally consistent over time. Hence, we measure the temporal consistency of an xAI method as the degree of randomness exhibited by explanations across different windows. Specifically, we adopt a definition of consistency similar to the concordance metric introduced in the Exathlon benchmark [START_REF] Jacob | Exathlon: A benchmark for explainable anomaly detection over time series[END_REF].

More precisely, given an MTS X ∈ R N ×T and a classification model F c , we compute consistency over the k highest scores of relevance attribution matrices resulting from the application of a classification model F c on each sub-sequence in X of length ℓ. Recall that the number of sub-sequences of length ℓ in X is given by the number of target TS variables we predict in X, namely |X| -ℓ. Hence, we define:

Consistency = - N i=1 ℓ j=1 p k (i, j) log 2 p k (i, j) (4)
p k (i, j) = |X|-ℓ t=1 1[R i,j ∈ R k (X 1,..,N t,(t+ℓ-1) , F c (X 1,..,N t,(t+ℓ-1) )[a])] k(|X| -ℓ)
(5) In equation ( 4), p k (i, j) corresponds to the relevance scores frequency of all sub-sequences at position i, j in the top k scores of each relevance attribution. Each relevance attribution matrix is computed with respect to a class a, where each sub-sequence belongs to. Consistency is thus the entropy of this distribution. In equation ( 5) 1[] is the indicator function.

The sum of all frequencies

N i=1 ℓ j=1 p k (i, j) is equal to 1.
Consistency is bounded between a fixed theoretical maximum log 2 (N ×ℓ) and a theoretical minimum log 2 k. As Consistency measures the entropy of explanations across time, the lower the value the more consistent the set of k-length explanations we obtain across time.

2) Explanation precision w.r.t. a Causal Graph: To measure whether highly scored relevance features correspond to true causes of the target variable, we introduce Precision and Recall for a fixed explanation length k, using as ground truth the causal graph of each dataset.

Given a MTS X ∈ R N ×T and a classification model F c applied on a number |X| -ℓ of sliding windows of X to predict the class of values in a target X i , we define:

P = |X|-ℓ t=1 |{R a,b ∈ R k (X 1,..,N t,(t+ℓ-1) |(a, b) ∈ C i t+ℓ (X)}| k(|X| -ℓ) (6) R = |X|-ℓ t=1 |{(a, b) ∈ C i t+ℓ (X)|R a,b ∈ R k (X 1,..,N t,(t+ℓ-1) }| |C i t+ℓ (X)|(|X| -ℓ) (7)
P measures on average how many relevant features are direct causes of the target attributes in X. Symmetrically, R reports the average number of direct causes that are relevant according to the explanation provided in the top k scores of each relevance attribution matrix.

3) Predictive Performance of Explanations: The predictive performance of explanations is often evaluated using the Insertion Deletion method [START_REF] Petsiuk | Rise: Randomized input sampling for explanation of black-box models[END_REF] that measures the impact of top k salient features on the classification performance (e.g., AUC), either by masking those or by masking all other features. While useful for comparing different xAI methods, this method does not suit binary causal relationships with a target MTS variable. While it is possible to define the causal strength of directed feature relationships [START_REF] Fitelson | 60029 Probabilistic measures of causal strength[END_REF], such information is not available in our ground truth.

In our study, we are interested in measuring the predictive performance of the final classification outcome of the known causes. In this respect, we rely on the remove-andretrain method (ROAR) [START_REF] Hooker | A benchmark for interpretability methods in deep neural networks[END_REF], which consists of re-training a model with reduced data dimensionality and testing its performance. In a nutshell, for each MTS and model we build and compare two different models: I) the Explanation-Only model (EO model), which is is trained with only the k most frequently salient features for p k , and II) the Cause-Only model (CO model) trained with only the features that are direct causes of the target variable.

To compare the performance of these two kinds of surrogate models, we measure the Area Under the Receiver Operating Characteristics Curve (ROC AUC) while performing FCCV. It is worth mentioning that we can not exploit the Insertion Deletion method [START_REF] Petsiuk | Rise: Randomized input sampling for explanation of black-box models[END_REF], measuring the AUC of model accuracy variation, when k changes. Since no order exists on the causes' importance in the causal graph ground truth, k must remain fixed to the number of true causes (k = |C i t (X)|).

V. EXPERIMENTAL RESULTS

In this Section, we report the performance results of the explanation methods presented in Section III according to the metrics we introduced in Section IV. We finally highlight the main conclusions drawn from our experiments.

We evaluate each metric on our 30% test data split. As relevance maps have to be computed w.r.t. a given class, we restrict our analysis to values of the target variable belonging to the positive class of each MTS instance. The number of relevance attribution maps over which metrics are computed varies across MTS instances, for a minimum of 17 in fMRI and CLIM and 47 for 7ts2h.

A. Temporal Consistency of xAI Methods

The series of experiments reported in this section aims to answer the following questions: (i) Is there an xAI method that is more consistent than others across windows? (ii) Is there a NN architecture that is favored by an xAI method when looking for temporally consistent explanations? (iii) How does the length k of explanations affect their consistency, with respect to a random and theoretical best baseline? Fig. 3 illustrates the temporal consistency of different xAI methods per NN architecture and dataset family, when varying the length k of explanations. The minimal (or optimal) consistency metric corresponding to a deterministic explanation is plotted in dashed font as a reference. Note that given a MTS X ∈ R N ×T , the largest value for k is thus given by N × T .

The Feature Occlusion (FO) method is evaluated only for k that are multiples of the number of TS variables. Note that the relevance attribution of all TS variables in a given timestamp is identical since the occlusion patch covers all features in a timestamp. In this sense, the FO method estimates temporal importance rather than feature importance. Hence, we consider only k values that correspond to selecting all features in the most salient timestamps.

We first focus on the temporal consistency of xAI methods (i). As we can observe in Fig. 3 with the exception of dCAM and XCM, explanation methods exhibit similar trends for the same NN architecture and dataset family. The two gradient methods IG and GS have nearly the same consistency, which we can attribute to the high similarity between the two methods. FA is similar but with better consistency at k < 5, showing that the method's higher saliency features tend to vary less.

The lower consistency metric value that DM shows around k = 4, 5, 6 is a consequence of the sparsity loss term optimized by the method. It forces an adjustable number of relevance coefficients to zero and the rest to 1 (and selects by itself this adjustable parameter). There is thus a reduced amount of randomness in its top k features as k gets close to this threshold. Method TS performs close to or better than FA, except for LSTMs. As for FO, it is more consistent than other methods due to ranking timestamps and not individual features, which decreases the entropy of the associated distribution. As expected, methods with coalition constraints (FO, DM, TS) are generally more consistent than univariate alternatives, except FA.

The randomness of the FP method could be attributed to the batch computation. To permute features, FP explains at the same time a batch of windows. The values of the permuted features come from the same batch. Since our batch size is 32, the distribution of permuted features varies in practice across batches.

The dCAM method shows a distinct behavior, especially at high k, due to the internals of its architecture. We observe that CAM coefficients regularly have low variance in the last time steps of the window. This side effect is likely due to the input padding, which consists of adding zero coefficients to the input of each convolutional layer. It derives that the feature relevance in the last timestamps of each window is consistently low.

As for the XCM method in 7ts2h, we observe that for about 60% of the models, the saliency maps produced consist of zero-valued coefficients. XCM essentially runs two parallel stacks of convolution layers with temporal (1D) and spatiotemporal (2D) kernels respectively. To build explanations, XCM only uses the spatiotemporal stack. We note that the temporal pipeline has wider kernels and three to four times more trainable weights. Obtaining zero-valued explanations denotes that the model has learned to use only this part. For the other 40% of the models, we observe a consistency at around four for k = 1. In fMRI and CLIM, the produced saliency maps are different from zero for more than 90% of the models.

Then, we investigate the effect of different NN architectures (ii). Our first observation is that up to k = 10, explanations of LSTM models are systematically more consistent than their counterparts computed on TCNN and Transformer (except FA on Transformer). For k = 1 in particular, for each xAI method and dataset family, the consistency metric on LSTM is at least one point below the consistency metric on the other types of models. The consistency at low k (LSTM, XCM, Transformer, DCAM) is ordered by the number of trainable parameters in our NN architectures. More complex models might distinguish more finely between different inputs, leading to a higher variation of explanations. A more extensive study of the impact of the model size on the temporal consistency of explanations is left as future work.

We finally study the consistency of explanations for different k values (iii). Our general observation is that the more features in the explanation, the less consistent the explanation is, across all NN architectures and xAI methods. For LSTM especially, the explanation methods get farther away from the minimal consistency as k increases. This indicates that the top few features are less randomly distributed than the rest ones, especially when it comes to the last 15 features, where consistency becomes close to a fully random explanation. We observe the ability of the xAI methods to report somewhat consistent most salient features while having very random small coefficients for others. Another noteworthy result is that for the LSTM model on 7ts2h with FA, explanations of length one have about the same entropy as a uniform distribution of support two. This means that the explanation methods are mostly concentrating on one or two features.

Summary: xAI methods produce more consistent explanations of LSTMs, which could be due to their small size. On the same NN architecture and dataset family, coalition-based method tend to be more consistent than univariate alternatives, with FA as an exception for low explanation sizes. dCAM and XCM have side effects leading to distinct behaviors.Finally, we observe that explanations quickly become random as we include features of lower relevance.

B. Causal Relevance of xAI Methods

The series of experiments reported in this section aims to identify the extent to which xAI methods are able to discover the true causes of the predicted variable in the causal graph. Similarly to the previous section, we are interested in the comparative performance i) of the different xAI methods, ii) of the different NN architectures, iii) and how performance metrics are affected as k grows, with respect to a random and theoretical best baseline.

We measure for each xAI method the precision and recall per explanation length k using as ground truth the respective causal graph. We first compute these metrics for individual MTS instances and predicted variables and average them over the time t. Then, we average across MTS instances and targets to plot in Fig. 4 a single curve per dataset family and NN architecture. We add two baselines: max and random, corresponding respectively to the maximal precision or recall that can be obtained, and the expectation of these metrics if the explanation was a uniformly randomly chosen feature set.

On TCNNs and LSTM, the xAI methods that achieves the best precision and recall are FA and TS on most datasets with close performances. DM seems well suited to the fMRI dataset. The sparsity constraint around k = 4, 5, 6 is effective in selecting the best features. Gradient-based methods IG and GS obtain near identical results. FP exhibits a performance either above or below IG and GS depending on the NN architecture and dataset family. dCAM explanations are extremely inaccurate: similarily to consistency, dCAM explanations miss the ground truth causes as they are usually situated within the last few timestamps of the saliency map that have systematically low coefficients. On the opposite, FO has a high precision and recall starting from the second marker on CLIM, since FO often identifies the two timestamps in which the ground truth causes are located. XCM produces close-torandom explanations. As seen in the previous section, XCM explains the less relevant spatiotemporal part of the network. The problem of 0-valued saliency maps is a particular case of this partial model explanation. Finally, we remark that in general, all xAI methods exhibit a similar standard deviation of the precision across all datasets and model combinations (differing by less than 0.1 to 0.05), decreasing as k increases. On top of it, we remark no clear ordering on all (dataset, model) pairs, except that GS and IG are nearly identical and have the highest variance on 7ts2h and CLIM.

Considering model differences (ii), we see that LSTM dominates over the others on all methods and dataset families. At k = 1, each xAI method dominates its counterparts in other models by at least 0.1. This is a surprising result as the AUCROC of the models used for explanation are not higher in LSTM compared to XCM and Transformer (see Table II). Nevertheless, we observe that the LSTM models obtain an AUCROC higher than 0.7 on more MTS instances. We conjecture that even if two architectures have similar performance, the architecture that can be applied reliably on many MTS instances will focus more on the true causes. The recurrent structure of the LSTM permits it to outperform other models since it is capable of better learning cause and effect with a short lag (1 or 2).

Regarding the effect of the explanation length (iii), we observe that for k = 1 LSTM paired with FA achieves a close to 0.9 precision for the most salient feature in 7ts2h, and above 0.75 for CLIM. The recall of explanations methods excluding dCAM and XCM shows that up to k = 8, the k-th feature is identifying true causes better than a random guessing. Afterward, added features bring little information until the recall reaches the random baseline. It becomes clear that only the top few features include the true causes of our target variable, while the rest of the explanation contains noise.

Finally, we observe a lower starting precision on fMRI than CLIM or 7ts2h. Table II reports a lower number of successful models for fMRI than for CLIM and 7ts2h. We reach a similar conclusion as in the analysis of per model type behavior: the ability of a model to achieve acceptable performances on a diverse set of MTS of similar dynamics is linked to model reliance on true causes.

Summary: Overall, methods other than dCAM and XCM succeed in recovering at least one of the true causes. Still, none of them are able to recover all causes, as features outside of the top 10 bring little information. FA and DM have a small lead on other xAI methods depending on the model and the dataset. Finally, models that perform well on more MTS instances, as LSTM does, have explanations closer to the causal graph.

C. Predictive Performance of xAI Methods

The last part of our analysis concerns the predictive performance of surrogate models trained only over the explaining features. We seek to know if the top salient features are sufficient predictors for training surrogate models and compare their performance to a baseline model where only the direct causes from ground truth are included as input. Thus, we focus on i) whether the features that are most frequently included in explanations can be used to predict the target TS variable as the original model, ii) how the performance of these surrogate Fig. 5. Critical difference (CD) diagram of the AUCROC of the EO, CO, and Full models (Nemenyi test CD=0.518). The scale denotes the average ranking of each method, and bold lines denote non significant differences. models compares to models exclusively built over the direct causes set of the target TS variable.

Both questions i) and ii) relate to the relative performance of the explanation-only models with FA, IG, GS, FP, DM, XCM, dCAM, TS (EO models), along with cause-only models (CO models) and Full (original) models. We build a critical difference diagram using the AutoRank library [START_REF] Herbold | Autorank: A python package for automated ranking of classifiers[END_REF] to test if there are any significant differences between the 10 models [START_REF] Demšar | Statistical comparisons of classifiers over multiple data sets[END_REF]. We rely on the non-parametric Friedman test [START_REF] Demšar | Statistical comparisons of classifiers over multiple data sets[END_REF] as an omnibus test, and reject the null hypothesis that all methods come from the same distribution (with a significance level of 5%). Next, we use the post-hoc Nemenyi test, in order to compare the methods in pairs. Two methods are significantly different when their average ranks differ more than a critical distance (CD) of 0.585 (at a significance level of 0.05). We conduct our test on the 536 data points that are formed by successfully trained models (Table II). Fig. 5 depicts the statistically significant models in a critical difference diagram.

We observe that i) EO models exhibit on average a poor predictive performance compared to Full models, and ii) that CO models outperform on average not only the EO models but also the Full models. xAI methods definitively fail to identify a minimal set of features necessary to the prediction task, unlike causal predictors leading on average to a higher model quality. Fig. 6 finally depicts how often EO models exhibit a worse performance than the CO models. Only in a few MTS instances, we observed that EO models are better than their CO counterparts. We speculate that this difference is attributed to two factors. The first is the independence relations not entailed by the causal graph (i.e. the faithfulness of the causal graph to the distribution it explains). Specifically, it is possible that some TS variables in our true causes have low causal strength. The second factor is specific to the CLIM dataset, where the causal ground truth is observed to be slightly inaccurate for some particular MTS instances.

Summary: EO models clearly underperform compared to the Full and CO models. This demonstrates not only the predictive power of direct causes w.r.t. the salient features but also w.r.t. to the full models trained on the original set of features. As a matter of fact, full models exploit correlations between all available features to gain predictive performance that are unlikely to generalize in different environments, especially in high-dimensional settings. Fig. 6. AUCROC of explanation-only restricted (EO) models compared to the causal-only restricted (CO) models, for each predicted variable and xAI method. Points above the diagonal represent datasets where the EO outperforms the CO model.

VI. CONCLUSION

The main conclusion drawn from our experiments is that relevance attribution methods do not produce consistent explanations across time and seldom discover the true causes of the predicted variables. As expected, the intrinsic measures underlying post-hoc methods are purely associational and usually expose potentially spurious and misleading correlations in input data used to train a NN model. For this reason, surrogate models build over salient features systematically underperform w.r.t. MTS classifiers build over the full feature space. On the contrary, sparse NN models build over causal explanations generalize much better than the full models.

We let as future work the study of how the size/complexity of the models influences the consistency of the xAI methods. We also remark on the need of adding more stochastically generated dataset families, ideally with controllable noise and diverse dynamics. This would allow experimenting on a link between the discovery of true causes and the robustness of the model type (to which extent the same NN architecture can learn on different MTS instances of shared dynamics). Dataset families with a large number of covariates could be added to study if NN can learn sparse causal information despite the curse of dimensionality, and if sparsity constraints could make the networks focus on causal variables. Another direction of improvement would be enlarging the causal evaluation of xAI methods to counterfactual explanations [START_REF] Ates | Counterfactual explanations for multivariate time series[END_REF], to forecasting and multi-class models interpretability methods [START_REF] Ozyegen | Evaluation of interpretability methods for multivariate time series forecasting[END_REF].

Fig. 1 .

 1 Fig. 1. (a) Synthetic MTS X, (b) Heatmap of a Relevance Attribution matrix computed on X. The highest value R 6 t is highlighted.

Fig. 2 .

 2 Fig. 2. Causal graphs underlying datasets of the three families: nodes are labeled by TS variables, directed edges represent direct causal relationships. Dashed lines indicate the presence of latent confounders (i.e., non-observable TS variables) that cause both linked observable TS variables).

Fig. 3 .

 3 Fig. 3. Mean temporal consistency when varying k per dataset family and NN architecture. Lines represent average consistency across all target TS variables in a dataset and areas the 95% confidence interval of predictions.

Fig. 4 .

 4 Fig. 4. Mean Precision (a) and Recall (b) as a function of k per NN architecture and family of datasets. Lines represent metrics at a given k averaged over all predicted features and instances and areas the 95% confidence intervals of predictions.

TABLE I DATASETS

 I CHARACTERISTICS: FMRI AND CLIM INSTANCES COMMIT TO DIFFERENT CAUSAL GRAPHS (IN-DEGREES ARE AGGREGATED).

	Dataset	Instances	Variables	Timestamps	Avg in-degree Max lag
	7ts2h [23]	10	7	4000	1.8	1
	fMRI [31]	17	5	200 to 5000	2.0 to 2.6	1
	CLIM [29]	200	5	250	2.0 to 4.4	2

1 https://captum.ai/api/

TABLE II AVERAGE

 II AUCROC AND NUMBER OF MODELS WITH AUCROC ≥ 0.7 PER NN ARCHITECTURE AND DATASET FAMILY.

		7ts2h	CLIM	fMRI
	LSTM	0.916 (51) 0.784 (116) 0.817 (45)
	XCM	0.901 (43) 0.801 (49)	0.826 (19)
	DCAM	0.867 (38) 0.786 (59)	0.793 (26)
	Transformer	0.903 (51) 0.792 (23)	0.906 (16)


Lino Da
email: damien.vacher@laept.univ-bpclermont.fr

DEFINITION OF A NEW LEVEL ONE TEST CASE MEASUREMENTS OF EQUILIBRIUM RADIATION FROM AN INDUCTIVELY COUPLED PLASMA IN THE NEAR-UV TO NEAR-IR SPECTRAL REGION FOR A TITAN-TYPE N 2 -CH 4 MIXTURE. PRELIMINARY RESULTS

A test case, using a similar methodology and experimental set-up than previous test case TC5-Level 1 [1;2], is presented in this paper. An inductively coupled plasma torch, working at atmospheric pressure, is used to create N 2 -CH 4 Titan-like plasma (98%N 2 -2%CH 4 ). The operating frequency and power are 64 MHz and 3 kW respectively. This kind of apparatus allows obtaining plasma in chemical and quasi-thermal equilibrium The spectral measurements cover the [300-800] nm range and are performed inside the induction coil. Each interesting spectrum is calibrated and compared to the line-by-line spectral code SPARTAN used for the simulation of the radiative emission of entry-type plasma.

Finally, a discussion is proposed about the nucleation phenomenon which is occurred in the ICP torch with the N 2 -CH 4 plasma. Preliminary studies show the synthesis of nanostructured carbon on the quartz tube.

INTRODUCTION

Different levels of test-cases have been defined in the frame of the ESA-CNES radiation working group. Such test-cases are useful for the evaluation and validation of radiation models used for atmospheric entry applications.

Test case n° 1 is devoted to the validation of models of line-by-line spectroscopic codes and of spectral data used for line-by-line radiation calculations. As the simplest case is the modelling of plasma in thermo-chemical equilibrium, an ICP torch, which allows obtaining such conditions, is used. Test case 1 was initially intended to validate measurements in air plasma. This new level of test case 1 proposes a spectroscopic study of plasma formed with a Titan-like atmosphere gas mixture, which will be useful for the validation of spectral databases (SPARTAN code).

The outline of the paper is as follows. In section 2, all the characteristics of the experimental set-up are given. In section 3, after a discussion concerning the sustained conditions of the plasma when methane is injected into a pure nitrogen plasma, experimental spectra in the range [320-840] nm are presented and an estimation of the axial temperature is given. In section 4, a discussion is proposed about the nucleation phenomenon which is occurred in the ICP torch with the N 2 -CH 4 plasma. Preliminary studies show the synthesis of nanostructured carbon on the inner surface of the quartz tube. The chemical equilibrium composition of a N 2 -CH 4 Titan-like plasma (98%N 2 -2%CH 4 ) is given through the calculation code based on the Gibbs free energy minimization. Finally, the perspectives of this preliminary study will be given in order to propose this new level one test case.

EXPERIMENTAL SET-UP

The ICP-T64 torch located at the L.A.E.P.T. (Thermal Plasmas and Electrical Arc Laboratory) in Clermont-Ferrand, France, is a classical ICP torch able to work with different kinds of plasma gas (air, argon, CO 2 , N 2 and gas mixtures). A seven-turn induction coil, cooled by air, is used to ignite (with a metallic wire) and sustain the N 2 -CH 4 plasma.

The main features of the experimental setup and operating conditions are reported in Fig. 1 and Tab. 1, respectively. The plasma is generated through the induction coil by a radio frequency (RF) of 64 MHz delivering a power up to 3 kW. The plasma is confined within a 28-mm quartz tube. The plasma gas is injected at a fixed rate of 0.2 g.s -1 .

The optical set-up, placed at 34.3 cm from the plasma axis, leads to a spatial resolution of 1mm. Spectral lines intensities are measured with a 0.50 m focal length Czerny-Turner monochromator connected to a CCD detector (1242 1152 pixels, each pixel having a width of 22.5 micrometers). An 1200 grooves/mm grating is used. Its apparatus function app λ ∆ , assimilated to a Gaussian profile with a Full Width at Half Maximum (FWHM), is calculated from the relation (Eq.1) :

1 - = ∆ D f FWHM p app λ (1)
where FWHM is expressed in pixels, D -1 represents the gratings dispersion (1.56 nm/mm) and f p defines the pixel dimension (22.5 µm). The calculated value of app λ ∆ gives 0.14 nm for an entrance slit equal to 100 µm. The measured emission spectra must be corrected from the spectral response of the optical device, which include lens, monochromator and CCD detector. So, a calibration in intensity, between 300 and 800 nm, is necessary to take into account all these effects. A tungsten lamp (PHILIPS WI14) is used to cover the range [300 -800] nm. The calibration procedure is realized in the same conditions that ones considered for the experimental spectral acquisitions.

EXPERIMENTAL SPECTRUM

Optical emission spectroscopy is the simplest non-intrusive method for plasma diagnostics. The analytical zone is limited between the fourth and the fifth induction coil where the luminous response is the greatest. Each spectrum corresponds to one acquisition of one second and the electronic noise is automatically subtracted.

One particularity of the Titan atmosphere plasma created in the ICP-T64 torch is the difficulty to have a sustainable plasma. At first, the plasma was ignited with pure nitrogen, and when the methane was injected the plasma generally was cut. It has been observed also that there was a loss of injected power, about 1 kW, when the methane was added. This effect is mainly due to the automatic adaptation of impedance but it is worth noting that there is no problem of extinction when the plasma changes to Titan atmosphere to pure nitrogen. The studies have to be improved to understand this unusual behavior of this kind of plasma. First explanations can be given by considering the different physical properties of the two plasmas (pure N 2 and Titan composition). As the viscosity doesn't seem to be an influent parameter, the two others may explain the difficulty to have a viable plasma when methane is added. Concerning the enthalpy (Fig. 5), no significant differences appear and for instant, no clear conclusion can be brought on the specific behaviour of the plasma when methane is added. More thorough studies must be realized.

Fig. 6 reports a rough spectrum recording from the plasma formed with a Titan-like atmosphere in the range [300 -800] nm. The CN molecular emission is predominant and the C 2 emission is also present. No atomic line is observed in the spectral domain. 8 shows the result of the application of Abel inversion. It must be noted that each intensity radial evolution of wavelength is interpolated by cubic splines and fitted by a polynomial of degree 9 before applying the Abel inversion.

It can be noted that the application of Abel inversion is easily practicable for the CN system whereas for the C 2 one, very big uncertainties are present. The presence of a plateau of intensity in the axis of the plasma leads to difficulties in order to apply the Abel inversion. In order to estimate the value of temperature, the spectra of the C 2 Swan system and the CN violet system has been considered and the SPARTAN code [START_REF] De | An adaptative Line-by-Line-Statistical Model for Fast and Accurate Spectral Simulations in Low-Pressure Plasmas[END_REF]used. Fig. 9 reports the best agreement obtained by the comparison between the experimental spectrum at the centre of the plasma for the C 2 Swan system and the synthetic one. As it is considered that the plasma is at or close to thermal and chemical equilibrium, two temperatures have been estimated in order to have the best agreement, the rotational temperature (T r ) and the vibrational one (T v ). The obtained values are T r =3200K and T v =3700K.

Fig. 10 reports the results concerning the CN violet system. When the experimental spectrum is observed, one can conclude that there is a strong absorption of the lines, especially for the first band head (0-0). So, in order to take into account this absorption, the following calculations are made : -Estimation of T v and T r followed by the calculation of emission ε n and absorption α n coefficients using the line-by-line code SPARTAN; -Calculation of the slab emitted intensity (Eq. 2)

(
)

)) exp( 1 / l I n n n × - × = α α ε (2)
-Convolution with a gaussian apparatus function simulating the slit. The agreement between simulation and experience is not ideal, as it can be seen in Fig. 10. The best gives T r =3200K and T v =3700K with a somehow large uncertainty of 500K. However, this first estimation of temperature, based on the molecular spectra, leads to acceptable results and is a good starting point for future investigations. It is worthy to note that owing to the large uncertainty on the temperatures, it is possible that a thermal equilibrium (T r =T v ) might be reached in the plasma.

DISCUSSION ABOUT THE NUCLEATION PHENOMENON

When the ICP torch works with a Titan like atmosphere, it is observed inside the quartz tube the formation and deposition of "black suit" on the injector (Fig. 11) and also on the walls of the quartz tube. It can be noted that the rotation of the deposit follows the movement of vortex of the gas flow. In order to determine what kind of particules was formed in the Titan plasma, an analysis of the suit has been made by using Scanning Electron Microscopy (SEM) with field effect (CASIMIR society). Fig. 12 shows a zoom of the particular aspect of the suit. If the zoom is further increased, spherical particles can be isolated with a diameter of around 35 nanometers. The studies which have been realized on carbon nano-particles synthesis [START_REF] Moreno | Carbon nanoparticles synthesis by gas phase non-equilibrium plasma[END_REF] shows similar images to that of Fig. 12. This particular aspect which represents the type of the carbon compounds nano-structured is called "ruffled paper structure". An X-ray analysis must be done to determine the true nature of the particules and to verify if it is pure carbon.

From the calculation of the plasma composition based on the Gibbs free energy minimization, one can give a first explanation about the presence of carbon in the plasma torch (Fig. 13). As it can be observed in the Fig. 13, there is the presence of a graphite phase [C(S)] until a temperature of 2700K. So, it confirms the formation of C on the walls of the quartz tube and on the injector which is placed just under the first inductor coil. By comparison with the Mars plasma which has been studied in a previous work and where no formation of carbon has occurred, the Fig. 14 shows also the fact that there is no formation of carbon specie. This absence is mainly due to the presence of O 2 .

CONCLUSION

The first investigations on the study of a Titan like atmosphere plasma formed with an inductively coupled plasma torch lead to interesting observations and results, especially the formation of carbon species. It was also shown that the addition of methane strongly disturbed the plasma and explanations about this change of behaviour have to be carried out.

Temperatures have been estimated but only at the centre of the plasma. Further works have to be done in order to determine the radial evolution of temperature and if a weak departure from equilibrium exists in this kind of plasma. For this, the addition of a very small percentage of argon in the plasma would be interesting as we could have atomic lines to determine the atomic excitation temperature.

A new optical diagnostic has to be added to the experimental set-up in order to determine the electronic density (Laser Interferometry). Moreover, it would be of interest to realize kinetic calculations to confirm that the chemical equilibrium is reached inside the inductor where the optical measurements are recorded (as it has been done for the study of the Mars plasma).

The aim of the future investigations is so to propose a complete level-one test case for the validation of radiative codes simulating equilibrium radiation from a Titan type plasma in equilibrium conditions in the near UV to near IR range.

Figure 1 .

 1 Figure 1. Plasma torch and detection experimental set-up.

Figure 2 . 2 Figure 3 . 2 Figure 4 .

 22324 Figure 2. Evolution of the electrical conductivity as a function of temperature.

Fig 2 - 3 Figure 5 .

 235 Fig 2-3-4 report the differences observed for the electrical conductivity, the viscosity and the thermal conductivity of the two plasmas between 2700K and 10000K.

Figure 6 .

 6 Figure 6. Experimental emission spectra for a Titan atmosphere plasma at atmospheric pressure.

Fig. 7

 7 Fig.7represents different intensity radial evolutions of CN and C 2 peaks respectively and Fig.8shows the result of the application of Abel inversion. It must be noted that each intensity radial evolution of wavelength is interpolated by cubic splines and fitted by a polynomial of degree 9 before applying the Abel inversion.It can be noted that the application of Abel inversion is easily practicable for the CN system whereas for the C 2 one, very big uncertainties are present. The presence of a plateau of intensity in the axis of the plasma leads to difficulties in order to apply the Abel inversion.

Figure 7 .Figure 8 .

 78 Figure 7. Emission profile measured over the plasma diameter for the first four band heads of the CN violet system (∼388 nm) and the first two band heads of the C 2 swan system (∼510 nm).

Figure 9 .

 9 Figure 9. Simulation of the radiative emission of the C 2 Swan system over the ICP torch axis.

Figure 10 .

 10 Figure 10. Simulation of the radiative emission of the C N violet system over the ICP torch axis for the ∆υ=0 transitions (best agreement obtained).

Figure 11 .

 11 Figure 11. Photographs of the suit deposited on the injector (left) and on the walls of the quartz tube (right).

Figure 12 .

 12 Figure 12. Photographs of the suit analysed by SEM with a scale of 200nm.

Figure 13 .Figure 14 .

 1314 Figure 13. Chemical composition (molar fraction) of a N 2 (98%)-CH 4 (2%) plasma at equilibrium at atmospheric pressure as a function of temperature.

Table 1 .

 1 Main features of the experimental set-up.

	Inductively coupled plasma
	Manufacturer/type: Défi Systèmes/ICP T64
	Power supply: 64 MHz, 3 kW
	Tuning: Automatic adaptation
	Inductor: Seven-turn air-cooled coil
	Plasma gas flow rate: 0.2 g.s -1
	Operating pressure: atmospheric pressure
	Torch: 28 mm internal diameter quartz tube
	Optical set-up
	Spectrometer: Chromex 500 IS
	500 mm focal length
	Czerny-turner mounting
	Entrance slit: e = 100 µm Gratings: 1200 grooves.mm -1
	Detector: CCD EEV 1152×1242 pixels
	Spatial resolution: 0.5 mm
	Apparatus function: 0.14 nm
	HF generator
	64 MHz
	DEFI Systèmes
		CCD
	Plasma gas
	N 2 -CH 4 (98%-2%)	ST 138 controller

ACKNOWLEDGMENTS

This research work has been carried out thanks to the financial support of the RayHen ANR and the support of the European Space Agency ESA in the scope of the AURORA program.


Greedy reduction of

Conclusion:

❑ Synthetic inactivation of the mismatch repair system resulted in an increased mutation rate and enabled the selection of competitive genome-reduced isolates. ❑ This establishes an engineering framework for efficient genome streamlining of chassis strains for synthetic biology applications.

INTRODUCTION

We investigate the consequences of genome reduction in Bacillus subtilis, which has many advantages for genome engineering including natural competence, efficient chromosomal replacement, existing repertoires of non-essential regions and genes [START_REF] Tanaka | Building the repertoire of dispensable chromosome regions in B. subtilis entails major refinement of cognate large-scale metabolic model[END_REF] , and availability of libraries of parts for fine-tuning gene circuit expression [START_REF] Guiziou | A part toolbox to tune genetic expression in Bacillus subtilis[END_REF] . As such, B. subtilis is an ideal candidate for the development of a well-characterized minimal chassis in which complex metabolic and regulatory circuits can be implemented and fine-tuned. First, we applied a greedy genome reduction strategy in which deletions in large dispensable regions were serially combined starting from the largest. At each step, the fittest clonal population was used for the next deletion cycle. This strategy generated a synthetic evolutionary tree of 298 strains with reduced genomes, the smallest genome (of 2.73 Mb) having lost 1.48 Mb in size [START_REF] Dervyn | Greedy reduction of Bacillus subtilis genome yields emergent phenotypes of high resistance to a DNA damaging agent and low evolvability[END_REF] .

Delivery:

❑ Library of 298 strains with genomes reduced to as small as 1.48 Mb in size. Method: ❑ Deletions were accumulated through a two-step deletion cycle.

GREEDY REDUCTION OF THE B. subtilis GENOME

❑ We built a library of 146 genome-reduced strains, each with a single deletion of a dispensable chromosome region, collectively covering ∼76% of the B. subtilis 168 genome 1 . We selected the 48 largest dispensable regions (from ∼131 kb down to ∼16 kb), covering 57% of the genome, for stepwise deletion (Figure 1). ❑ Deletions were accumulated through a two-step cycle (Figure 2):

(i) Deletion of a chromosomal region by homologous replacement with a selectable cassette.

(ii) Eviction of the cassette via spontaneous pop-out recombination to generate marker-free deletions, preparing genome-reduced strains for further cycles. ❑ From MGP052 p/o (p/o for pop out) two separate lineages were pursued accumulating 10 and 17 deletions, respectively. Finally, we obtained MGP184 (∼2.87 Mb) and MGP254 (∼2.73 Mb), two strains respectively reduced by ∼1.34 and ∼1.48 Mb relative to the B. subtilis 168 strain.

LOSS OF FUNCTION AFFECTS FITNESS OF GENOME-REDUCED STRAINS

∼ 24% showed increased growth rates, while ∼ 40% exhibited decreased growth rates compared to the Master strain (MS).

∼1% of these strains displayed significantly increased biomass production, whereas 72% exhibited significantly decreased biomass production, possibly due to challenges in carbon source consumption or increased cell maintenance.

❑ Cell fitness decreased with the accumulation of deletions in the genome. ❑ Accumulation of deletions globally led in cumulated synthetic-sick interactions within each lineage.

❑ When comparing MGP181 and MGP254 with MS (using a fold change cut-off of 1.5), ∼ 83% and 72% of core genome genes (present in all genome-reduced strains) showed no significant differential expression.

Etienne Dervyn 1, † , Anne-Gaëlle Planson 1, † , Kosei Tanaka

GENOME-REDUCED STRAINS EXHIBIT LOW SPONTANEOUS MUTAGENESIS

❑ Low mutation frequency in MGP229 p/o (Table 1) relies entirely on a functional mismatch repair (MMR) system. ❑ The MMR system might be more efficient in genome-reduced strains, drastically reducing spontaneous mutagenesis and likely limiting evolution. Genome-reduced strains have limited potential for short-term adaptive evolution. MS cat ΔmutS 2.4 x 10 -6 ± 1.0 x 10 -6 MGP229 p/o ΔmutS 3.5 x 10 -6 ± 1.0 x 10 -6

GENOME-REDUCED STRAINS ARE HIGHLY RESISTANT TO THE DNA DAMAGING MITOMYCIN C

These findings highlight the transcriptional network's robustness against extensive genome reduction and suggest that synthetic-sick interactions may stem from deregulation of only a few genes.

❑ Transcriptome analysis revealed a common feature in our genome-reduced strains: upregulation of genes involved in cell responses to oxidative stresses.

CONCLUSION

Design and construction of efficient cell factories shall be envisioned as an engineering framework coupling reductive and evolutive cycles with high-throughput screening. Reductive cycles encompass successive deletions of unwanted genomic regions in a context of low evolvability. Evolutive cycles include a condition-dependent increase of spontaneous mutagenesis for experimental evolution or the integration in a chromosomal landing-pad of a library of DNA fragments from the Bacillus species pangenome. This blueprint will allow selecting efficient cell factories at each cycle of deletion, before entering the next cycle and pave the way for innovative construction of synthetic chassis. ❑ Over 272 strains carrying a single deletion of dispensable regions 1 (Figure 3):

The DNA-damaging agent mitomycin C (MMC) induces double-stranded DNA breaks.

LexA, which controls the SOS response to DNA damage, is constitutively active in genomereduced strains.

The upregulation of the LexA regulon may contribute to the higher MMC resistance (Figure 4).

Figure 1 :

 1 Figure 1: A synthetic evolutionary tree of 298 strains with genomes reduced by up to 1.48 Mb.

Figure 2 :

 2 Figure 2: The deletion cycle.A deletion mutant strain (Δn) is transformed with DNA carrying a deletion of rank n + 1 linked to a selectable marker M. Upon its integration by homologous recombination, the deletion is integrated into the genome. Marker eviction occurs spontaneously by pop-out recombination and is selected by restoration of an antibiotic resistance.

Figure 3 :

 3 Figure 3 : High-throughput phenotyping of the 298 genome-reduced strains. Growth rate is negatively correlated to the deleted-genome size. Dots represent the single interval deletion ( ) and multiple interval deletion ( , , , ) strains from MS ( ). Edges represent the lineages: from MS to MGP052 ( ), from MGP052 to MGP184 ( ) and to MGP254 ( ). The 95% confidence intervals are given as vertical bars.

Figure 4 :

 4 Figure 4: Genome-reduced strain sensitivity to MMC. Cultures were plated on LB agar plates with either 0 or 80 ng/mL MMC and incubated at 37°C for 72 hours.

  1, † , Victor Chubukov 2 , Cyprien Gérin 3 , Sandra Derozier 3 , François Lecointe 1 , Uwe Sauer 2 , Ken-Ichi Yoshida 4 , Pierre Nicolas 3 , Philippe Noirot 1 and Matthieu Jules 1

	HIGHLIGHTS
	Characteristics:
	❑ Growth was progressively impaired with genome reduction.
	❑ Accumulation of synthetic-sick interactions.
	❑ Transcriptional signatures of a constitutive thiol-induced oxidative stress.
	❑ Decreased frequency of spontaneous mutation (up to 100-fold).

Table 1 :

 1 Mutation rates :