Oana Cocarascu

Elena Cabrio

Serena Villata

Francesca Toni

Dataset Independent Baselines for Relation Prediction in Argument Mining

Keywords: Argument Mining, Relation Prediction, Machine Learning Methods Computational Models of Argument H. Prakken

Argument(ation) Mining (AM) is the research area which aims at extracting argument components and predicting argumentative relations (i.e., support and attack) from text. In particular, numerous approaches have been proposed in the literature to predict the relations holding between arguments, and applicationspecific annotated resources were built for this purpose. Despite the fact that these resources were created to experiment on the same task, the definition of a single relation prediction method to be successfully applied to a significant portion of these datasets is an open research problem in AM. This means that none of the methods proposed in the literature can be easily ported from one resource to another. In this paper, we address this problem by proposing a set of dataset independent strong neural baselines which obtain homogeneous results on all the datasets proposed in the literature for the argumentative relation prediction task in AM. Thus, our baselines can be employed by the AM community to compare more effectively how well a method performs on the argumentative relation prediction task.

Introduction

Argument(ation) Mining (AM) is "the general task of analyzing discourse on the pragmatics level and applying a certain argumentation theory to model and automatically analyze the data at hand" [START_REF] Habernal | Argumentation mining in user-generated web discourse[END_REF]. Two tasks are crucial in AM [START_REF] Lippi | Argumentation mining: State of the art and emerging trends[END_REF][START_REF] Cabrio | Five years of argument mining: a data-driven analysis[END_REF][START_REF] Lawrence | Argument mining: A survey[END_REF]: 1) argument component detection within the input natural language text, aiming at the identification of the textual boundaries of the arguments and their classification (claim, premise); and 2) relation prediction, aiming at identifying (support, attack) relations between argumentative components, possibly identified in the first stage. In this paper we focus on the second task. Despite the high volume of approaches tackling the relation prediction task with satisfying results (see [START_REF] Cabrio | Five years of argument mining: a data-driven analysis[END_REF] for an extensive list), a problem arises: these solutions heavily rely on the peculiar features of the dataset taken into account for the experimental setting and are hardly portable from one application domain to another. On the one side, this issue can be explained by the huge number of heterogeneous application domains where argumentative text may be analysed (e.g., online reviews, blogs, political debates, legal cases). On the other side, it represents a drawback for the comparison of the different approaches proposed in the literature, which are often presented as solutions addressing the relation prediction task from a dataset independent point of view. A side drawback for the AM community is therefore a lack of large annotated resources for this task, as most available resources cannot be successfully reused, being highly context-based. Even the employment of pretrained language models (e.g., BERT [START_REF] Devlin | BERT: pre-training of deep bidirectional transformers for language understanding[END_REF]) does not address this issue.

In this paper, we tackle this issue by proposing a set of strong cross-dataset baselines based on different neural architectures. Our baselines are shown to perform homogeneously over all the datasets proposed in the literature for the relation prediction task in AM, differently from individual methods proposed in the literature. Our contribution is to bestow the AM community with a set of strong cross-dataset baselines to compare with in order to demonstrate how well a relation prediction method for AM performs.

We focus on two types of argumentative relations: attack and support, given that the majority of datasets target only these two types of relations. We define neural baselines to address the corresponding binary classification problem, analysing, to the best of our knowledge, all available datasets for this task, ranging from persuasive essays to user-generated content, to political speeches. Given two arguments, we are interested in determining the argumentative relation between the first, called child argument, and the second, called parent argument, using a neural model. For example, the child argument People know video game violence is fake may attack the parent argument Youth playing violent games exhibit more aggression. In our baselines, each of the two arguments is represented using embeddings as well as other features. We propose three neural network architectures for the classification task, two concerned with the way child and parent are passed through the network (concat model and mix model), and an attention-based model. We also explore BERT as an alternative to our baselines: although this is used successfully to boost performances for other tasks in Natural Language Processing, it is generally not competitive for relation prediction with the datasets we consider.

We conduct experiments with a number of datasets, chosen either because they were specially created for relation prediction in AM or because they can be easily transformed to be used for this task. These are: Essays (essay) [START_REF] Stab | Parsing argumentation structures in persuasive essays[END_REF], Microtexts (micro) [START_REF] Peldszus | Joint prediction in MST-style discourse parsing for argumentation mining[END_REF], Nixon-Kennedy (nk) [START_REF] Menini | Never retreat, never retract: Argumentation analysis for political speeches[END_REF], Debatepedia (db) [START_REF] Cabrio | Node: A benchmark of natural language arguments[END_REF], IBM (ibm) [START_REF] Bar-Haim | Stance classification of contextdependent claims[END_REF], ComArg (com) [START_REF] Boltužić | Back up your stance: Recognizing arguments in online discussions[END_REF], Web-content (web) [START_REF] Carstens | Towards relation based argumentation mining[END_REF], CDCP (cdcp) [START_REF] Park | A corpus of eRulemaking user comments for measuring evaluability of arguments[END_REF], UKP (ukp) [START_REF] Stab | Cross-topic argument mining from heterogeneous sources[END_REF], AIFdb (aif) [START_REF] Bex | On logical specifications of the argument interchange format[END_REF][START_REF] Chesñevar | Towards an argument interchange format[END_REF][START_REF] Iyad | The argument interchange format[END_REF][START_REF] Reed | AIF+: dialogue in the argument interchange format[END_REF]. Datasets' statistics can be found in Table 1 1 .

Neural baselines for relation prediction

We use four types of features: word embeddings, sentiment features, syntactic features, computed for both child and parent, and textual entailment from child to parent. We refer to the last three types of features as standard features. Word embeddings are distributed representations of texts in an n-dimensional space. Textual entailment represents the class (amongst entailment, contradiction, or neutral) obtained using AllenNLP2 , a textual entailment model based on decomposable attention [START_REF] Parikh | A decomposable attention model for natural language inference[END_REF]. The features related to sentiment are based on manipulation of SentiWordNet [START_REF] Esuli | Sentiwordnet: A publicly available lexical resource for opinion mining[END_REF] and the sentiment of the entire (child and parent) texts analysed using the VADER sentiment analyser [START_REF] Hutto | VADER: A parsimonious rule-based model for sentiment analysis of social media text[END_REF]. Every WordNet synset [START_REF] Miller | Wordnet: A lexical database for english[END_REF] can be associated to three scores describing how objective, positive, and negative it is. For every word in the (child and parent) texts, we select the first synset and compute its positive score and its negative score. In summary, the features related to sentiment for a text t that consists of n words, W i = 1 . . .w n , are the following: (i) sentiment score (∑ w i pos score(w i )neg score(w i )), (ii) number of positive/negative/neutral words in t, (iii) sentiment polarity class and score of t. Syntactic features consist of text statistics (e.g., number of words) and word statistics with respect to part-of-speech tags (i.e., number of words, nouns, verbs, first person singular, etc.) and lexical diversity (i.e., number of unique words divided by the total number of words in text t).

We describe the three neural architectures we propose for determining the argumentative relation (of attack or support) holding between child and parent. For all, we report only configurations of the architectures and number/size of the hidden layers which performed the best3 . For our models, we use GRUs [START_REF] Cho | Learning phrase representations using RNN Encoder-Decoder for statistical machine translation[END_REF] as they take less time to train and are more efficient.

Concat model (C).

In this model, each of the child and parent embeddings is passed through a GRU. We concatenate the standard features of the child and of the parent. The merged standard vector is then concatenated with the outputs of the GRUs. The resulting vector is passed through 2 dense layers (of 256 neurons and 64 neurons, with sigmoid as activation function), and then to softmax to determine the argumentative relation.

Mix model (M).

In this model, we first concatenate the child and parent embeddings and then pass them through a GRU, differently from the concat model where we pass each embedding vector through a GRU first. We concatenate the standard features that we obtain for the child and for the parent. The merged standard vector is then concatenated with the output of the GRU. From this stage, the network resembles the concat model: the resulting vector is passed through 2 dense layers (of 256 neurons and 64 neurons, with sigmoid as activation function), to be then finally passed to softmax. Attention model (A). Inspired by the demonstrated effectiveness of attention-based models [START_REF] Yang | Hierarchical attention networks for document classification[END_REF][START_REF] Vaswani | Attention is all you need[END_REF], we combine the GRU-based model with attention mechanisms. Each of the child and parent embeddings is passed through a GRU and we compute attention in two directions. We concatenate the standard features of the child and of the parent. The merged standard vector is then concatenated with the outputs of the GRUs. The resulting vector is passed through a single dense layers (128 neurons, with sigmoid as activation function), that is then passed to softmax.

Experimental results

Non-neural baselines. For training we have used the larger datasets, aif, essay, ibm and web. We resampled the minority class from the essay dataset and used our models on the oversampled dataset. We did not use for training the ukp dataset as the parent is a topic instead of an argument. The models were then tested on the remaining datasets, with the average being computed on testing datasets. We report the F 1 performance of the attack class (A) and the support class (S) for the non-neural baselines in Table 2. We used Random Forests (RF) [START_REF] Breiman | Random forests[END_REF] with 15 trees in the forest and gini impurity criterion and SVM with linear kernel using LIBSVM [START_REF] Chang | LIBSVM: A for support vector machines[END_REF], obtained as a result of performing a grid search, as it is the most commonly used algorithm in the works that experiment on the datasets we considered [START_REF] Bar-Haim | Stance classification of contextdependent claims[END_REF][START_REF] Boltužić | Back up your stance: Recognizing arguments in online discussions[END_REF][START_REF] Carstens | Using argumentation to improve classification in natural language problems[END_REF][START_REF] Menini | Never retreat, never retract: Argumentation analysis for political speeches[END_REF][START_REF] Niculae | Argument mining with structured SVMs and RNNs[END_REF]. On top of the standard features used for our neural models, for the baselines we added the following features: TF-IDF, number of common nouns, verbs and adjectives between the two texts as in [START_REF] Menini | Never retreat, never retract: Argumentation analysis for political speeches[END_REF], a different sentiment score nr pos-nr neg nr pos+nr neg+1 as in [START_REF] Bar-Haim | Stance classification of contextdependent claims[END_REF], with all features being normalized.

Neural baselines with non-contextualised word embeddings. Table 3 shows the best baselines for relation prediction in AM. We experimented with GloVe (300-dimensional) embeddings [START_REF] Pennington | GloVe: Global vectors for word representation[END_REF], using pre-trained word representations in all our models. We used 100 as the sequence size as we noticed that there are few instances with more than 100 words. We used a batch size of 32 and trained for 10 epochs (as a higher number of epochs led to overfitting). We report the results using embeddings and syntactic features and the results with all the features presented in Section 2. We also conducted a feature ablation experiment (with embeddings being always used) and observed that syntactic features contribute the most to performance, with the other types of features bringing small improvements when used together only with embeddings. In addition, we have run experiments using two datasets for training to test whether combining two datasets improves performance. During training, we used one of the large datasets (aif, essay, ibm, web) and one of the remaining datasets (represented as blanks in the table ). Amongst the proposed architectures, the attention model generally performs better. Using only a single dataset for training, the model that performs the best is the mix model using all features and trained on the essay dataset. The best results are obtained when using another dataset along one of the larger datasets for training. This is because combining data from two domains we are able to learn better the types of argumentative relations. When using syntactic features, adding micro, cdcp, and ukp does not improve the results compared to using a single dataset for training. Indeed, cdcp has only one type of relation (i.e. support) resulting in an imbalanced dataset, and in ukp, the parent argument is a topic, which does not improve the prediction task. When using all features, micro, com, ukp, and nk do not contribute to an increase in performance. The best performing model is the attention mechanism trained on the web and essay datasets using syntactic features (0.544 macro average F 1 ). 4. Experimental results with F 1 for Attack and for Support relations. XB stands for the number of BERT layers used (X=3,4) and YD stands for the number of dense layers (Y=1,2) used before the final layer that predicts the class. The blanks indicate the training datasets. The Average (Avg) and the Macro (Macr) Avg do not include the results for the training datasets.

beddings [START_REF] Devlin | BERT: pre-training of deep bidirectional transformers for language understanding[END_REF] analyse the entire sentence before assigning embeddings to individual words. We employ BERT embeddings to test whether they bring any improvements to the classification task. While for GloVe vectors we do not need the original, trained model in order to use the embeddings, for the BERT embeddings we require the pre-trained language models that we can then fine tune using the datasets of the downstream task. We try different combinations: using 3 or 4 BERT layers and using 1 dense layer (of 64 neurons) or 2 dense layers (of 128 and 32 neurons, respectively) before the final layer that determines the class. Table 4 shows the results with BERT embeddings instead of GloVe, using feature ablation (syntactic vs all features) and two datasets for training to test whether this can improve performance. The best results are obtained using 4 BERT layers and 2 dense layers (0.537 macro average F 1 ). However, this best BERT baseline does not outperform the best results with the attention model and GloVe embeddings.

Related work

In terms of results on the datasets we have conducted experiments on, most works perform a cross-validation evaluation or, in the case of datasets consisting of several topics, the models proposed are trained on some of the topics and tested on the remaining topics. For essay, an Integer Linear Programming model was used to achieve 0.947 F 1 for support and 0.413 F 1 for attack on the testing dataset using cross-validation to select the model [START_REF] Stab | Parsing argumentation structures in persuasive essays[END_REF]. Using SVM, 0.946 F 1 for support and 0.456 F 1 for attack were obtained [START_REF] Stab | Parsing argumentation structures in persuasive essays[END_REF]. Using a modification of the Integer Linear Programming model to accommodate the lack of some features used for essay but not present in micro, 0.855 F 1 was obtained for support and 0.628 F 1 for attack. On micro, an evidence graph model was used to achieve 0.71 F 1 using cross-validation [START_REF] Peldszus | Joint prediction in MST-style discourse parsing for argumentation mining[END_REF]. On nk, 0.77 F 1 for attack and 0.75 F 1 for support were obtained using SVM and cross-validation [START_REF] Menini | Never retreat, never retract: Argumentation analysis for political speeches[END_REF]. SVM accuracy results on the testing dataset using coverage (i.e. number of claims identified over the number of total claims) were reported in [START_REF] Bar-Haim | Stance classification of contextdependent claims[END_REF] as follows: 0.849 accuracy for 10% coverage, 0.740 accuracy for 60% coverage, 0.632 accuracy for 100% coverage. RF were evaluated on web and aif using cross-validation, achieving 0.717 F 1 and 0.831 F 1 , respectively [START_REF] Carstens | Using argumentation to improve classification in natural language problems[END_REF]. Structured SVMs were evaluated in a cross-validation setting on cdcp and ukp using various types of factor graphs, full and strict [START_REF] Niculae | Argument mining with structured SVMs and RNNs[END_REF]. On cdcp, F 1 was 0.493 on the full graph and 0.50 on the strict graph, whereas on ukp, F 1 was 0.689 on the full graph and 0.671 on the strict graph. No results on the two-class datasets were reported for db, com, and ukp. The results on ukp treat either supporting and attacking arguments as a single category or consider three types of relations: support, attack, neither. The latter type of reporting results on three classes is also given on the com. Some other works have started investigating the dataset independence in AM. [START_REF] Opitz | Dissecting content and context in argumentative relation analysis[END_REF] showed how models may overlook textual content when provided with the context surrounding the span by relying on contextual markers for predicting relations and tested their method on the essay dataset. [START_REF] Lin | Lexicon guided attentive neural network model for argument mining[END_REF] integrated (claim and other domain) lexicon information into neural networks with attention tested on ukp. [START_REF] Kuribayashi | An empirical study of span representations in argumentation structure parsing[END_REF] experimented with span representations, originally developed for other tasks, on the essay dataset. Other works have used contextualised word embeddings for relation prediction in AM [START_REF] Durmus | Determining relative argument specificity and stance for complex argumentative structures[END_REF][START_REF] Reimers | Classification and clustering of arguments with contextualized word embeddings[END_REF]. More recently, [START_REF] Ein-Dor | Corpus wide argument mining -A working solution[END_REF] proposed and tested on ukp an argument retrieval system.

Discussion and Conclusion

Dataset independence is one of the biggest challenges in AM. An AM model for relation prediction trained on every individual dataset considered in this paper would perform better than any general baseline on that dataset. We believe AM model would require leveraging a diverse corpus to be of use in a real-world system. Most works have previously focused on a moderate-sized corpus distributed across a small set of topics [START_REF] Ein-Dor | Corpus wide argument mining -A working solution[END_REF]. This paper is a step towards the applicability of AM techniques across datasets. Our baselines perform homogeneously in terms of average over all existing datasets for relation prediction in AM while using generic features. We propose as baseline the model that performed the best, with the baseline using attention mechanism with GloVe embeddings and syntactic features trained on the web and essay datasets (0.544 macro average F 1 ). The results for the attack class are generally worse than those for support as the datasets that are used in training ibm) have fewer instances for the attack class than for support (see Table 1). The datasets differ at granularity: some consist of pairs of sentences (e.g., ibm) whereas others include pair of multiple-sentence arguments (e.g., nk). Additionally, the argumentative on relations can be domain-specific and their semantic nature may vary between corpora (e.g., com). We considered the unified task of determining support or attack between any two texts.

Embeddings represent the differentiating feature for the models we experimented with. Whilst word embeddings are often used as the first data processing layer in a deep learning model, we employed TF-IDF features for the non-neural models that we considered as baselines. Other works that address the task of relation prediction make use of features specific to the single dataset of interest, making it difficult to test those models on other datasets. For instance, for the essay dataset, [START_REF] Stab | Parsing argumentation structures in persuasive essays[END_REF] use structural features such as number of preceding and following tokens in the covering sentence, number of components in paragraph, number of preceding and following components in paragraph, relative position of the argument component in paragraph. For the other datasets, [START_REF] Stab | Cross-topic argument mining from heterogeneous sources[END_REF] use topic similarity features (as the parent argument is a topic), [START_REF] Menini | Never retreat, never retract: Argumentation analysis for political speeches[END_REF] use the position of the topic and similarity with other related/unrelated pair from the dataset, keyword embeddings of topics from the dataset. We have used only general purpose features that are meaningful for all datasets addressing the relational AM task. Surprisingly, BERT embeddings (achieving state-of-the-art performances in many tasks [START_REF] Devlin | BERT: pre-training of deep bidirectional transformers for language understanding[END_REF]) do not bring improvements here, compared to non-contextualised word embeddings.

To conclude, several resources have been built recently for the task of argumentative relation prediction, covering different topics like political speeches, Wikipedia articles, persuasive essays. Given the heterogeneity of these different kinds of text, it is hard to compare cross-dataset the different proposed approaches. We addressed this nonportability issue by making a broad comparison of different deep learning methods using both non-contextualised and contextualised word embeddings for a large set of datasets for the argumentative relation prediction task, an important and still widely open problem. We proposed a set of strong dataset-independent baselines based on several neural architectures and have shown that our models perform homogeneously over all existing datasets for relation prediction in AM.

Table 1 .

 1 Datasets' statistics.

		essays micro	nk	db	ibm	com	web	cdcp	ukp	aif
	# attacks	497	108	378	141	1069	296	1301	0	5935 9854
	# supports	4841	263	353	179	1325	462	1329 1220 4759 7543

Table 2 .

 2 Experimental results for non-neural baselines with F 1 for Attack and for Support. The blanks indicate the training dataset. The Average (Avg) and the Macro (Macr) Avg exclude the results for the training datasets.

			essay micro db ibm com web cdcp ukp nk	aif	Avg Macr Avg
	non-neural baselines	RF RF SVM SVM	F 1 A 0.32 0.24 0.40 F 1 S 0.57 0.74 0.64 F 1 A 0.57 0.40 0.45 0.53 0.43 0.33 0.38 0.67 0.59 0.85 0.53 0.59 0.54 0.636 -0.39 0.55 0.44 0.381 -0.60 0.52 0.57 0.509 F 1 S 0.44 0.47 0.52 0.41 0.57 0.51 0.45 0.50 0.38 0.472 F 1 A 0.34 0.36 0.33 0.29 0.38 -0.42 0.42 0.40 0.368 F 1 S 0.71 0.67 0.65 0.67 0.59 0.84 0.57 0.56 0.49 0.639 F 1 A 0.49 0.35 0.39 0.39 0.38 -0.56 0.57 0.520 0.456 F 1 S 0.50 0.54 0.52 0.59 0.60 0.67 0.46 0.47 0.500 0.539	0.508 0.490 0.503 0.498

Table 3 .

 3 Neural baselines with contextualised word embeddings. Contextualised word embeddings such as the Bidirectional Encoder Representations from Transformers (BERT) em-Experimental results with F 1 for Attack and for Support for the Concat, Mix, and Attention architectures, with GloVE embeddings. The blanks indicate the training datasets. The Average (Avg) and the Macro (Macr) Avg do not include the results for the training datasets.

				essay micro db ibm com web cdcp ukp nk	aif	Avg Macr Avg
	embed. + syntactic	C G A G A G	F 1 A F 1 S F 1 A F 1 S F 1 A F 1 S	0.35 0.43 0.48 0.31 0.45 0.71 0.68 0.58 0.70 0.54 0.77 0.47 -0.58 0.37 0.58 0.53 0.53 -0.61 0.59 0.55 0.537 0.43 0.433 0.50 0.619 0.61 0.60 0.42 0.50 0.72 0.43 0.38 0.47 0.516 0.36 0.48 0.43 0.39 -0.52 0.45 0.51 0.449 0.75 0.66 0.62 0.68 0.79 0.52 0.56 0.54 0.640	0.526 0.526 0.544
	all features	M G A G A G	F 1 A F 1 S F 1 A F 1 S F 1 A F 1 S	0.37 0.43 0.43 0.40 0.46 0.71 0.64 0.61 0.70 0.55 0.78 0.11 0.78 0.51 0.599 -0.71 -0.46 0.466 0.36 0.54 0.50 0.51 -0.59 0.59 0.55 0.520 0.67 0.63 0.49 0.51 0.74 0.47 0.41 0.49 0.551 0.43 0.54 0.49 0.46 -0.59 0.63 0.63 0.539 0.68 0.55 0.57 0.56 0.65 0.46 0.38 0.47 0.540	0.532 0.535 0.539
				essay micro db ibm com web cdcp ukp nk	aif	Avg Macr Avg
	BERT + syntactic	4B F 1 A 2D F 1 S 4B F 1 A 1D F 1 S 4B F 1 A 0.50 0.36 0.46 0.55 0.47 0.53 0.50 0.61 0.57 0.59 0.49 0.69 0.47 0.48 0.46 0.545 -0.56 0.48 0.45 0.506 0.36 0.48 0.40 0.45 0.42 -0.53 0.37 0.430 0.69 0.67 0.61 0.62 0.57 0.79 0.50 0.50 0.619 0.50 -0.52 0.47 0.50 0.473 2D F 1 S 0.61 0.62 0.59 0.61 0.74 0.52 0.50 0.61 0.600	0.526 0.525 0.537
	all features	4B F 1 A 1D F 1 S 3B F 1 A 0.48 0.34 0.48 0.53 0.50 0.54 0.51 0.59 0.56 0.55 0.47 0.67 0.45 0.48 0.49 0.533 -0.59 0.51 0.49 0.524 0.45 -0.45 0.50 0.54 0.463 2D F 1 S 0.57 0.65 0.60 0.64 0.73 0.55 0.52 0.54 0.600	0.529 0.532
	Table

For more details about the individual datasets, we refer the reader to the relevant publications.

https://allennlp.org

We also experimented with 1 and 2 hidden layers, and hidden layer sizes of 32, 64, 128, and 256, trying all possible combinations towards best configurations. We did not consider a higher number of hidden layers due to the small size of the data.

O. Cocarascu et al. / Dataset Independent Baselines for Relation Prediction in Argument Mining


Marta Nuñez-Garcia

Nicolas Cedilnik

Shuman Jia

Hubert Cochet

Marco Lorenzi

Maxime Sermesant

Estimation of imaging biomarker's progression in post-infarct patients using cross-sectional data

Keywords: Disease progression modeling, Cross-sectional data, Ventricular arrhythmia, Post-infarct cardiac remodeling

Many uncertainties remain about the relation between postinfarct scars and ventricular arrhythmia. Most post-infarct patients suffer scar-related arrhythmia several years after the infarct event suggesting that scar remodeling is a process that might require years until the affected tissue becomes arrhythmogenic. In clinical practice, a simple timebased rule is often used to assess risk and stratify patients. In other cases, left ventricular ejection fraction (LVEF) impairment is also taken into account but it is known to be suboptimal. More information is needed to better stratify patients and prescribe appropriate individualized treatments. In this paper we propose to use probabilistic disease progression modeling to obtain an image-based data-driven description of the infarct maturation process. Our approach includes monotonic constraints in order to impose a regular behaviour on the biomarkers' trajectories. 49 post-MI patients underwent Computed Tomography (CT) and Late Gadolinium Enhanced Cardiac Magnetic Resonance (LGE-CMR) scans. Image-derived biomarkers were computed such as LVEF, LGE-CMR scar volume, fat volume, and size of areas with a different degree of left ventricular wall narrowing, from moderate to severe. We show that the model is able to estimate a plausible progression of post-infarct scar maturation. According to our results there is a progressive thinning process observable only with CT imaging; intramural fat appears in a late stage; LGE-CMR scar volume almost does not change and LVEF slightly changes during the scar maturation process.

Introduction

Scar-related substrate after myocardial infarction (MI) induces most lifethreatening ventricular arrhythmias [START_REF] Martin | Scar Tissue: Never Too Old to Remodel[END_REF]. Common recommendations for primary prevention of ventricular arrhythmia and sudden cardiac death advise cardioverter defibrillator (ICD) implantation in post-MI patients with left ventricular ejection fraction (LVEF) lower than 35% and symptoms of heart failure [START_REF] Al-Khatib | AHA/ACC/HRS guideline for management of patients with ventricular arrhythmias and the prevention of sudden cardiac death: Executive summary: A report of the American College of Cardiology/American Heart Association Task Force on Clinical Practice Guidelines and the Heart Rhythm Society[END_REF].

However, a risk stratification rule solely based on LVEF lacks sensitivity and specificity [START_REF] Goldberger | American Heart Association/American College of Cardiology Foundation/Heart Rhythm Society scientific statement on noninvasive risk stratification techniques for identifying patients at risk for sudden cardiac death: a scientific statement from the American Heart Association Council on Clinical Cardiology Committee on Electrocardiography and Arrhythmias and Council on Epidemiology and Prevention[END_REF]. Additionally, several studies have reported long periods between MIs and arrhythmia [START_REF] Martin | Characteristics of scar-related ventricular tachycardia circuits using ultra-high-density mapping: a multi-center study[END_REF] suggesting scar remodeling or maturation is a dynamic process, and it may take many years until the affected tissue becomes arrhythmogenic. Therefore, a deeper investigation of the nature of the arrhythmogenic substrate and its potential evolution over time is required. New biomarkers and more information about the right time point to evaluate them are needed.

Non-invasive imaging techniques have the potential to improve our knowledge about the nature of scar-related arrhythmogenic substrate and its dynamic remodeling process. Late Gadolinium Enhancement Cardiovascular Magnetic Resonance (LGE-CMR) imaging is the reference method to classify myocardial tissue. Many methods have been proposed to detect and quantify scar in the left ventricular (LV) wall using LGE-CMR images, most of them relying on a threshold-based rule. After segmenting the LV wall, tissue is typically divided into healthy tissue (darkest voxels), dense scar (brightest voxels) and border zone (image intensities between the other 2 thresholds). Scar remodeling using longitudinal LGE-CMR data have been recently described by Jáuregui et al. [START_REF] Jáuregui | Follow-up after myocardial infarction to explore the stability of arrhythmogenic substrate: the FOOT-PRINT study[END_REF]. After scanning post-MI patients at 7 days, 6 months and 4 years after the infarction, the authors showed that dense scar and border zone mass decreased over time suggesting the existence of a long-term scar healing process. bayesi Due to the relatively low spatial resolution of LGE-CMR images and several contraindications (patients carrying ICDs, claustrophia, etc.) Computed Tomography (CT) imaging has been recently introduced to accurately characterize scar. Furthermore, several studies have shown the existence of a progressive thinning of the LV wall in scarred areas [START_REF] Anversa | Cellular basis of ventricular remodeling after myocardial infarction[END_REF]. Wall thickness and intramural fat mapping from CT are promising imaging biomarkers that could be useful to improve risk stratification in patients considered for ICD implantation and to better personalize ablation therapies.

Disease progression modeling (DPM) provides a data-driven description of the natural evolution of a given pathology and it aims at revealing long-term pathological trajectories from short-term clinical data [START_REF] Lorenzi | Probabilistic disease progression modeling to characterize diagnostic uncertainty: application to staging and prediction in Alzheimer's disease[END_REF]. Reformulating DPM within a probabilistic setting has the potential of predicting a plausible evolution of the different biomarkers considered.

In this paper, we used CT and LGE-CMR imaging biomarkers from crosssectional (i.e. one single time point per patient) post-MI data, and a probabilistic DPM to characterize scar maturation. The model shows a plausible evolution of the different biomarkers taken into account and allows us to estimate a relative temporal timeline of disease progression.

Methods

Data processing

Forty-nine post-MI patients underwent CT and LGE-CMR imaging at the Centre Hospitalier Universitaire (CHU) de Bordeaux (France). A schematic representation of the proposed method can be seen in Figure 1. The LV wall was manually segmented from the LGE-CMR images by experts. As typically done in clinical practice, a threshold-based approach was used to classify voxels into healthy tissue (image intensity below the 40% of the maximum intensity in the LV wall), border zone (between 40% and 60% of the maximum intensity) and dense scar (above 60% of the maximum intensity).

The LV wall was automatically segmented from CT images using a deep learning approach based on two successive U-net networks [START_REF] Cedilnik | Fully Automated Electrophysiological Model Personalisation Framework from CT Imaging[END_REF][START_REF] Jia | Automatically segmenting the left atrium from cardiac images using successive 3D U-nets and a contour loss[END_REF][START_REF] Ronneberger | U-net: Convolutional networks for biomedical image segmentation[END_REF]. The high spatial resolution of CT images requires high memory resource while at the same time the ventricles take only a fraction of the entire CT volume. For that reason, we input to the first network a low-resolution version of the CT data and the output (coarse segmentation) was used to locate the ventricles. Only the region around them was kept and a cropped high-resolution version of the image was fed to the second 3D U-net. The output masks were post-processed (including up-sampling to the original CT image resolution) to obtain clean and non-overlapping masks. The model was trained using 450 CT scans with available expert segmentations of the LV endocardium and epicardium, and the right ventricular epicardium. 50 CT scans were used for validation with a loss function defined as the opposite of the Dice score.

To estimate the LV wall thickness a previously described method based on solving a partial differential equation using the endocardium and epicardium masks was used [START_REF] Yezzi | An Eulerian PDE approach for computing tissue thickness[END_REF]. This method assigns a thickness value to each voxel of the LV wall. Wall thickness measurements were then projected onto a mid-wall 3D surface mesh obtained from the LV wall segmentation where areas of different thickness bins were computed. We considered 4 wall thickness bins representing different degrees of wall thinning, from moderate to severe: thickness between 5

Disease progression modeling

The statistical framework used in this study (detailed in [START_REF] Lorenzi | Probabilistic disease progression modeling to characterize diagnostic uncertainty: application to staging and prediction in Alzheimer's disease[END_REF]) formulates disease progression modeling based on Bayesian Gaussian Process (GP) regression [START_REF] Rasmussen | Gaussian processes in machine learning[END_REF] by modeling individual time transformations encoding the information of the associated pathological stage, and introducing monotonicity constraints to impose a plausible behaviour on the biomarkers' trajectories from normal to pathological stages. Briefly, biomarkers' evolution are modeled as monotonic GPs, while individuals are assigned to their specific time point relative to the regression time axis, according to the severity of the associated biomarker measurements. The time axis corresponds to the patient's scar age, i.e. delay between the infarction event and the imaging study. Monotonic constraints were defined as follows: increasing for fat volume and all degrees of wall thinning [START_REF] Anversa | Cellular basis of ventricular remodeling after myocardial infarction[END_REF], and decreasing for LVEF, dense scar volume, and border zone volume [START_REF] Jáuregui | Follow-up after myocardial infarction to explore the stability of arrhythmogenic substrate: the FOOT-PRINT study[END_REF].

Results

The predicted biomarkers temporal trajectories are shown in Figure 2. It can be observed that the model predicts a plausible progressive thinning of the LV wall and an increase of fat volume in a late stage. On the contrary, the trajectories corresponding to border zone volume, dense scar volume, and LVEF remain almost constant.

All biomarkers trajectories can be seen together in Figure 3 (left). Figure 4 and Table 2 show the estimated distributions of the maximum change time for all biomarkers suggesting when the biggest change for each biomarker occurs. These parameters allow us to estimate the timeline for remodeling changes by inferring which biomarkers may change before the others. The progressive wall thinning as seen by CT thickness seems clear showing, moreover, an initial fast thinning followed by a more softened, spaced in time, narrowing. According to the model, fat appears only in a late stage. The high variance of dense scar and border zone volumes indicates that no clear evolution pattern is seen for these features.

Importantly, due to the lack of longitudinal data (i.e. only one time point per patient was available), the time axis can only be interpreted as relative time and it does not represent real years. We believe however that this methodology can be useful to describe the dynamics across biomarkers, while distinguishing the most informative features to state the pathological stage of an individual during the course of the scar maturation process. Additionally, our method may enable to position a given patient among other patients at a similar state on the scar maturation process. It may be crucial to determine when the scarred tissue starts to become severely thin (thickness < 2mm) because the total area of severe thinning has been shown to be a good predictor for ventricular arrhythmia. In this cohort the area of severe thinning was significantly greater in patients with arrhythmia (17 patients), p = 0.0019 in a two-sample t-test. Other significant features (at the 5% level) were: the area of thickness between 3 and 2 mm (p = 0.0028); area of thickness between 4 and 3 mm (p = 0.037); dense scar volume (p = 0.010); LVEF (p = 0.028); and fat volume (p = 0.035). Progression time course

Conclusions

This work suggests that probabilistic disease progression modeling can be used to estimate the long-term progression of imaging biomarkers related to post-MI scar maturation using purely cross-sectional data. According to the model, it exists a progressive LV wall thinning process in post-MI patients that lasts for years. CT may outperform LGE-CMR imaging in characterizing post-infarct scar remodeling as well as in predicting disease stage. This prediction may be useful for example to assess individual cardiac risk to advise ICD insertion (currently determined by the time passed since the infarction, and suboptimal metrics such as LVEF). It can also be useful to adapt treatments (e.g. ablation therapy) to the patient's stage according to the model.

Fig. 1 :

 1 Fig. 1: Data processing pipeline. LGE-CMR images are manually segmented and the LV wall tissue is classified using a threshold-based approach. CT images are automatically segmented using 2 successive U-nets, and wall thickness is computed afterwards. Then, 8 imaging biomarkers are calculated and used to estimate the probabilistic disease progression model. LGE-CMR = Late Gadolinium Enhanced Cardiac Magnetic Resonance; CT = Computed Tomography; LVEF = Left Ventricular Ejection Fraction.

Fig. 2 :

 2 Fig. 2: Modeled biomarker progressions. Each red line displays the estimated trajectory for a given patient and a given biomarker. Dots correspond to individual samples colored according to disease status: not VA (blue) and VA (orange). Note that since only one time point per patient was available the time axis can only be interpreted as relative time and it does not represent real years. LVEF = Left Ventricular Ejection Fraction; VA = Ventricular arrhythmia.

Fig. 3 :

 3 Fig. 3: On the left, all estimated biomarker trajectories; on the right, corresponding derivatives showing the temporal point of maximum change (i.e. point of maximum slope in the trajectories) and the relative magnitude of the change.

Fig. 4 :

 4 Fig. 4: Distributions of the biomarkers' maximum change time.

Table 1 :

 1 Population characteristics. Values are mean ± STD [Min, Max]. STD = Standard deviation. LAD = Left Anterior Descending; LCX = Left Circumflex; RCA = Right Coronary Artery.

	Parameter	Mean ± STD [Min, Max]
	Patient's age (years)	65.53 ± 11.62 [38, 91]
	Scar age (years)	7.39 ± 9.06 [0.3, 40]
	LV end-diastolic volume (mL)	137.47 ± 38.00 [90.61, 242.14]
	LV end-diastolic area (cm 2 )	181.10 ± 32.33 [123.40, 250.69]
	LVEF (%)	44.04 ± 12.96 [23, 75]
	Dense scar volume (mL)	8.24 ± 8.44 [0.72, 50.07]
	Border zone volume (mL)	13.95 ± 8.68 [2.45, 34.69]
	Fat volume (mL)	0.69 ± 1.66 [0, 9.13]
	Area of wall thickness between 5mm -4mm (cm 2 ) 17.81 ± 11.07 [1.62, 45.44]
	Area of wall thickness between 4mm -3mm (cm 2 )	10.69 ± 7.44 [1.37, 27.99]
	Area of wall thickness between 3mm -2mm (cm 2 )	7.51 ± 6.37 [1.56, 25.44]
	Area of wall thickness below 2mm (cm 2 )	7.16 ± 4.07 [3.59, 23.43]
	Gender	5 (10.20%) female patients
	LAD scar	20 (40.81%) patients
	LCX scar	17 (34.69%) patients
	RCA scar	23 (46.94%) patients
	and 4 mm; thickness between 4 and 3 mm; thickness between 3 and 2 mm; and
	thickness below 2mm.
	Eight biomarkers were included in the analysis. From LGE-CMR: LVEF;
	volume of dense scar; and volume of border zone. From CT: fat volume; and areas
	of wall thinning considering the 4 thickness bins mentioned before. Population
	characteristics are shown in Table 1 including scar age (i.e. lapse of time between
	MI and imaging studies).

Table 2 :

 2 Distributions of the maximum change time in ascending order. STD = Standard deviation.

	Biomarker	Time (mean ± STD)
	4mm <thickness <5mm	2.1335 ± 4.4195
	3mm <thickness <4mm	4.06589 ± 1.7665
	2mm <thickness <3mm	16.3780 ± 3.5544
	LVEF	18.06976 ± 15.06029
	Dense scar volume	23.2351 ± 21.7677
	thickness <2mm	25.9893 ± 3.6027
	Border zone volume	30.9844 ± 23.8233
	Fat volume	35.7065 ± 5.5985

Acknowledgements. Part of this work was funded by the ERC starting grant EC-STATIC (715093), the IHU LIRYC (ANR-10-IAHU-04), the Equipex MUSIC (ANR-11-EQPX-0030) and the ANR ERACoSysMed SysAFib projects. This work was also supported by the French government, through the 3IA Côte d'Azur Investments in the Future project managed by the National Research Agency (ANR) with the reference number ANR-19-P3IA-0002, and by the ANR JCJC project Fed-BioMed 19-CE45-0006-01. We would like to thank all patients who agreed to make available their clinical data for research.


Evaggelos Kritsikis
email: evaggelos.kritsikis@grenoble.cnrs.fr

Jean-Christophe Toussaint

Olivier Fruchart

Fast computation of magnetostatic fields by Non-uniform Fast Fourier Transforms

The bottleneck of micromagnetic simulations is the computation of the long-ranged magnetostatic fields. This can be tackled on regular N -node grids with Fast Fourier Transforms in time N log N , whereas the geometrically more versatile finite element methods (FEM) are bounded to N 4/3 in the best case. We report the implementation of a Non-uniform Fast Fourier Transform algorithm which brings a N log N convergence to FEM, with no loss of accuracy in the results.

The power of computers steadily increases over the years while the size of devices used in fundamental science or technology is shrinking. Today we have reached a cross-over where numerical simulations are capable of describing in detail the physics of nanodevices, for which they thus play a leading role in their understanding and designing. In numerical micromagnetics for spin electronics, a bottleneck is the computation of the magnetostatic interactions which by nature are long-ranged. These interactions can be expressed in terms of either magnetostatic field H or scalar pseudo-potential φ such that H = -∇φ. The latter is convenient since it boils the problem down to a single scalar unknown. To deal with magnetostatic interactions, essentially two distinct approaches have been implemented, depending on the type of mesh used :

1. for Finite Difference (ie. translation invariant) meshes, a Green approach with Fast Fourier Transforms, called FD-FFT. Computation time is moderate (N log N with N the number of nodes), but the curved boundaries that occur often in experimental devices are not ideally described ;

2. for Finite Element (much more general) meshes, a Finite Element Method coupled to a Boundary Element Method, called FEM-BEM. This can faithfully describe curved boundaries but computation time is higher, at least N 3/2 in 2D (resp. N 4/3 in 3D).

In this Letter we report the implementation of a new magnetostatic code which combines the advantages of both cited approaches : it uses a FEM mesh thus describing curved boundaries as well as FEM-BEM, albeit with computation time N log N . It is based on a algorithm reported recently for computing non-periodic Fast-Fourier Transforms (NFFT) 1 . Our code, called FEM-NFFT, proves to be significantly faster than FEM-BEM with no loss of precision. As a first step the implementation was done for a 2D geometry (which pertains to 3D systems with one direction of translational invariance, i.e. cylinder-like) for the proof of concept. The gain is expected to be even greater in more realistic 3D calculations. This demonstrates the potential of FEM-NFFT for micromagnetism and thus spin-electronic devices.

Let us recall the principle and features of FD-FFT and FEM-BEM before presenting our approach and results. We consider a system Ω with boundary ∂Ω, displaying a known magnetization distribution M(r).

Finite Difference micromagnetic codes use a translation invariant grid. On such a grid, FFTs can be used to compute convolutions in time N log N . This motivates a Green approach for magnetostatics : φ is calculated as a convolution of the Green function

G = -(1/2π) log r in 2D [resp. G = 1/(4πr) in 3D] with the magnetic charges, volumic ρ = -∇ • M and surfacic σ = M • n : φ(r) = Ω ρ(r ′ )•G(r-r ′ )dr ′ + ∂Ω σ(r ′ )•G(r-r ′ )dr ′ (1)
The N log N speed explains the wide and lasting use of FD in micromagnetic simulations 2 . However, most devices have curved boundaries, either by design or as a result of experimental imperfections. Simulating these cases with FD requires the use of saw-tooth boundaries to describe the magnetic material. This geometrical approximation may induce inadequate descriptions 3 . Finite Element micromagnetic codes use, on the contrary, complex-shaped meshes with triangles (resp. tetrahedrons) as 2D (resp. 3D) unit cells. They consequently suffer much less from the above-mentioned limitations. However, without translational invariance of the mesh, FFTs are not available. Bearing in mind that direct summation of Eq.(1) on a FEM mesh, called FEM-direct, would cost N 2 time, one understands why the Green approach is thought incompatible with FEM codes.

To deal with magnetostatics, FEM codes thus go back to the Poisson equation -∆φ = ρ, with the usual regularity condition that the field should decay at infinity. Applying FEM, the equations are translated into a linear system, which is solved by standard iterative methods [START_REF] Braess | Finite Elements[END_REF] in time N 3/2 in 2D (resp. N 4/3 in 3D) Not only is this asymptotically slower than the N log N time required for FD-FFT, but N here takes higher values, because the mesh must extend well beyond Ω in order to tackle the regularity condition at infinity. This induces an additional slowdown, and also creates finite-size artifacts.

To avoid meshing outside Ω, the main approach consists in coupling FEM with a Boundary Element Method, resulting in the so-called FEM-BEM [START_REF] Fredkin | [END_REF] . The asymptotic complexity of the Poisson solver is unchanged but the BEM step introduces another time limitation N 2 ∂ where N ∂ is the number of boundary nodes. In the most favorable case consisting of compact systems N ∂ ≈ N 1/2 in 2D (resp. N ∂ ≈ N 2/3 in 3D). However for flat geometries, of particular relevance to applications, N ∂ ≈ N , in which case the time limitation for FEM-BEM may be pretty severe.

Our innovation is to revert, within the FEM framework, to a Green approach. The convolution (1) is discretized in a way typical for FEM, and computed using a fairly recent mathematical method called NFFT (Nonuniform Fast Fourier Transform). NFFTs allows one to compute discrete convolutions in time N log N without the equispaced data requirement of FFTs.

More in details, we seek φ at the nodes (r i ), i=1...N of the mesh. A linear interpolation inside each element, of known magnetization values M(r i ), is used to evaluate charges ρ and σ at points r j , j=1...M defined as the quadrature points for the integrals in (1) 6 .

Consequently, (1) is rewritten the following way :

φ(r i ) = M j=1 ρ j G(r i -r j ) ω j det J(r j ) + M j=1 σ j G(r i -r j ) ω j det J(r j ) (2)
where -ρ j = ρ(r j ) if r j is in the interior of Ω, 0 otherwise, -σ j = σ(r j ) if r j is on ∂Ω, 0 otherwise, -ω j is the weight of r j in the quadrature scheme, -J(r j ) is the Jacobian of the affine transformation mapping the unit element on that containing r j .

We then use NFFTs to compute (2). Although the seminal paper 7 dates back to 1993, the NFFT method remains little-known even in the mathematical community. A presentation of the method can be found Ref. 8 . Here we sketch the basic strategy and show that computation time does not exceed N log N , ie. that of the classical FFT.

Let us look at our non equispaced data as a sum of Dirac functions. The goal is to find the spectrum of this data function. The idea is to convey initial information over a regular grid, so that an FFT can be used. We therefore choose a regular grid X and let the data diffuse to the X nodes through convolution with a Gaussian function (or more generally a smooth localized function). The Gaussian is localized in space, so we can consider that each piece of data diffuses only to a fixed number of the nearest X nodes. Computation time of this diffusion step is thus proportional to N .

The question arises how to choose the period of the X grid. In our case, the answer depends on how smooth the Green function G is. A necessary preliminary step before executing the NFFT is therefore to smooth G around the origin ; the price to pay is an afterwards correction in the smoothing zone. It can be shown that a grid of size p 2 N in 2D (resp. p 3 N in 3D) is convenient, where p is the degree of smoothness chosen for G.

We then perform, according to the initial idea, a FFT on the X grid, in time proportional to N log N . Based on the convolution theorem, what we get is the Fourier coefficients of the data function multiplied by those of the Gaussian. Therefore, we finally divide these numbers by the Fourier coefficients of the Gaussian to get the desired spectrum. The number of divisions is proportional to N . As a whole, the NFFT is expected to behave asymptotically like N log N , as all extra steps behave like N .

We implemented FEM-NFFT to 2D test cases where an analytical solution φ a is available, so that errors can be readily estimated. Interpolation and quadrature routines are written in C++ and the NFFT package used 9 is in C99. For G we have chosen a smoothness degree of 2. On each test case, we provide computation times and error estimates for FEM-direct, the classic FEM-BEM and our NFFT-based method. The computed error is the normalized root mean square (M s L) -1 ( N i=1 |φ(r i )φ a (r i )| 2 /N ) 1/2 , where M s is the saturation magnetization and L the system diameter set at unity, . Computations are done on an Intel P4 2GHz with 1GB RAM running Fedora 5.

The first test case is a disk uniformly magnetized along the x-axis (a cylinder in 3D space). The analytical solution is φ(x) = M s x/2. The second case is the so-called magic cylinder, a circular annulus of radii R 1 , R 2 where the angle between magnetization and the x-axis equals twice the polar angle. The name stems from the uniform magnetic field thus induced in the inner region. The analytical solution inside the annulus is φ(r, θ) = M s r cos θ log(r/R 2 ) in polar coordinates.

Tables 1 and2 display the numerical results for the two cases, respectively. It can be readily seen that FEM-NFFT provides results very similar to FEM-direct. The error induced by the NFFTs is thus negligible. Compared to FEM-BEM, errors are comparable for non-uniformly magnetized systems (see Table 2) whereas for uniform distributions, Green approaches, to which FEM-NFFT belong, are more accurate by one order of magnitude (see Table 1). This is because they can treate apart volu-mic and surfacic charge contributions. Concerning computational time, FEM-direct is as expected the quickest, gaining a factor around 5 over FEM-BEM for the finer meshes. Nodes required in 3D cases of interest commonly count up to 10 5 , around which number the time advantage of FEM-NFFT over classical methods such as FEM-BEM is expected to reach one order of magnitude. To conclude, we have successfully implemented a Nonuniform Fast Fourier Transform (NFFT) algorithm to compute magnetostatic fields for micromagnetic simulations based on Finite Element methods (FEM). The new approach, called FEM-NFFT, combines the advantages previously found separately in Finite Difference methods Tab. II -Test case 2 : error (in ppm) and computation time (in seconds) of FEM-direct, FEM-BEM and FEM-NFFT for different mesh sizes. (computation time scaling like N log N ) and FEM (faithful description of curved boundaries). Thus FEM-NFFT promises a leap in the attractiveness of micromagnetic simulations of spin electronic devices.

Fig. 1 -

 1 Fig. 1 -Test case 1 : the uniformly magnetized disk.Mesh used for N = 400 (left) and magnetization distribution (right).

Fig. 2 -

 2 Fig. 2 -Test case 2 : the magic cylinder. Mesh used for N = 1192 (left) and magnetization distribution (right).

G. S. D. Potts and M. Tasche, Fast Fourier transforms for nonequispaced data : A tutorial (Birkhäuser, Boston,

2001), pp. 247-270.2 B. Kevorkian, Ph.D. thesis, Université Joseph Fourier, Grenoble (1998). 3 Z. G. C. J. García-Cervera and E. Weinan, J. Comput.


Jean-Claude Tourneur
email: jc.tourneur@nb.sympatico.ca

Claire Cole

Joël Meunier
email: joel.meunier@univ-tours.fr

The first of five moults of Forficula auricularia L. (Dermaptera: Forficulidae)

Most studies about the European earwig Forficula auricularia L. (Dermaptera: Forficulidae) assume that four larval moults occur from hatching to adult. In this note, we use observational and quantitative approaches to reveal the presence of an additional intermediate moult occurring immediately after egg hatching. We then briefly discuss how this finding may open up new possibilities for follow-up studies to improve our general understanding of the taxonomy, biology, and evolution of this species.

All insects pass through successive moulting events to reach adulthood. These moults typically allow general body growth and the development (or loss) of specific appendix and organs that become essential (or useless) at a certain period of their life [START_REF] Chapman | Postembryonic development[END_REF]. This is the case, for instance, of the genitalia and sexual ornaments that often appear only after the last moult into adulthood [START_REF] Belles | Origin and evolution of insect metamorphosis[END_REF][START_REF] Rolff | Complete metamorphosis of insects[END_REF]. The number and temporality of moulting events during insect development are under strong selection pressure, because they typically reflect age-specific physiological constraints and newly moulted individuals are often more vulnerable to pathogens and predation due to their resulting soft cuticle [START_REF] Field | Infection of the Norway lobster Nephrops norvegicus by a Hematodinium-like species of dinoflagellate on the west coast of Scotland[END_REF][START_REF] Corteel | Molt stage and cuticle damage influence white spot syndrome virus immersion infection in penaeid shrimp[END_REF].

Somewhat surprisingly, the number of moulting events that occur during the development of the European earwig Forficula auricularia L. (Dermaptera: Forficulidae) remains unclear. This insect is present worldwide and exhibits a relatively long lifespan, with males and females living up to two years. Females typically produce a clutch of 30 to 60 eggs in early winter [START_REF] Tourneur | Egg laying in a northeastern North American (Montréal, Québec) population of Forficula auricularia L. (Dermaptera: Forficulidae)[END_REF][START_REF] Meunier | One clutch or two clutches? Fitness correlates of coexisting alternative female life-histories in the European earwig[END_REF][START_REF] Tourneur | Variations in seasonal (not mean) temperatures drive rapid adaptations to novel environments at a continent scale[END_REF], provide them with extensive forms of care through winter, and then remain with the newly hatched juveniles for several weeks during which they continue providing post-hatching care [START_REF] Lamb | Parental behavior in the dermaptera with special reference to Forficula auricularia (Dermaptera: Forficulidae)[END_REF][START_REF] Kölliker | Benefits and costs of earwig (Forficula auricularia) family life[END_REF][START_REF] Koch | Mother and offspring fitness in an insect with maternal care: phenotypic trade-offs between egg number, egg mass and egg care[END_REF][START_REF] Diehl | Surrounding pathogens shape maternal egg care but not egg production in the European earwig[END_REF][START_REF] Van Meyel | Love them all: mothers provide care to foreign eggs in the European earwig Forficula auricularia[END_REF]. Whereas it is generally assumed that the number of larval stages present in this species is four (L1 to L4) and thus that larva pass through four moulting events during their lifetime (L1 to L2, L2 to L3, L3 to L4, and then L4 to adults; [START_REF] Jones | The European earwig and its control[END_REF][START_REF] Fulton | The European earwig[END_REF][START_REF] Crumb | The European earwig[END_REF][START_REF] Lhoste | Les stades larvaires et la division des articles antennaires chez Forficula auricularia L. (Dermaptera.)[END_REF][START_REF] Henson | On the malpighian tubules of Forficula auricularia (Dermaptera)[END_REF][START_REF] Henson | The growth and form of the head and antennae in the earwig (F. auricularia Linn[END_REF][START_REF] Behura | The biology of the European earwig, Forficula auricularia Linn[END_REF]), the literature contains contrasting support for this claim. In particular, Chapman (1917a) mentioned the presence of six "ecdyses" while reporting only four moults in his subsequent paper (Chapman 1917b). Moreover, the occurrence of a moult immediately after egg hatching (a process commonly known as the intermediate moult; [START_REF] Chapman | Postembryonic development[END_REF] has been reported in a few other dermapteran species [START_REF] Albouy | Dermaptères ou perce-oreilles[END_REF]) such as Marava arachidis (Yersin) (Dermaptera: Spongiphoridae), where the first developmental stage only lasts the time of hatching [START_REF] Herter | Zur fortpflanzungsbiologie eines lebendgebärenden ohrwurmes (prolabia arachidis yersin)[END_REF]). Finally, a few studies report that F. auricularia larva are already well developed at the end of egg development [START_REF] Fulton | The European earwig[END_REF][START_REF] Behura | The biology of the European earwig, Forficula auricularia Linn[END_REF], suggesting that an intermediate moult could be possible just before or immediately after egg hatching [START_REF] Chapman | The insects: Structure and function[END_REF]. Whether this moult indeed occurs in F. auricularia remained unknown [START_REF] Lhoste | Les stades larvaires et la division des articles antennaires chez Forficula auricularia L. (Dermaptera.)[END_REF].

In this study, we used observational and quantitative approaches to demonstrate that an intermediate moult occurs immediately after egg hatching in F. auricularia. Our observations and measurements were based on 33 F. auricularia families originating from earwig females field sampled in August 2019 in Ripon,Québec,Canada (75° 14' 77" W,45° 81' 06" N,altitude 195 m) and then maintained under standard conditions [START_REF] Meunier | One clutch or two clutches? Fitness correlates of coexisting alternative female life-histories in the European earwig[END_REF]. Of the sister species A and B currently included under the species name F. auricularia [START_REF] Wirth | Molecular and reproductive characterization of sibling species in the European earwig (Forficula auricularia)[END_REF], the population studied here belongs to subspecies A [START_REF] Tourneur | Factors affecting the egg-laying pattern of Forficula auricularia (Dermaptera: Forficulidae) in three climatologically different zones of North America[END_REF]. From field sampling to the day of oviposition, these females were fed with fresh carrots (Apiaceae) and an artificial diet (dry power of egg yolk 1/3, bee-collected pollen 1/3, and cricket 1/3) weekly. At oviposition, each female was isolated in a new Petri dish (diameter 10 cm) with a thin layer of humid sand. At that time, males and food sources were also removed to prevent egg cannibalism by males and because females stop foraging activity during the period of egg care [START_REF] Kölliker | Benefits and costs of earwig (Forficula auricularia) family life[END_REF]. After hatching, the larvae were maintained under the same standard conditions to conduct morphological measurements until they reached adulthood (see below).

Our observations of the intermediate moult began at the end of egg development. At that time, the eyes, mandibles, and ruptor ovi of the larvae present in the eggs are already visible (Fig. 1). The ruptor ovi (also called egg burster) is a tooth-like projection that an insect embryo or larva typically uses to break the embryonic envelope and egg chorion when hatching [START_REF] Cooper | Ruptor ovi, the number of moults in development, and method of exit from masoned nests[END_REF][START_REF] Cai | Morphology of the immature stages of Panorpa qinlingensis (Mecoptera: Panorpidae) with notes on its biology[END_REF]. In F. auricularia, the ruptor ovi is located on the larval frons (Fig. 1A). We observed that this first moult begins with the rupture of the frontal and coronal sutures and occurs immediately after egg hatching. During this phase, the larva swallows a large quantity of air, which greatly distends its abdomen and thus facilitates the exuviation. As the exuviation takes place, the ruptor ovi initially present at the front position (of the larval stage thereafter called L0) then slides between the next stage's mandibles to the base of its abdomen (Fig. 2A to F). We propose naming the first stage "L0", both because such a pharate first-stage larva is not always considered as a larval instar in insects [START_REF] Chapman | Postembryonic development[END_REF] and to remain consistent with existing literature referring to L1 to L4 larva in the European earwig (e.g., [START_REF] Jones | The European earwig and its control[END_REF][START_REF] Fulton | The European earwig[END_REF][START_REF] Crumb | The European earwig[END_REF][START_REF] Lhoste | Les stades larvaires et la division des articles antennaires chez Forficula auricularia L. (Dermaptera.)[END_REF][START_REF] Henson | On the malpighian tubules of Forficula auricularia (Dermaptera)[END_REF][START_REF] Henson | The growth and form of the head and antennae in the earwig (F. auricularia Linn[END_REF][START_REF] Behura | The biology of the European earwig, Forficula auricularia Linn[END_REF][START_REF] Vogelweith | Age, pathogen exposure, but not maternal care shape offspring immunity in an insect with facultative family life[END_REF].

The fact that the exuvia remains stuck at the base of the chorion facilitates the exuviation of the antennae, legs, and abdominal cerci. Once the larva is released (the frontal ovi ruptor is no longer present), it regurgitates its abdominal air and thus allows its abdomen, the intersegmental membranes of which have been distended, to return to normal size and shape. The sclerotisation of the cephalic capsule of the resulting L1 larva then begins. Our description of this additional moult in the European earwig fits with the description of intermediate moults reported in other groups of insects, such as Orthoptera and Hemiptera [START_REF] Chapman | Postembryonic development[END_REF].

To confirm that the larvae visible in the eggs just before hatching belong to a different developmental stage (L0) than the larvae resulting from egg hatching (L1), we measured and compared the interocular distance of larvae from each of the newly defined five developmental instars. These measurements were conducted to the 0.01 mm using a binocular scope (Nikon-SMZ-2T; Nikon Canada Inc., Mississauga, Ontario, Canada) on 18 to 51 larvae per instar (Fig. 3).

The eyespots are the first cephalic structures observable in the first developmental stage, and their measurement is possible through the larvae's translucent chorion (Fig. 1A). Our results confirm that the interocular distances of the measured larvae differed amongst the five developmental stages (linear model; F154,4 = 623.6, P < 0.0001) with mean values (± SE) of 0.389 ± 0.005 mm in the newly defined L0 larva (in the eggs), 0.495 ± 0.004 mm in L1, 0.563 ± 0.006 mm in L2, 0.681 ± 0.013 mm in L3, and 0.834 ± 0.005 mm in L4 (Fig. 3). Each interocular distance differed significantly from another, including between the interocular distances of L0 larvae measured through the eggshell and of L1 larva measured at egg hatching (pairwise comparisons based on the estimated marginal means of the linear model presented above, with P-values corrected for multiple testing using Tukey methods; all P < 0.0001; Fig. 3). Finally, these interocular distances (Fig. 1) fitted a straight line across the five newly defined larval stages (R 2 = 0.92; y ~ 10.81 * x + 26.56; x being the larval stage), a result in line with predictions of Dyar's Law [START_REF] Dyar | The number of molts of lepidopterous larvae[END_REF]. This law, which predicts that the head width of a larva should follow a regular geometrical progression in its successive instars, was originally formulated to define the number of moults across insects [START_REF] Dyar | The number of molts of lepidopterous larvae[END_REF][START_REF] Shull | Application of Dyar's law to life stages of sigara mathesoni (Heteroptera: Corixidae)[END_REF] and therefore confirms the presence of an intermediate moult in F. auricularia.

Overall, our data reveal the presence of an intermediate moult in F. auricularia (and thus a total of five, not four, moults, as commonly stated in the literature), which occurs almost immediately after egg hatching. This finding may open new opportunities for follow-up studies to improve our general understanding of the taxonomy, biology, and evolution of this species.

For example, number of moults is an important criterion in taxonomy [START_REF] Chapman | Postembryonic development[END_REF], and these results call for future work to explore whether this new information could be used to discriminate amongst the sister species A and B currently included under the species name, F. auricularia [START_REF] Wirth | Molecular and reproductive characterization of sibling species in the European earwig (Forficula auricularia)[END_REF]. Similarly, natural variation in environmental constraints might affect the success, timing, and duration of this intermediate moult, as they do for any moulting event [START_REF] Chapman | The insects: Structure and function[END_REF]; our findings thus suggest investigations into the possible impacts of these constraints and of their potential stability in the invasion success of the European earwig worldwide [START_REF] Quarrell | The invasion biology of the invasive earwig, Forficula auricularia in Australasian ecosystems[END_REF][START_REF] Tourneur | Variations in seasonal (not mean) temperatures drive rapid adaptations to novel environments at a continent scale[END_REF]. Finally, maternal care towards the eggs is both extended and essential in this species [START_REF] Boos | Maternal care provides antifungal protection to eggs in the European earwig[END_REF][START_REF] Koch | Mother and offspring fitness in an insect with maternal care: phenotypic trade-offs between egg number, egg mass and egg care[END_REF][START_REF] Van Meyel | Love them all: mothers provide care to foreign eggs in the European earwig Forficula auricularia[END_REF]; our results prompt further tests to shed light on this behaviour's possible impact on the success of this first moulting event (e.g., Fig. 1B reports the case of a failed intermediate moult resulting in the death of a sclerotinised L1 larva in the eggs).

Fig. 2 .

 2 Fig. 2. Migration of the ruptor ovi after hatching and during first exuviation (from A to F). A circle locates the ruptor ovi in each picture.

Fig. 3 .

 3 Fig. 3. Interocular distances measured in each of the five larval developmental stages. Boxplots depict median (bold bar) and interquartile range (light bar), with whiskers extending to 1.5 times the interquartile range and dots representing experimental values. Pairwise comparisons between developmental stages were conducted using estimated marginal means, and the Pvalues were corrected for multiple testing using Tukey methods. Different letters correspond to P < 0.0001. The number of measured larvae is indicated in parentheses at the bottom of each bar. The dashed line illustrates the period of egg hatching.